Fully three-dimensional propagation model of horizontal borehole hydraulic fractures in strata under the effect of bedding planes

Abstract

The bedding plane effect will occur when hydraulic fractures propagate to the bedding plane in sedimentary strata, resulting in the “≠,” “工,” or “/” shaped fracture morphology. Based on previous physical experiments results, this article analyzed the mecroscopic propagation mechanism of tensile failure and the mechanical conditions for main hydraulic fracture and the bedding plane fracture propagating, proposing the criteria for hydraulic fracture to penetrate through the bedding plane. A fully three-dimensional model of hydraulic fracture morphology in horizontal borehole hydraulic fracturing is established with the vertical water flow, water leak-off, and bedding plane effect taken into consideration. Basic equations of continuity, pressure decline, hydraulic fracture morphology, and others are solved. After that, true triaxial hydraulic fracturing experiments with samples containing bedding planes are conducted to verify the aperture, length, width, and height of hydraulic fractures in this model. The model is proved to be accurate and reliable.

Introduction

The structural transformation of coal and rock matrix relates to many difficult problems including stimulated extraction of oil and shale gas, controlling the hard roof of coals, prevention of rockburst, stimulated extraction of gas in low penetration coal seams, and prevention of outburst of gas and coal. The morphology of hydraulic fracture is a key factor affecting the effect of hydraulic fracturing and, therefore, must be controlled (Huang, 2010; Huang et al., 2011). Reservoirs of oil, shale gas, and coalbed methane and the strata in coal mine are usually sedimentary rock with horizontal and layered structure. Visible coal–coal, rock–rock, and coal–rock interfaces, which are referred to as bedding planes, exist in the above sedimentary strata. Field tests show that interface effect will occur when the hydraulic fracture propagates to the bedding plane (Blair et al., 1989; Chen and Chen, 1997; Daneshy, 1978). The bedding plane has a significant effect on the morphology of hydraulic fractures (Chen and Chen, 1997). Study of the interface effect on hydraulic fracturing is of great importance to control hydraulic fractures propagation in specific rock strata.

Field investigation shows that bedding planes are tightly closed with low aperture, smooth wall, and low roughness. The fillers within the bedding plane are little and the interlayer spacing is large. The failure mode of bedding plane fracture is tensile failure, which belongs to the type I failure mode in fracture mechanics.

When the hydraulic fracture intersects with the bedding plane perpendicularly, the fracturing water will first flow within the bedding plane. Later the hydraulic fracture will penetrate through the bedding plane and extend along its original direction (Blair et al., 1989). Weak interface will prevent hydraulic fractures propagating along its original direction, while strong interface will not. However, the propagation law of hydraulic fractures will not change as the property of the two formations of the weak interface changes (Daneshy, 1978). Hydraulic fractures will penetrate natural fissures and propagate along its original direction only when the stress difference is high and the approaching angle is large. In most cases, hydraulic fracture will stop propagating or turn after intersecting with natural fissures (Blanton, 1982, 1986; Thiercelin et al., 1987). Quantitative researches show that the horizontal stress difference and the approaching angle have great impact on the open and the shear failure mechanics of natural fissures (Zhou et al., 2007, 2008, 2010). Hydraulic fractures can penetrate the interface directly or propagate along the bedding plane, which depends on the difference of the mechanical property of the production formation and the interlayer (Zhou et al., 2008), the interface property, the vertical compressive stress difference (Anderson, 1981; Biot et al., 1983; Chen and Chen, 1997; Warpinski et al., 1980, 1982), and the horizontal stress difference (Chen et al., 2008). The stress difference that prevents the propagation of hydraulic fractures is affected by the water pressure, fracture geometry, and fluid loss (Teufel and Clark, 1984; Warpinski et al., 1980).

Complex numerical model of hydraulic fracturing in natural reservoir has been developed by researchers (Weng, 2015). The stress intensity factor of hydraulic fractures propagating in layered strata has also been calculated (Zhao and Chen, 2010). Criteria for fracture penetrating through the bedding plane are developed (Biot et al., 1983), simple criteria to predict whether the hydraulic fracture will penetrate nonorthogonal natural fissures are derived (Gu et al., 2012), and criteria for hydraulic fracture to penetrate through natural fissures in the three-dimensional space are established (Cheng et al., 2014). Currently, there are two methods to calculate the size of two-dimensional fracture: the PKN model and the KGD model. The PKN model is based on the theory of vertical plane strain, while the KGD model is derived from horizontal plane strain (Chen et al., 2008). The fracture height is assumed to be constant in both PKN and KGD model. However, the fracture extends in the height direction to some extent. Therefore, four types of pseudo-three-dimensional models are developed. In these models, the hydraulic fracture is assumed to propagate with the shape of ellipse, and the variation of fracture height is taken into consideration. To simulate the three-dimensional fracture geometry, two methods are used: one is adopting the fracture propagation criterion (K = K_IC) to introduce the fracture height into the model, the other is describing the fracture extension in the vertical direction by the KGD model while describing the horizontal propagation by the PKN model. The common features of the pseudo-three-dimensional models are as follows: the three-dimensional flow in the fracture is simplified into two- or one-dimensional flow because the fracture length is larger than its height. The vector water flow parallel to the fracture surfaces is not considered and the water loss is ignored. The difference between elastic modulus and fracture toughness of different formations is not considered. The calculation model for hydraulic fracture morphology in layered strata is studied. Calculation methods for fracture strength factor are discussed (Chen et al., 1994). From the perspective of fracture toughness, the vertical extension height of hydraulic fractures is calculated and then the mechanism of hydraulic fracturing is systematically studied. Based on the basic equations of elasticity, a three-dimensional model of hydraulic fractures in layered medium is derived and numerically solved. Then a fully three-dimensional fracture propagation model is proposed to simulate the hydraulic fracturing more closely. In the fully three-dimensional model, the vertical and horizontal extensions are separately simulated. In this model, the geometrical morphology, three-dimensional deformation, and two-dimensional water flow within hydraulic fractures are more closely reflected. Internationally, a simple and practical set of models are established to predict fully three-dimensional hydraulic fracture (Bouteca, 1984), and a coupling model of three-dimensional hydraulic fracture propagation and two-dimensional proppant transportation is developed (Lee and Jantz, 1989). In China, a fully three-dimensional model of hydraulic fractures under nonuniform conditions is established, and the calculation method of fracture width is given out (Chen and Chen, 1997; Chen et al., 1995). The present study suggests that bedding planes have great influence on the propagation of artificial hydraulic fractures. However, the study of propagation law of hydraulic fractures under three-dimensional conditions is not sufficient. Further studies should be conducted on this problem.

