Fracture Propagation Behavior of Jointed Rocks in Hydraulic Fracturing

Jointed rocks are typical examples of heterogeneous materials with joints.*e existence of joints influences the physical properties of rock mass and propagation of fractures, which can affect production operations in engineering. A series of simulations is performed to understand the failure patterns and fracture propagation behavior of jointed rocks in hydraulic fracturing. *ree points, that is, dip-angle joint, joint strength, and complex joints, are considered in the simulations. Results demonstrate three basic kinds of hydraulic fractures on jointed rock, namely, along the joint, across the joint, and partly along the joint and partly across the joint.*emaximum principal stress is the control factor of fracture propagation in global scale, and the joint plane is the control factor of fracture propagation in local scale. In the propagation path, when the dip angle is small or the joint is weak, the fracture propagates along the joint; otherwise, the fracture propagates across the joint. In the multijoint interconnection models, hydraulic fractures propagate along joints in the tensile stress regions near the propagating fracture tip without dip angle limitation. Subsequently, the fractures connect with one another to form a complex fracture network based on the joint morphology.


Introduction
Joint is a structural plane of discontinuities formed by sedimentation.Jointed rock masses are commonly encountered in civil, mining, and petroleum engineering because they are universally distributed over the earth's surface.As a control factor of rock mass stability [1], joints influence rock failure patterns and fracture propagation behavior.Hydraulic fracturing is commonly used in petroleum engineering, mining, and geotechnical industries, and it is carried out as an important technique to further improve production efficiency.e first experimental hydraulic fracturing treatment in the United States was conducted in the Hugoton Gas Field in Grant County, Kansas [2], and this technology has been widely used worldwide since 1947.In recent years, hydraulic fracturing was also used for jointed rock masses of sandstones, mudstones, shales, and coal seams to form complex fracture networks and improve the permeability of unconventional reservoirs [3][4][5].Hydraulic fracture propagation on jointed rocks is complex due to the different morphologies and composition of joints.Extensive experimental and numerical studies on the actual kinematics of hydraulic fracture propagation of jointed rocks were performed by researchers all over the world [6][7][8][9][10][11][12].Discontinuities dominate the geometry, deformation modulus, strength, failure behavior, and permeability of rocks, and even the local magnitude and direction of in situ stress fields [1].A study on the mechanical influence of joints and the growth process of hydraulic fractures of jointed rocks is essential, and it is even helpful in determining efficient fracturing scenarios or production designs in field operations.
e main purpose of this study is to numerically investigate how joints affect the initiation and propagation of hydraulic fractures in different conditions.Numerical simulations have been proven advantageous in modeling.
us, the hydraulic fracturing process of parallel and complex joint models are simulated by using the Rock Failure Process Analysis (RFPA) 2D 2.0-Flow.e RFPA is a code to simulate fracture propagation following the continuum damage principle [13,14] of heterogeneous rock materials.

Brief Introduction of the RFPA Code
Rocks are typical examples of heterogeneous natural geological materials [24].
ese solid media are often described as homogeneous and continuous materials, and their failure processes are explained by elasticity-plasticity principles in classical mechanics.However, recent studies on the mechanical and unique characteristics of rock materials have shown the limitations of the classical method.Subsequently, rock failure is assumed to be a nonlinear process from mesoscopic weakening to macroscopic damaging [25][26][27].

Nonlinear Macroscopic Model Described by Linear Mesoscopic Elements.
e RFPA code is based on finite element and statistical damage theory, and it can simulate the failure process of heterogeneous and permeable geomaterials.e heterogeneity of materials and the random distribution of defects of these materials are also considered by the RFPA code.e mesoscopic elements in the numerical model are assumed to be isotropic and homogeneous, and their mechanical properties (e.g., Young's modulus, Poisson's ratio, and strength properties, among others) are also assumed to be linear.e mesoscopic elements are statistically distributed (e.g., normal, Poisson, and Weibull distributions) to describe the mechanical properties of the nonlinear macromodel (Figure 1).
In RFPA, four-node isoparametric elements are used as the basic elements (Figure 2).In the simulation, the mesoscopic elements are defined as damaged if they satisfy the strength criterion.As stress increases, numerous elements fail.Moreover, when these elements connect with one another, the heterogeneous materials reach a failure process.
e statistical distribution density f(α) of the mechanical parameter of the mesoscopic elements is defined by (1).e heterogeneity index m [12,28] in the RFPA code is used to describe the heterogeneity of solid materials.A high m value indicates the presence of highly homogeneous materials, whereas a small m value denotes the existence of highly inhomogeneous materials (Figure 3).In the equation [29], α represents the mechanical parameters (e.g., Young's modulus, Poisson's ratio, tensional strength, and compressive strength) of the elements and α 0 denotes the scale (1)

Primitive Phase Transition.
e mesoscopic elements of an RFPA-based model are called primitives.ree primitive types (matrix, air, and contact) are considered for describing the deformation and failure processes of solid materials (Figure 4).When the strain is located between the maximum tensile and compressive strains, the element is called a matrix primitive.A matrix primitive contains elastic and residual phases.e matrix primitive changes into air primitive when the deformation exceeds the maximum tensile strain, and it changes into contact primitive when the deformation exceeds the maximum compressive strain.Changes occur in the mechanical properties of elements during primitive transformation, and this change is called phase transition.

Main Governing Equations of the RFPA-Flow Code.
e geomaterial (rock) in the RFPA-Flow code is assumed to be fully saturated with fluid flow in accordance with Darcy's law.In addition, the coupled process between stress/strain and fluid flow in the deforming rock mass is governed by Biot's consolidation theory.e rock medium is an elasticbrittle material with residual strength, and the acoustic emission (AE) of rock also proves it.us, the mechanical behavior of rock material in the loading and unloading process conforms with elastic damage theory.e maximum tensile strength criterion and Mohr-Coulomb criterion are used as the damage judgment of the mesoscopic element in the RFPA code.e fundamental equations [30] of fluid-solid coupling in the RFPA code are as follows: equilibrium equation : where ρ is the volume density; σ ij is the total stress; ε ij is the total strain and ε v is the volume strain; U i is the displacement; p is pore pressure; λ is the Lame coefficient; δ ij is the Kronecher constant; G is the modulus of shear deformation; Q is Biot's constant; K is the permeability coefficient; and α is the coefficient of pore pressure which is determined by the experiment.
Under the same loading conditions, fluid flow behavior differs between the propagating fractures and the preexisting fractures.In the failure process of rock mass, the permeability of a rock changes all the time.us, the propagation of preexisting fractures, initiation of new fractures, coupling relationship of flow-stress-damage in these two kinds of fractures, and variations of rock permeability in the failure process are considered in the RFPA code.
In the elastic state, the relationship of the stress and permeability coefficient is described by where K 0 is the initial permeability coefficient; σ ′ is the effective stress; and b is the coupling parameter.
In describing the dependency of permeability on the stress and the damage, the following coupling equation [31] is introduced: Advances in Materials Science and Engineering where α is the coefficient of pore pressure; ξ(ξ ≥ 1.0) is a mutation coefficient of permeability accounting for the increase in permeability once the element reaches the damage state; and β is the coupling parameter reflecting the influence of stress on the permeability coefficient.e value of the coefficients ξ, α, and β are determined experimentally, and they vary with the stress state and damage evolution of the rock.In the RFPA code, 0 ≤ α < 1 and ξ � 1 are assumed for the mesoscopic element in the elastic stage (as indicated by D � 0 in ( 7)).Once damage occurs, the permeability of the mesoscopic element increases significantly (i.e., see ( 8) and ( 11)).Furthermore, ξ � 5 is assumed for the damaged element (denoted by 0 < D < 1), and α � 1, ξ � 100 is assumed for the fully damaged element (denoted by D � 1).As illustrated in Figure 5, the elastic-brittle damage constitutive relation of an element under uniaxial tensile and compressive stress is used to govern the failure of this mesoscopic element.When the stress in an element satisfies the strength criteria, the element begins to accumulate damage and the elastic modulus of this mesoscopic element degrades gradually as damage progresses.e elastic modulus of the damaged element is defined as where D represents the damage variable, and E and E 0 are the elastic modulus of the damaged and undamaged elements.e element and its damage are assumed to be isotropic and elastic.erefore, as the stress increases, an increasing number of elements fail.Moreover, as the failed elements connect with one other, the whole failure process of the heterogeneous materials is reached.
When the tensile stress in the mesoscopic element reaches its tensile strength f t , and the damage variable D of the element under uniaxial tension can be expressed as where f tr is the residual tensile strength; ε t0 is the maximum tensile strain at the elastic limit; and ε tu is the ultimate tensile strain of the element at which this element is damaged completely (Figure 5).An ultimate strain coefficient η is employed to describe the relationship of ε t0 and ε tu , such that ε tu � ηε t0 .
As the damage variable D of the element changes in the damage process, the permeability coefficient of the mesoscopic element is expressed as Similarly, the damage variable D and permeability coefficient K of the mesoscopic element in the compressive state (Figure 5) in the damage process of geomaterials can be obtained on the basis of the Mohr-Coulomb strength criterion.
When the shear stress in the element satisfies the Mohr-Coulomb failure criterion, and the damage variable D of the element under uniaxial compression can be expressed as where ϕ is the friction angle; f c is the compression strength; f cr is the residual compression strength; and ε c0 is the maximum compression strain.
As the damage variable D of the element changes in the damage process, the permeability coefficient of the mesoscopic element is expressed as