Previous experiments have been conducted with samples containing prefabricated bedding planes using the true triaxial hydraulic fracturing experiment system. On the basis of previous experiments, the basic form, three-dimensional conceptual model and three behavior types the “≠” shaped hydraulic fracture are revealed (Huang and Liu, 2017). In this work, based on the previous experiments, a fully three-dimensional propagation model for multihydraulic fractures is derived with the water loss, vertical water flow, and bedding plane effect taken into consideration.

Development of the mechanical model of hydraulic fracturing in horizontal strata borehole

Fully three-dimensional propagation model of hydraulic fractures

At present, good effects are usually achieved in exploitation of coalbed methane and oil by hydraulic fracturing with horizontal borehole (Figure 1(a)). When hydraulic fractures generated in horizontal borehole propagate to the bedding plane, it may propagate both within the bedding plane with the shape of ellipse and penetrate the bedding plane vertically, forming “≠” shaped fractures (Figure 1(b)). To describe the morphology of the main hydraulic fracture and the hydraulic fracture within the bedding plane (which is referred to as bedding plane fracture) more precisely, a three-dimensional coordinate system (Figure 1(b)) is established. The center of the unsealed section of the borehole at the bottom of the borehole is taken as the origin of the coordinate. In the vertical direction, the origin of the coordinate is located in the center between the two bedding planes. The x-axis of the coordinate is along the axial of the borehole, the y-axis is along the height of the main hydraulic fracture, and the z-axis is along the vertical direction. The model consists of three layers, i.e. the upper, middle, and bottom layer (Figure 2(a) and (b)), which represent the fracturing layer, roof layer, and floor layer.

Figure 1.

Fully three-dimensional model for hydraulic fractures morphology in hydraulic fracturing in horizontal borehole in rock strata. (a) Horizontal fracturing borehole, (b) three-dimensional propagation morphology of hydraulic fractures in strata (Huang and Liu, 2017), (c) subdivision of the hydraulic fracture by unit plane along the x-axis, and (d) parameters to describe the hydraulic fracture in cross section.

Figure 2.

Mechanical model of the growth of the main hydraulic fracture and the bedding plane fracture. (a) Geometry of the main hydraulic fracture before intersecting with the bedding plane, (b) geometry of hydraulic fractures after the main hydraulic fracture intersects with the bedding plane, (c) stress state of the zone A at the tip of the main hydraulic fracture, and (d) stress state of the zone B at the tip of the bedding plane fracture.

To study the aperture, propagation direction, and range of hydraulic fractures, subdivide the hydraulic fracture into n cross sections by a series of planes that are perpendicular to the x-axis (Figure 1). The parameters to describe the hydraulic fracture in each plane include the aperture of the hydraulic fracture, $W_z(x,z,t)$ ; the aperture of the bedding plane, $W_c(x,y,t)$ ; the height of main hydraulic fracture, $H(x,t)$ ; the length of the bedding plane, $L(x,t)$ (Figure 1(d)).

Propagation criteria for the main hydraulic fracture tip without bedding plane

Before the main hydraulic fracture intersects with the bedding plane (Figure 2(a)), the fracture tip is subject to principal stresses ${\sigma }_1$ , ${\sigma }_2$ , ${\sigma }_3$ ; water pressure $p_z$ ; and pore pressure $p_0$ (Figure 2(c)). The rock in front of the fracture tip is assumed to be homogeneous and isotropous. Microcracks zone, which is induced by the splitting effect of water pressure, exists around the tip of the hydraulic fracture. Under the combined action of the effective water pressure and the far-field stress, microcracks will grow and interconnect along the predominant direction, forming macrofractures (Huang and Liu, 2017). Therefore, the propagation law of the main hydraulic fracture tip is as follows: when the stress intensity factor of the main hydraulic fracture tip is larger than its fracture toughness, the main hydraulic fracture begins to propagate.

Propagation criteria for hydraulic fracture tip with bedding plane

At the intersection of the hydraulic fracture with the bedding plane, fractures exist both in the direction of the main hydraulic fracture and in the bedding plane (Figure 2(b)). The stress state of the main hydraulic fracture tip is almost the same before and after intersection with the bedding plane. The tip of the bedding plane fracture is subject to the principal stress ${\sigma }_1$ , the water pressure $p_c$ , the pore pressure $p_c$ , and the tensile strength ${\sigma }_T$ (Figure 2(d)). The propagation criteria for the bedding plane fracture are as follows: when the resultant force of the effective water pressure and the normal stress at the fracture tip are larger than the tensile strength of the bedding plane, the fracture along the bedding plane grows forward.

Criteria for hydraulic fracture penetrating through the bedding plane

True triaxial experiments of hydraulic fracturing have proved that three types of propagation behaviors and shapes of hydraulic fracture may occur after the main hydraulic fracture intersects with the bedding plane: propagating along the bedding plane (with hydraulic fracture of “工”), penetrating through the bedding plane (with hydraulic fracture shape of “/”), or both the above two cases occur at the same time (with hydraulic fracture of “≠”) (Huang and Liu, 2017) (Figure 3).

Figure 3.

Propagation behavior models of hydraulic fractures after intersecting with the bedding plane. (a) Model I: hydraulic fracture propagates along the bedding plane, (b) Model II: hydraulic fractures propagate along the bedding plane and penetrate through the bedding plane at the same time, and (c) Model III: hydraulic fracture penetrates through the bedding plane.

1. Criteria for hydraulic fracture penetrating the bedding plane

The condition that the hydraulic fracture penetrates the bedding plane is

$$K_s>K_1$$

(1)

The propagation of the main hydraulic fracture belongs to the Model I fracture type

$$K_s\text{=}\left(p-\Delta p_z-p_{{\hspace{0.17em}}_{0z}}-{\sigma }_{nz}\right)\sqrt{\pi a}$$

(2)

Substituting equation (2) into equation (1) leads to

$$\left(p-\Delta p_z-p_{0z}-{\sigma }_3\right)\sqrt{\frac{\pi H}{2}}>K_{\text{I}}$$

(3)

Equation (3) can also be written as

where $K_s$ is the effective stress intensity factor at the fracture tip caused by the effective water pressure ( $\text{MPa}\cdot \sqrt{\text{m}}$ ), $K{}_1$ is the fracture toughness of the fracturing formation ( $\text{MPa}\cdot \sqrt{\text{m}}$ ), ${\sigma }_{nz}$ is the pressure difference between inside and outside the hydraulic fracture (MPa), p is the water pressure at the breakthrough point of the main hydraulic fracture and the bedding plane (MPa), $p_{0z}$ is the rock pore water pressure (MPa), $\Delta p_z$ is the pressure loss of water flowing within the main hydraulic fracture (MPa) and H is the total fracture height (m).

2. Criteria for hydraulic fracture propagating along the bedding plane

When the resultant force of the effective water pressure and the normal stress at the tip of the bedding plane fracture is larger than the tensile strength of the bedding plane, the hydraulic fracture propagates along the bedding plane (Figure 3(a))

$$p-\Delta p_c-p_{\text{0}c}-{\sigma }_{nc}>{\sigma }_t$$

where $p_{0c}$ is the original pore pressure (MPa), $\Delta p_{\text{c}}$ is the water pressure loss during the water flow from the fracturing borehole to the bedding plane (MPa), ${\sigma }_{nc}$ is the normal stress in the bedding plane (MPa), and ${\sigma }_t$ is the tensile strength of the bedding plane (MPa).

The stress state of the bedding plane can be written as (Figure 4)