Analysis of Hydraulic Fracture Propagation after Fracture Connection to Joint
Hydraulic fracture propagation, which changes path as it approaches a natural structural plane, is controlled by in situ stress, joint strength, and dip angle between fracture and joint, among others [32].With the existence of a structural As the hydraulic fracture approaches the joint and connects with the initial stage of the fracture-joint intersection, the fluid is injected and flows into the joint.e pressure in the joint is where p net is the net pressure.ree types of hydraulic fracture propagations occur as fluid flows into the joint (Figure 6(b)).

Across the Joint.
When the hydraulic fracture propagates across the joint and intersects with the joint, the pressure in the intersection point is where T is the tensile strength of the rock matrix.Equations ( 12) and ( 13) are integrated into (14), and the net pressure in the fracture is obtained as 3.2.Along the Joint.In this situation, the hydraulic fracture connects with the joint and creates a fracture-joint overlap.en, as the hydraulic fracture propagates along the joint, the fluid is injected.e pressure in the joint is where T b is the tensile strength of the joint.Equations ( 12) and ( 13) are integrated into (16), and the net pressure in the fracture is obtained as After the fracture propagates along the joint and overlaps the joint plane, it propagates along these overlapped joints.
en, the fracture continues to propagate in the maximum principal stress direction in the joint tip.

Partly along the Joint and Partly across the Joint.
In the propagation path along the joint plane (Figure 6(b)), the hydraulic fracture across the joint is in a weak point.e pressure in this weak point is where T w is the tensile strength of the weak point.Equations ( 12) and ( 13) are integrated into (18), and the net pressure in the fracture is obtained as