$$\{\begin{array}{c}{\sigma }_{n1}=\frac{{\sigma }_{11}+{\sigma }_{31}}{2}-\frac{{\sigma }_{11}\text{-}{\sigma }_{31}}{2}\text{cos}2\theta \\ {\sigma }_{n2}=\frac{{\sigma }_{12}+{\sigma }_{32}}{2}-\frac{{\sigma }_{12}\text{-}{\sigma }_{32}}{2}\text{cos}2\theta \\ S_1=\frac{{\sigma }_{11}\text{-}{\sigma }_{31}}{2}\text{sin}2\theta \\ S_2=\frac{{\sigma }_{12}\text{-}{\sigma }_{32}}{2}\text{sin}2\theta \end{array}$$

(6)

Figure 4.

Stress state of the bedding plane during hydraulic fracturing.

In the direction perpendicular to the bedding plane, the forces reach a balance. Thus, the stress state of the area near the bedding plane can be written as

$$\{\begin{array}{c}{\sigma }_{nc}={\sigma }_{n1}=\frac{{\sigma }_{11}+{\sigma }_{31}}{2}-\frac{{\sigma }_{11}\text{-}{\sigma }_{31}}{2}\text{cos}2\theta ={\sigma }_{n2}=\frac{{\sigma }_{12}+{\sigma }_{32}}{2}-\frac{{\sigma }_{12}\text{-}{\sigma }_{32}}{2}\text{cos}2\theta \\ S_1 =\frac{{\sigma }_{11}\text{-}{\sigma }_{31}}{2}\text{sin}2\theta \\ S_2 =\frac{{\sigma }_{12}\text{-}{\sigma }_{32}}{2}\text{sin}2\theta \end{array}$$

(7)where ${\sigma }_{nc}$ is the normal stress perpendicular to the bedding plane (MPa), ${\sigma }_{11}$ is the maximum principal stress in the roof layer (MPa), ${\sigma }_{12}$ is the maximum principal stress in the fracturing layer, ${\sigma }_{31}$ is the minimum principal stress in the roof layer (MPa), ${\sigma }_{32}$ is the minimum principal stress in the fracturing layer (MPa), $S_1$ is the horizontal stress in the fracturing layer (MPa), $S_2$ is the horizontal stress in the roof layer (MPa), and $\theta$ is the angle between the bedding plane and the minimum principal stress ${\sigma }_1$ (°).

The propagation criteria for hydraulic fracture intersecting with the bedding plane can be written as follows.

When

$$\Delta p_c+p_{0c}+{\sigma }_{nc}+{\sigma }_t>\frac{\sqrt{2}{K}_1}{\sqrt{\pi H}}+{\sigma }_3+\Delta p_z+p_{0z}$$

(8)

The hydraulic fracture preferentially penetrates the bedding plane.

When

$$\Delta p_c+p_{0c}+{\sigma }_{nc}+{\sigma }_t<\frac{\sqrt{2}{K}_1}{\sqrt{\pi H}}+{\sigma }_3+\Delta p_z+p_{0z}$$

(9)

The hydraulic fracture preferentially propagates along the bedding plane

When

$$\Delta p_c+p_{0c}+{\sigma }_{nc}+{\sigma }_t=\frac{\sqrt{2}{K}_1}{\sqrt{\pi H}}+{\sigma }_3+\Delta p_z+p_{0z}$$

(10)

The hydraulic fracture both penetrate the bedding plane and propagate along the bedding plane at the same time.

Fully three-dimensional propagation model for hydraulic fractures

The height of the hydraulic fracture is not equal to the thickness of the fracturing layer. When hydraulic fractures propagate to a certain range, growing fractures exist in the roof and floor layer (Figure 1(b)). The main hydraulic fracture will propagate along the direction perpendicular to the minimum principal stress ${\sigma }_3$ , i.e. the direction of the fracture length. Because the minimum principal stress is usually along the horizontal direction, the main hydraulic fracture is therefore assumed to be vertical. Figure 4(a) shows that as the main hydraulic fracture propagates along the axis of the fracturing borehole (the x-axis), the fracture height grows too in the vertical direction (the z-axis). Take any unit cross section of the main hydraulic fracture, the fracture height is represented by $H\left(x,t\right)$ . With the propagation of the bedding plane fracture along the axis of the borehole (the x-axis), the fracture width grows along the y-axis within the bedding plane. Take any unit cross section of the bedding plane fracture, the horizontal width of the hydraulic is presented by $B\left(x,t\right)$ . The above analysis reveals that, in this model, the hydraulic fracture propagates in all of the length, width, and height direction. Correspondingly, water flows in the above three directions exist. Due to low aperture of the hydraulic fracture, the water flow in the fracture height direction is negligible. Thus, the water flow within the fracture can be assumed to be two-dimensional and parallel to the fracture surfaces. Moreover, the water loss is taken into consideration in this model to more closely simulate hydraulic fracturing in rock mass containing bedding plane.

Based on previous true triaxial physical simulation experiments of hydraulic fracturing, a fully three-dimensional propagation model for hydraulic fractures is developed (Figure 5). In this model, principles of elastic mechanics, fluid mechanics, and fracture mechanics are used and the vertical water flow, water leak-off, and interface effect are taken into consideration. Basic equations of this model include the continuity equation, equations of water pressure decline, equations of the length and aperture of the main hydraulic fracture and the bedding plane fracture, and equations of the height of the main hydraulic fracture and the width of the bedding plane fracture.

Figure 5.

Schematic diagram of the fully propagation model for hydraulic fractures in rock strata containing bedding planes.

Assumptions in the model

The rock and the bedding plane are homogeneous and linear elastic. The calculation of the physical dimension of the hydraulic fracture is based on the linear elastic theory.

The shape of the vertical cross section of the hydraulic fracture is ellipse.

The water flow within the fracture is laminar flow.

The impact of thermal field is not taken into consideration.

The water injection rate is constant.

Equation of continuity

The total injected water consists of the water to generate the main hydraulic fracture and the bedding plane fracture and the infiltrated water in the main hydraulic fracture and the bedding plane fracture. According to the principle of mass conservation, the total injected water equals the sum of the water generating the main hydraulic fracture and the bedding plane fracture, i.e.

$$\begin{array}{c}Q\cdot t\text{=} q_z\cdot t+q_c\cdot t\\ V_{f_z}\left(t\right)+V_{f_c}\left(t\right)+V_{L_z}\left(t\right)+V_{L_c}\left(t\right)\end{array}$$

(11)where Q is the water injection rate (m³/min), t is the injection time (min), $q_z\left(x,t\right)$ is the water flow at length of x in the main hydraulic fracture at t (m³/min), $q_c\left(x,t\right)$ is the water flow at length of x in the bedding plane fracture at t (m³/min), $V_{f_c}\left(t\right)$ is the total volume of the main hydraulic fracture at t (m³), $V_{f_c}\left(t\right)$ is the total volume of the bedding plane fracture at t (m³), $V_{L_z}(t)$ is the leak-off water volume of the main hydraulic fracture at t (m³), and $V_{L_c}(t)$ is the leak-off water volume of the bedding plane fracture at t (m³).

Along the length of the main hydraulic fracture and the bedding plane fracture, the increment of water flow, $q\left(x,t\right)$ , is the sum of the leak-off water volume, $\lambda \left(x,t\right)$ , and the increment of fracture volume with the propagation of hydraulic fractures, $dA\left(x,t\right)/dt$ , i.e.

$$\{\begin{array}{c}-\frac{d{q}_z\left(x,t\right)}{dx}={\lambda }_z\left(x,t\right)+\frac{d{A}_z\left(x,t\right)}{dt}\\ -\frac{d{q}_c\left(x,t\right)}{dx}={\lambda }_c\left(x,t\right)+\frac{d{A}_c\left(x,t\right)}{dt}\end{array}$$

(12)in equation (12)

$$\{\begin{array}{c}{\lambda }_z\left(x,t\right)\text{=}2H\left(x\right)\frac{C_z}{\sqrt{t-{\tau }_z\left(x\right)}}\\ {\lambda }_c\left(x,t\right)\text{=}2B\left(x\right)\frac{C_c}{\sqrt{t-{\tau }_c\left(x\right)}}\end{array}$$

(13)where t is the current water injection time (min), $A_z\left(x,t\right)$ is the water loss velocity at length x in the main hydraulic fracture at t (m/min), $A_c\left(x,t\right)$ is the water loss velocity at length x in the bedding plane fracture at t (m/min), ${\lambda }_z\left(x,t\right)$ is the cross sectional area of the main hydraulic fracture at length x and time t (m²), ${\lambda }_c\left(x,t\right)$ is the cross sectional area of the bedding plane fracture at length x and time t (m²), $C_z$ is the overall filtration coefficient of the main hydraulic fracture ( $\text{m/}\sqrt{\text{min}}$ ), $C_c$ is the overall filtration coefficient of the bedding plane fracture ( $\text{m}\text{/}\sqrt{\text{min}}$ ), ${\tau }_z\left(x\right)$ is the time the fracturing water flows from the borehole to length x of the main hydraulic fracture (min), and ${\tau }_c\left(x\right)$ is the time the fracturing water flows from the borehole to length x of the bedding plane fracture (min).

On the basis of previous research, the distribution of the flow rate along the fracture length is

$$\{\begin{array}{c}{q}_z\left(x\right)=\frac{Q_z}{2}\left[1-\frac{2}{\pi }\sin \left(\frac{x}{L_z}\right)\right]\\ q_c\left(x\right)=\frac{Q-Q_z}{2}\left[1-\frac{2}{\pi }\sin \left(\frac{x}{L_c}\right)\right]\end{array}$$

(14)where $L_z$ is half length of the main hydraulic fracture (m), $L_c$ is the half length of the bedding plane fracture (m), and $Q_z$ is the water volume for the propagation and leak-off of the main hydraulic fracture.

Pressure decline equations

The pressure decline equation of water flow in the length direction of the main hydraulic fracture (Wang and Zhang, 1998) is

$$\frac{dp}{dx}=-\frac{64\mu q_z\left(x\right)}{\pi H\left(x\right)W_{0z}{\hspace{0.17em}}^3}$$

(15)

Similarly, the pressure decline equation of water flow in the length direction of the bedding plane fracture is

$$\frac{dp}{dx}=-\frac{64\mu q_c\left(x\right)}{\pi B\left(x\right)W_{0c}{\hspace{0.17em}}^3}$$

(16)

Therefore

$$\frac{dp\left(x\right)}{dx}=F_2\left[W_0\left(x\right),H\left(x\right),q\left(x\right)\right]$$

(17)

For vertical hydraulic fracture, the maximum aperture of the ellipse-shaped cross section (Figure 5) (Wang and Zhang, 1998) is

$$W_z\left(x,t\right)=\frac{(1-v_z)}{G_z}H\Delta P_{zl}$$

(18)

Meanwhile, the following equation holds

$$G_z=\frac{E_z}{2\left(1+v_z\right)}$$

(19)

Substituting equation (19) into equation (18) and then substituting the result into equation (15), the water pressure decline equation can be written as (Wang and Zhang, 1998)

$$\frac{{\partial }^4\Delta P_{zl}}{\partial x}=-\frac{256}{\pi }\frac{\mu q_z{G}_z{\hspace{0.17em}}^3}{{\left(1-v_z\right)}^3{H}^4}$$

(20)

Calculate the integral of equation (20) from x=0 to x=L, then the pressure decline equation of the vertical hydraulic fracture in its length direction is

$$\Delta P_{zl}={\left(\frac{256}{\pi }\frac{\mu q_z{L}_z{G}_z{\hspace{0.17em}}^3}{{\left(1-v_z\right)}^3{H}^4}\right)}^{\frac{1}{4}}$$

(21)

Similarly, the pressure decline equation of vertical hydraulic fracture in its height direction can be derived as

$$\Delta P_{zh}={\left(\frac{256}{\pi }\frac{\mu q_{z}H{G}_z{\hspace{0.17em}}^3}{{\left(1-v_z\right)}^3{L}_z{\hspace{0.17em}}^4}\right)}^{\frac{1}{4}}$$

(22)

The pressure decline equation of the bedding plane fracture in its length direction is

$$\Delta P_{cl}={\left(\frac{256}{\pi }\frac{\mu q_c{L}_c{G}_c{\hspace{0.17em}}^3}{{\left(1-v_c\right)}^3{B}^4}\right)}^{\frac{1}{4}}$$

(23)

The pressure decline equation of the bedding plane fracture in its width direction is

$$\Delta P_{ch}={\left(\frac{256}{\pi }\frac{\mu q_{c}B{G}_c{\hspace{0.17em}}^3}{{\left(1-v_c\right)}^3{L}_c{\hspace{0.17em}}^4}\right)}^{\frac{1}{4}}$$

(24)where $W_z$ is the aperture of the main hydraulic fracture (m), $\mu$ is the viscosity of liquid (MPa s), $E_z$ is the elastic modulus of rock (MPa), $G_z$ is the shear modulus of rock (MPa), $v_z$ is the Poisson’s ratio of rock, $E_c$ is the elastic modulus of the bedding plane (MPa), $G_c$ is the shear modulus of the bedding plane (MPa), and $v_c$ is the Poisson’s ratio of the bedding plane.

Equations of the length and aperture of the main hydraulic fracture and the bedding plane fracture

The aperture of the three-dimensional fracture is equivalent to the fracture aperture in the vertical cross section of the fracture (Figure 1(c)), which is calculated by two-dimensional methods. The fracture aperture $W\left(x\right)$ is a function of the static pressure p and the fracture height $H\left(x\right)$ . To calculate the length and aperture of hydraulic fractures at the same time, the PKN model is modified by taking the water leak-off into consideration. Then the length of the main hydraulic fracture, $L_z$ , can be written as

$$L_z=\frac{q_z\sqrt{t}}{\pi H{C}_z}$$

(25)

And the width of the bedding plane fracture, $L_c$ , can be written as

$$L_c=\frac{q_c\sqrt{t}}{\pi B{C}_c}$$

(26)

The width of the main hydraulic fracture at its initiation location is

$$W_z\left(0,t\right)=4{\left[\frac{2\left(1-{\nu }_z\right)\mu q_z{\hspace{0.17em}}^2}{{\pi }^3{G}_z{C}_{z}H}\right]}^{\frac{1}{4}}{t}^{\frac{1}{8}}$$

(27)

The equation of the width of the main hydraulic fracture at any point can be derived as

$$W_z\left(x,t\right)=W_z\left(0,t\right){\left\{\frac{x}{L_z}{\sin }^{-1}\frac{x}{L_z}+{\left[1-{\left(\frac{x}{L_z}\right)}^2\right]}^{\frac{1}{2}}-\frac{\pi }{2}\frac{x}{L_z}\right\}}^{\frac{1}{4}}$$

(28)

The equation of the width of the bedding plane fracture at the initiation location can be written as

$$W_c\left(0,t\right)=4{\left[\frac{2\left(1-{\nu }_c\right)\mu q_c{\hspace{0.17em}}^2}{{\pi }^3{G}_c{C}_{c}B}\right]}^{\frac{1}{4}}{t}^{\frac{1}{8}}$$

(29)

The equation of the width at any point of the bedding plane fracture is

$$W_c\left(x,t\right)=W_c\left(0,t\right){\left\{\frac{x}{L_c}{\sin }^{-1}\frac{x}{L_c}+{\left[1-{\left(\frac{x}{L_c}\right)}^2\right]}^{\frac{1}{2}}-\frac{\pi }{2}\frac{x}{L_c}\right\}}^{\frac{1}{4}}$$

(30)where $W_c$ is the maximum aperture of the bedding plane fracture.

Equations of the height of the main hydraulic fracture and the width of the bedding plane fracture

According to the theory of fracture mechanics, the fracture propagates only when the stress intensity factor reaches the fracture toughness of the rock mass. Based on the calculation formula of stress intensity factor, the following equations can be obtained.

①. When hydraulic fractures propagate in the fracturing layer, i.e. $-\frac{H_0}{2}\le z\le \frac{H_0}{2}$ (Wang and Zhang, 1998)

$$W_{0z}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_z}{K}_1\sqrt{H\left(x\right)}$$

(31)

②. When hydraulic fractures extend beyond the fracturing layer, i.e. $\left|z\right|>\frac{H_0}{2}$ (Wang and Zhang, 1998)

$$W_{0z}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_z}{K}_2\sqrt{H\left(x\right)} +\frac{4\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\pi E_z}\left(S_2-S_1\right)H_0\left[\ln H\left(x\right)+\sqrt{H^2\left(x\right)-H_0{\hspace{0.17em}}^2}-\ln H_0\right]$$

(32)

Namely

$$W_0\left(x\right)\text{=}{F}_1^{\prime }H\left(x\right)$$

(33)

③. When hydraulic fractures propagate into the bedding plane, i.e. $\left|z\right|\text{=}\frac{H_0}{2}$

$$W_{0c}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_c{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_c}{K}_c\sqrt{B\left(x\right)}$$

(34)where $K_2$ is the fracture toughness of the roof and floor layer ( $\text{MPa}\cdot \sqrt{\text{m}}$ ) and $K_c$ is the fracture toughness of the bedding plane ( $\text{MPa}\cdot \sqrt{\text{m}}$ ).

Based on the above analysis, the fully three-dimensional propagation model for hydraulic fractures are as follows.

(1) Equations of continuity

The continuity equation of the main hydraulic fracture is

$$q(x\prime )=q(0)-\frac{V_2(x\prime )-V_1(x\prime )}{dt}-\frac{2H\left(x\right)Cn{a}_1}{\sqrt{t_1}}g(\tau \prime /t_1)\text{-}2S_{p}H\left(x\right)\frac{nx\prime }{t_1}$$

(35)

The continuity equation of the hydraulic fracture within the bedding plane

$$q(x\prime )=q(0)-\frac{V_2(x\prime )-V_1(x\prime )}{dt}-\frac{2\left(H_0+B\left(x\right)\right)+Cn{a}_1}{\sqrt{t_1}}g(\tau \prime /t_1)\text{-}2S_p\left(H_0+B\left(x\right)\right)\frac{nx\prime }{t_1}$$

(36)

(2) Equations of pressure decline

The pressure decline equation of the main hydraulic fracture is

$$\frac{dp}{dx}=-\frac{64\mu q_z\left(x\right)}{\pi H\left(x\right)W_{0z}{\hspace{0.17em}}^3}$$

(37)

The pressure decline equation of the hydraulic fracture within the bedding plane

$$\frac{dp}{dx}=-\frac{64\mu q_c\left(x\right)}{\pi B\left(x\right)W_{0c}{\hspace{0.17em}}^3}$$

(38)

(3) Equations of fracture length

The equation of the length of the main hydraulic fracture

$$W_{0z}\left(0,t\right)=4{\left[\frac{2\left(1-{\nu }_z\right)\mu q_z{\hspace{0.17em}}^2}{{\pi }^3{G}_z{C}_{z}H}\right]}^{\frac{1}{4}}{t}^{\frac{1}{8}}$$

(39)

The equation of the length of the hydraulic fracture within the bedding plane

$$W_{0c}\left(0,t\right)=4{\left[\frac{2\left(1-{\nu }_c\right)\mu q_c{\hspace{0.17em}}^2}{{\pi }^3{G}_c{C}_{c}B}\right]}^{\frac{1}{4}}{t}^{\frac{1}{8}}$$

(40)

(4) Equations of the height of the main hydraulic fracture is as follows

$$\{\begin{array}{cc}{W}_{0Z}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_z}{K}_1\sqrt{H\left(x\right)}& -\frac{H_0}{2}\le z\le \frac{H_0}{2}\\ W_{0Z}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_z}{K}_2\sqrt{H\left(x\right)}+\frac{4\left(1-{\nu }_z{\hspace{0.17em}}^2\right)}{\pi E_z}\left(S_2-S_1\right)H_0\left[\ln H\left(x\right)+\sqrt{H^2\left(x\right)-H_0{\hspace{0.17em}}^2}-\ln H_0\right]& |z|>\frac{H_0}{2}\end{array}$$

(41)

(5) Equation of the width of the hydraulic fracture within the bedding plane

$$W_{0c}\left(x\right)=\frac{2\sqrt{2}\left(1-{\nu }_c{\hspace{0.17em}}^2\right)}{\sqrt{\pi }{E}_c}{K}_c\sqrt{B\left(x\right)}z\text{=}\pm \frac{H_0}{2}$$

(42)

Verifications and analysis of the model

The fully three-dimensional propagation model is verified by physical experiments with true triaxial hydraulic fracturing experimental system (Huang, 2009; Huang et al., 2011). Samples with maximum size of 500 mm × 500 mm × 500 mm can be tested in this experimental system. Table 1 shows the experimental program of hydraulic fracturing experiments. The water injection rate is set as 500 ml/min.

Table 1.

Experimental program.

Sample order	Stress state (MPa)	Difference of principal stresses (MPa)σ1−σ3	Drying time of the bedding plane (h)	Effective water injection time (min)
21-2	σ1 = 6σ2 = 5.5σ3 = 5	1	21	9
42-1	σ1= 6σ2 = 5.5σ3 = 3	3	42	5

Cement mortar samples are used in this experiment to simulate rock. The main parameters in the experiments are as follows: the dynamic viscosity of water, $\mu$ , is 1.683 × 10⁻¹¹ N·s/m² (20°C); the total water injection rate is 0.0005 m³/min; the flow rate of the main hydraulic fracture, $q_z$ , is 0.0003 m³/min; the flow rate in the fracture within the bedding plane, $q_c$ , is 0.0002 m³/min; the total fluid loss coefficient of the main hydraulic fracture, $C_z$ , is 0.01 $\text{m}/\sqrt{\min }$ ; the total fluid loss coefficient of the bedding plane fracture, $C_c$ , is 0.013 $\text{m}/\sqrt{\min }$ ; the Poisson’s ratio of cement mortar, ${\nu }_z$ , is 0.18; the elastic modulus, $E_c$ , is 7.5 GPa; the shear modulus, $G_z$ , is 3.177 GPa; the Poisson’s ratio of the bedding plane, ${\nu }_c$ , is 0.4; the elastic modulus of the bedding plane, $E_z$ , is 0.59 GPa; the shear modulus of the bedding plane, $G_c$ , is 0.21 GPa; the fracture toughness of the cement mortar, $K_2$ , is 0.464 $\text{MPa}\cdot \sqrt{\text{m}}$ ; the fracture toughness of the bedding plane, $K_c$ , is 0.116 $\text{MPa}\cdot \sqrt{\text{m}}$ ; the main hydraulic fracture height is 0.225 m; the height of the bedding plane fracture is 0.111 m; the aperture of the main hydraulic fracture is 0.054 mm; the aperture of the fracture within the bedding plane is 0.090 mm (sample 42-1); and the horizontal principal stress difference, $S_2-S_1$ , is 1.5 MPa.

The length and height of the main hydraulic fracture

Substitute $q_z$ and $C_z$ into equation (25) to calculate the actual length of the main hydraulic fracture. In the sample 42-1, the effective water injection time is 9 min. When the fracture height is 225 mm, the calculated fracture length is 255 mm. Compared with the actual fracture length of 275 mm (Figure 6(a)), the calculated error is 7.3%. In the sample 21-2, the effective water injection time is 5 min when the fracture height is 185 mm, the calculated fracture length is 231 mm, the error of which is 15.6% compared with the actual value of 199 mm (Figure 6(b)).

Figure 6.

The geometry of the main hydraulic fracture (Huang and Liu, 2017). (a) Sample 21-2 and (b) sample 42-1.

The length and width of the bedding plane fracture

Substitute $q_c$ and $c_c$ into equation (26) to calculate the length of the bedding plane fracture. In the sample 42-1, the effective water injection time is 9 min. When the fracture width is 111 mm, the calculated fracture length is 265 mm with deviation of 1.5%, compared with the actual fracture length of 261 mm. In the sample 21-2, the effective water injection time is 5 min. When the fracture width reaches 108 mm, the calculated fracture length is 203 mm while the actual length is 192 mm (Figure 7(b)), the deviation of which is 5.7%.

Figure 7.

Geometry of the bedding plane fracture (Huang and Liu, 2017). (a) Sample 21-2 and (b) sample 42-1.

The aperture and the height of the main hydraulic fracture

The aperture of the main hydraulic fracture is calculated to be 0.054 mm by substituting the fracture height of 225 mm (Figure 6(a)) and other parameters into equation (27). Then substitute the aperture of the main hydraulic fracturing into equation (32) and obtain the theoretical fracture height of 180 mm. Compared with the actual fracture height, the deviation of the calculated fracture height is 20%.

Substitute the fracture height of 185 mm (Figure 6(b)) and other parameters into equation (27) and the aperture of the main hydraulic fracture is calculated to be 0.053 mm. Then substitute the aperture into the equation of fracture height (equation (32)), the theoretical fracture height is calculated to be 180 mm, the deviation of which is 2.7%.

The aperture and fracture width of the bedding plane fracture

Substitute the fracture width of 111 mm (Figure 7(a)) and other parameters into equation (29), the aperture of the main hydraulic fracture is calculated to be 0.09 mm. Then substitute the calculated fracture aperture into the equation of fracture height (equation 34) and obtain the theoretical fracture width of 117 mm. There exists an error of 5.1% between the calculated value and the actual value of the fracture width.

Substitute the fracture width of 108 mm (Figure 7(b)) and other parameters into equation (29), the aperture of the main hydraulic fracture is calculated to be 0.084 mm. Then substitute the calculated fracture aperture into equation (34) and the theoretical fracture width is calculated as 102 mm, with a deviation of 5.6%.

The deviations of the above calculated results are all less than 20%, which meets the requirements of engineering applications. Therefore, the mechanical model can be used in engineering.

Conclusions

Based on previous physical experiments, a mechanical model for hydraulic fracturing in horizontal borehole is established. A fully three-dimensional model for hydraulic fracture propagation morphology is developed.

On the basis of propagation behaviors of the main hydraulic fracture and the bedding plane fracture in physical experiments, the fracture failure mechanism of the main hydraulic fracture as well as the tensile failure mechanism of the bedding plane fracture is obtained, with the effect of pore pressure and the growth of microcracks taken into consideration. The mechanical conditions for the propagation of the main hydraulic and the bedding plane fracture are also derived. The criteria for hydraulic fractures penetrating through the bedding plane are given out.

On the basis of previous true triaxial physical simulation experiments of hydraulic fracturing, employing theories of elastic mechanics, fluid mechanics, fracture mechanics, and others, a fully three-dimensional propagation model of hydraulic fractures is established. In this mechanical model, the vertical water flow, water leak-off, and the effect of bedding plane are taken into consideration. The following basic equations of the model are solved: the equation of continuity, the equation of water pressure decline, the equation of the height and aperture of the main hydraulic fracture and the bedding plane fracture, and the equation of the height of the main hydraulic fracture and the width of the bedding plane fracture.

Hydraulic fracturing experiments in samples containing bedding planes are used to verify the basic equations of the length, aperture, width, and height of hydraulic fractures. All the deviations of the calculated values are less than 20%, which meet the requirements of the engineering application.

The model has its limits so some assumptions should be included: The rock and the bedding plane are homogeneous and linear elastic; the calculation of the physical dimension of the hydraulic fracture is based on the linear elastic theory; the shape of the vertical cross section of the hydraulic fracture is ellipse; the water flow within the fracture is laminar flow; the impact of thermal field is not taken into consideration; the water injection rate is constant.