Numerical Modeling
e simulation of hydraulic fracture propagation on jointed rocks is divided into three groups.e first group is the simulation of dip angle models, as shown in Figure 7(a).
e dip angle α is the angle between the joint plane and the maximum principal stress direction.e α values of the seven models are 0 °, 15 °, 30 °, 45 °, 60 °, 75 °, and 90 °(f c � 20 MPa).e second group is the simulation of different joint strength models.In this group, the dip angle of the joint is 60 °, and the joint strengths of the five models are 10, 20, 30, 40, and 80 MPa.As shown in Figures 7(b) and 7(c), the complex joint in the third group consists of two multijoint models and two masonry-shaped joint models after the maximum principal stress directions are changed.In the models, a wellbore with a radius of 32 mm is assumed to be located at the center of the 640 mm × 640 mm specimen discretized into 102,400 elements.e process of defining joints: first draw lines on the existing rock material model and then assign the joint material properties to the mesoscopic elements on the lines.e thickness of the joint is the size of one element, 2 mm.In the parallel joint models, the distance between the adjacent planes is 40 mm.In the multijoint models, the two dip angle joints and the two Advances in Materials Science and Engineering strength joints are positioned horizontally and vertically around the well.In the masonry-shaped joint models, the distances of the horizontally and vertically adjacent parallel joint planes are 40 and 80 mm, respectively.e confining pressures of σ H � 15 MPa and σ h � 10 MPa (in complex joint models: σ H � 15 MPa, σ h � 10 MPa; σ H � 10 MPa, σ h � 15 MPa) are imposed on the outward boundaries of the models.An increasing water pressure, with an initial value of 12 MPa followed by 0.1 MPa increments per step, is considered for the wellbore.e rock and joint parameters employed in the simulations are shown in Table 1.

Fracture Propagation of Parallel Joint Models with Different Dip Angles
As shown in Figure 8, it depicts the hydraulic fractures propagation of different dip angle specimens.Hydraulic fracture initiation and propagation are significantly changed when the dip angle is increased.Figure 8 shows that when the dip angle is set between 0 °and 15 °, the fracture initiated in the wellbore propagates in a straight path along the joint to form a well-connected pathway for fluid flow because of the small fracture toughness of the joint material.In the process of hydraulic fracturing, a tensile stress region (highlighted region) emerges around the propagating fracture tip, which causes the fracture to propagate along the weakened joint plane.In this situation, hydraulic fracture propagation is mainly controlled by the joint plane.
When the dip angle is slightly increased (30 °to 45 °), the main fracture and some secondary ones propagate along the parallel joint planes in the tensile stress region to form multiple parallel straight fractures.Within this dip-angle range, the enlarged dip angle develops more secondary fractures.Fracture propagation is partially controlled by the joint plane.
When the dip angle is increased to 60 °, the fractures propagate along the joint plane.Moreover, a secondary fracture is induced in the maximum principal stress direction and connects with other secondary parallel fractures.When fractures further propagate, a complex fracture network is formed.In this situation, fracture propagation is controlled by the joint plane and the maximum principal stress.
For large dip angles (75 °to 90 °), the hydraulic fracture ceases to propagate along the joint planes.In the propagation process, when the fracture meets discontinuities (joints) of strongly contrasting mechanical properties, it propagates across multiple joints and forms a jagged fracture with multiple branches in the maximum principal stress direction.In the numerical models, the mechanical properties of mesoscopic elements are presented as a Weibull distribution [33] owing to the heterogeneity of jointed rocks.
en, fracture initiation and propagation occur on the weak elements, that is, these weak elements fail under high stress when pressure increases, and they may connect with one another and with other numerous isolated failure elements to form macroscopic fractures [33,34].
e termination of fracturing at the frequently contacting joint interfaces limits horizontal flow and instead produces complex flow paths.As a result, some secondary tensile fractures propagate along the joint near the terminated main fracture tip.Furthermore, the main fracture length is shorter than those depicted by previous simulations.In this situation, the fracture propagation is mainly controlled by the maximum principal stress.
Figure 9 presents the linear increase in breakdown pressure as the dip angle is increased.e smaller the dip angle, the easier the specimen is fractured.However, a small dip angle does not necessarily lead to good production operations in engineering, that is, a single fracture unrealistically portrays the efficient economic development of unconventional reservoirs, and thus, multiple fractures are   Advances in Materials Science and Engineering needed for reservoir stimulation [35][36][37].e simulation of this study suggests that the dip angles of 30 °to 60 °for parallel jointed rocks can efficiently propagate a main fracture and induce secondary ones to form a complex fracture network.

Fracture Propagation in Joint Strength Models
To understand the hydraulic fracture propagation behavior of the joint strength specimens, a dip angle of 60 °is considered on the basis of the first group simulation, that is, the hydraulic fracture propagation of the 60 °dip model is simultaneously controlled by the joint plane and the maximum principal stress.
Figure 10 shows that as joint strength increases, the main fracture shifts from along-the-joint to across-the-joint propagation.e influence of joints on fracture evolution is weakened gradually.In this situation, the initiation and propagation of hydraulic fractures are controlled initially by the joint plane, followed by the maximum principal stress.
A joint strength of 10 MPa implies that the joint is much weaker than the rock material.In this situation, fracture propagation is mainly controlled by the weak joint plane.When the joint strength is set to 20 MPa and hydraulic pressure is increased, secondary fractures are induced along and across the joints, and they connect with one another upon intersection.In this situation, fracture evolution is controlled by the joint and the maximum principal stress.When the joint strength is increased to 30 MPa, the joint influence on fracture growth is significantly weakened, and the fracture propagates partly along the joint and partly across the joint to form a horizontally crooked fracture.When the joint strength is set to 40 MPa, the main fracture propagates horizontally and a jagged fracture with multiple branches is formed.A joint strength of 80 MPa is similar to that of a rock matrix.In this situation, the joint influence on fracture evolution is negligible, and the maximum principal stress mainly controls fracture propagation.When the fluid flows into the fracture, the fracture continues to propagate and later attains a relatively smooth interface.

Numerical Simulation of Complex
Joint Models e previous simulations in this study show that the joint plane and the maximum principal stress both control the propagation of hydraulic fractures in different situations.
e multijoint model with different dip angles and strength joints (Figure 7  e simulation of the multijoint models indicates that the maximum principal stress controls fracture propagation in global scale, whereas the joint plane controls fracture propagation in local scale.However, all these models consider only independent joints.Subsequently, masonryshaped joint models are numerically simulated to capture the fracture propagation processes of interconnecting jointed rocks. As shown in Figure 12(a), the hydraulic fracture propagates along the straight joint plane under the confining pressures of σ H � 15 MPa and σ h � 10 MPa. e breakdown pressure and the failure mode are similar to those of the dip angle model of 0 °(Figure 8).When the confining pressure is changed to σ H � 10 MPa and σ h � 15 MPa (Figure 12(b)), the main fracture propagates vertically, while secondary hydraulic fractures initiate and propagate along the horizontal and vertical joints under tensile stress conditions (Figures 12(b) and 12(d)).When the hydraulic pressure is further increased, the propagating fractures connect with one another and merge to form a masonry-shaped fracture network.

Comparison of Numerical and Experimental Results
To prove the reliability of the RFPA code in the numerical simulation of rocks during hydraulic fracturing, the numerical model of shale was also considered.e dimension of the shale model was 300 mm × 300 mm, and the radius of the well was 12 mm.e model was discretized into 90,000 elements.e confining pressure of σ H � 20 MPa and σ h � 17 MPa was imposed on the outward boundaries of the en, the simulation results were compared with experimental results [38] (Figure 13).Plane-strain calculation was used for the simulation.e shale and joint parameters employed in the simulation are shown in Table 2.
Figure 13 shows the main hydraulic fracture initiating and propagating perpendicularly to the joint under the control of the maximum principal stress σ H .Given the weak cementation effect and small fracture toughness of the joint material, some secondary fractures propagate along the joint planes when the main fracture propagates.
is finding suggests that hydraulic fracture propagation is controlled by the joint plane in local scale.However, because the intensity of fracture toughness is large in the perpendicular direction of the joint, hydraulic fracture propagation is prevented after it comes across multiple joints.en, a jagged main fracture with turning and intersection points is formed.In the hydraulic fracturing process of the shale, a complex fracture network is formed when the hydraulic fractures and structural planes connected with one another.
However, some differences are observed between the 2D numerical simulation and the 3D experiment.For instance, the number of fractured joints and the degree of fracturing are lower in the simulation compared with those in the experiment.Nonetheless, the initiation and propagation of Advances in Materials Science and Engineering fractures are essentially consistent between the simulated and experimental results.

Conclusions
In this study, the RFPA code is applied to simulate the hydraulic fracturing of heterogeneous jointed rocks.However, reality is often much more complex than simulations to which numerical models are applied.Nevertheless, the   As depicted by the multijoint interconnection models, the main fracture propagates in the maximum principal stress direction.When multiple joints interconnect, secondary hydraulic fractures propagate along the joints within tensile stress regions without dip angle limitation.When hydraulic pressure is increased, the fractures connect with one another to form a complex fracture network based on the joint morphology.e intensity of fracture toughness is large in the perpendicular direction of the joint, and this prevents the propagation of hydraulic fractures.When hydraulic fractures frequently appear across the joints, an apparent decrease in fracture length is observed.During the main fracture propagation process, some joints appear with secondary fractures because of the weak cementation effect and limited fracture toughness of the joint material.As hydraulic fractures connect with structural planes, a complex fracture network can be formed.

Data Availability
e data used to support the findings of this study are available from the corresponding author upon request.

Figure 1 :
Figure 1: Heterogeneous distribution of mechanical properties from linear mesoscopic mode to nonlinear macroscopic mode.

Figure 5 :
Figure 5: Elastic-brittle damage constitutive law of mesoscopic elements subjected to uniaxial compressive and tensile stress.

4
Advances in Materials Science and Engineering plane, hydraulic fractures either cease to propagate or change direction to find the path of least resistance and energy damage.As depicted by Figure 6(a), the normal and shear stresses of the joint under in situ stress on the basis of Mohr-Coulomb strength theory are

Figure 6 :
Figure 6: Reaction of joint on hydraulic fracture.
(b)), and the masonry-shaped joint model with horizontal and vertical joints (Figure 7(c)), are simulated to capture the hydraulic fracture propagation on jointed rocks.

Figure 11 (
a) presents the simulated hydraulic fracture initiation and propagation of the multijoint model in the horizontal direction under the confining pressures of σ H � 15 MPa and σ h � 10 MPa. e fracture deterministically selects the path of least resistance and energy damage through the rock, and it propagates to the proximal tip of the joint after the hydraulic fracture is initiated in the right side of the well.When fluid flows into the fracture, the fracture propagates along the small-dip (45 °) joint.en, the fracture turns horizontally under maximum principal stress condition in the distal tip of the joint.At the other side of the well, the fracture propagates across the large-dip (75°)joint and forms a jagged fracture with branches.When the confining pressure is changed to σ H � 10 MPa and σ h � 15 MPa (Figure11(b)), the hydraulic fracture initiates and propagates vertically and spreads along the weak joint (60 °, f c � 10 MPa) and across the high-strength joint (60 °, f c � 80 MPa).
present study provides several interesting insights into fracture initiation and propagation mechanisms associated with hydraulic fracturing.A number of conclusions are derived from the numerical simulations.Fracture propagation is controlled by the maximum principal stress in global scale and the joint plane in local scale.In the maximum principal stress direction, hydraulic fractures propagate along small-dip joints (α ≤ 60 °in this paper) or weak joints (f c ≤ 20 MPa in this paper) and across

Figure 13 :
Figure 13: Fracture propagation in the numerical and experimental [38] results of the shale model.(a) Numerical results.(b) Experimental results.

Table 1 :
Physicomechanical parameters of rock and joint employed in simulations.

Table 2 :
Physicomechanical parameters of shale and joint employed in simulation.
Advances in Materials Science and Engineering