Project Topics
Get Complete Project Material File(s) Now! »
Hydromechanical couplings
Now looking into the hydraulic fracturing from a phenomenological point of view we will see that hydraulic fracturing involves several coupled phenomena [2, 31, 139]:
the preferential uid ow in the fracture depends on the fracture faces.
mechanical deformation of the surrounding porous medium is induced by the uid pressure on fracture lips.
the deformations in the porous matrix induced by the tectonic shifts, the regional stress state or the residual stresses caused by the exploitation of the reservoir (subsidence).
the fracture and the porous medium exchange uid.
the fracture propagates, and therefore, hydraulic fracture is a moving boundary problem. Therefore, fully modelling the hydraulic fracturing process requires solving a coupled system of governing equations consisting of (1) elasticity equations that determine the relationship between the fracture opening and the uid pressure, (2) non-linear partial dierential equations for uid ow (usually lubrication theory) that relate the uid ow in the fracture to the fracture opening and the uid pressure gradient, (3) a fracture propagation criterion (usually given by assuming LEFM is valid) that allows for quasi-static fracture growth when the stress intensity factor is equal to the fracture toughness, and (4) diusion of fracturing uid in the rock formation (see Figure 1.4).
Propagation of uid-driven fractures has been the subject of numerous scientic papers [27, 2]. Most of these eorts have been driven by the oil and gas industry, and have been oriented towards the development of numerical models to predict the propagation of hydraulic fractures in the often complex geological settings under which extraction of hydrocarbons takes place.
In recent years, many researchers have recognized the existence of the complexity of hydrauli cfracture extension. Warpinski [186, 187] discovered that the main fracture and the branch fractures extended simultaneously through the eld tests and put forward the concept of fracture propagation zone. Afterwards, through physical simulation experiments [17, 18] and [36, 37] found that hydraulic fracture presented three kinds of extension path when it intersected with natural fractures: crossing natural fractures, extending along natural fractures, or the two cases occurring simultaneously. Although due to the insucient and imperfect understanding of the fracture network forming mechanism in shale reservoirs, there is always blindness in the fracturing design of shale reservoirs (Figure 1.5a).
Based on fracture extension characteristic in shale reservoirs, [188] classied hydraulic fractures into four major categories: the single plane biwing fracture, complex multiple fracture, complex multiple fracture with open natural fractures, and complex fracture network, as shown in Figure 1.5b.
Generalized Gri th-Irwin theory
In examining the stability of cracks, it is paramount to determine the relationships between the elastic stresses and the input energy rate or strain energy release rate in crack extension. The elastic energy release rate may be regarded as the force tending to open the crack [100]. Gri th’s theory shows that fracture occurs when the energy stored in the structure overcomes the surface energy of the material. He states that as the crack grows by a small amount da (a being the crack half length) at each end, energy ows into the crack tips where it is consumed in the fracture process, i.e., in overcoming the forces binding the material. Simply put, crack propagation will occur if the energy released is su cient to provide all the energy that is required for crack growth. If the crack grows stably, then dW = d (1.8).
where = 4a is the surface energy consumed in the creation of the cracks per unit thickness, and is the surface energy density, i.e., the energy required to create a unit traction-free crack surface. The quantity W is the work that can be extracted from the loading system during the crack formation. If the boundaries of the system are xed, then W equals the elastic strain energy released during the crack formation. If the boundaries are free to move as in the problem studied by Gri th [74], W equals to the change in the potential energy of the system.
For stable crack growth, the energy balance condition in Equation 1.8 may be rewritten as @ (W+=0 or @W = @ or G = R (1.9) @a @a @a where G = @W = energy release rate @a R = @ = crack resistance @a.
As a rst approximation, it can be assumed that the energy required to produce a crack is the same for each increment of crack growth, i.e., R is a constant. Thus, the above crack propagation criterion can be rewritten as G = Gc (1.10).
Cohesive Zone Model{CZM
LEFM has proven to be a useful tool for solving fracture problems, although it su ers from some limitations stemming from the strong assumptions that could be mainly listed as follows: If the size of the process zone approaches any relevant lengthscale of the medium (structure typical size, initial crack length, distance between crack tips), the response of the structure is subjected to a size e ect that LEFM fails to predict, In 3D, Gri th’s criterion stating that crack propagates when G = Gc is very hard to verify along a front, Crack may not initiate from a sound structure: a pre-crack has to be put in. Thus, the sensitivity of the response to the shape and orientation of the pre-crack has to be studied, and the failure load cannot always be accurately reproduced.
In order to adress some of these issues, cohesive zone model was introduced originally by Barenblatt [12] and Dugdale [57]. Later on, Hillerborg [83] used this model for concrete under the name of ctitious crack model. The cohesive process zone model is a general model which, in principle, is applicable to materials other than concrete or cementitious composites; e.g., crazing in polymers has been modelled using a cohesive surface methodology [176], and Schwalbe and collaborators modelled successfully the e ect of strength mismatch in welded joints using a cohesive zone model [116]. The cohesive zone model has also been applied successfully by [141, 63] to rocks to investigate hydraulic fracture propagation.
The cohesive zone model consists in introducing cohesive traction forces on the fracture walls, that obey a softening traction-displacement relation [141]. The cohesive traction is generally derived from a potential and is directly related to the displacement jump (see Figure 1.10). Once the critical stress c is reached, the damage process starts, irreversibly.
In this work the cohesive zone model relies on the principle of minimization of the total potential energy of the system: u; ; =JuK E p (u; ) min (1.41).
Finite Element Method{FEM
The FEM is the most widely employed numerical method for studying fracture mechanics problems (Figure 1.12). The formulation of the FEM is based on a variational statement of the governing physics. The domain is discretized into subdomains (elements) which are interconnected through common discrete points (nodes). The primary unknown eld variables are nodal values. The formulation reduces the problem to the solution of a system of algebraic equations in terms of the nodal variables (for dynamic problems, the result is a system of ordinary di erential equa-tions). Finite element systems tend to be relatively banded and symmetric for most problems [115]. FEM is capable of modelling complex geometries, loading conditions and heterogeneous material distributions [124]. [23] used FEM for modeling of fractures in orthotropic materials while [154] used FEM to model fractured induced anisotropy in poroelastic media. In particular, where they modeled fractures as very thin, highly permeable and compliant porous layers, allowing them to compute the wave velocities and quality factors at the macroscale as a function of frequency and propagation angle. Nevertheless, the classical displacement-based FEM is not able to describe the strain localisation properly since the di erential motion equations change type and lead to an ill-posed boundary value problem. A number of techniques have been implemented into the standard FEM to facilitate the computational simulation of crack propagation problems, such as element erosion methods and XFEMs [124, 168].
Finite Di erence Method{FDM
FDM is a continuum-based method similar to FEM that di ers in using a grid of nodes instead of elements for approximating the unknown elds such as displacement eld. However, the conven-tional FDM su ers from the use of regular grid system for the description of material heterogeneity, complex boundary conditions and fractures [59, 95]. To overcome these shortcomings, the general FDM has been improved particularly thanks to the development of nite volume methods, which make it capable of using irregular quadrilateral, triangular and Voronoi grids (Figure 1.15b). [3] used FDM to model fracture in bimaterial (Figure 4.1a). He divided the nodal points into inner points and ctitious boundary points, and the location of the crack tip is assumed to be at the cen-ter of the mesh and never on a mesh point. The commercial FLAC code is the most common FDM tool for stress analysis in geomechanics problems. [105] developed and implemented an algorithm based on LEFM approach in FLAC 2D code. According to the algorithm, each element comprises a microcrack with a random length that propagates when the critical value is satis ed by SIFs. [179] developed a constitutive model based on FDM to simulate dynamic fracturing in coal. Despite all these improvements, FDM still su ers from inability to model fracture propagation appropriately due to its continuum nature where the entire domain is employed for calculation [126]. There is no recorded evidence of extensive use of FDM in modeling crack propagation in anisotropic media nor hydrofracturing.
Boundary Element Method{BEM
The Boundary Element Method (BEM) has emerged as a powerful alternative to FEM par-ticularly in cases where better accuracy is required due to problems such as stress concentration. Also, if the medium extends to in nity, no arti cial boundaries such as those needed in FDM or FEM are required because BEM automatically satis es far- eld conditions. The most important feature of BEM, however, is that the solution is approximated at the boundaries, while equilib-rium and compatibility are exactly satis ed in the interior of the medium. In FDM and FEM, the approximations are made inside the medium [20]. The advantage of limiting the discretization to the boundaries is that the problem is reduced by one order: from three-dimensional (3D) to a 2D surface problem at the boundary and from 2D to a line problem. Hence, it only requires discretization of the surface rather than the volume (Figure 1.16).
This is in contrast to continuum methods (although we list BEM here amongst continuum methods), where the entire medium has to be discretized. The method is very attractive for those cases where the volume to boundary surface ratio is large. The technique used in BEM consists in essence of transforming the governing di erential equations, which apply to the entire medium, to integral equations that only consider boundary values [180, 26]. In a boundary value problem, some parameters such as stresses and displacements are known, while others are not, which then are part of the solution. There are two approaches to solve for the unknown parameters. In the rst approach (direct BEM), the unknowns are solved directly, and once they are found, stresses and displacements at any point in the continuum can be obtained directly from the solution. In the second approach (indirect BEM), the solution is given in terms of some ‘ ctitious’ quantities, typically stresses or displacements. The ctitious quantities are found rst, and the stresses and displacements at any point in the medium are expressed in terms of these ctitious quantities [20]. The di culties of the standard direct BEM in dealing with fracture problems such as the coincidence of crack nodes, gave rise to new techniques such as sub-region boundary element method (SBEM), displacement discontinuity methods (DDM), dual-boundary element method (DBEM) and dual-reciprocity boundary element method. The SBEM and DBEM are direct BEMs, while the DDM is an indirect BEM (Figure 1.17).
There has been extensive use of BEM in anisotropic fracture modeling [137, 101, 102, 169]. Pan et al. [137] proposed a new formulation of the BEM to calculate SIFs for cracked 2D anisotropic materials. They presented a new approach to collocate the displacement and traction integral Sourena MOOSAVI 38 After [126].
equations on the outside boundary of the problem (no-crack boundary) only and on one side of the crack surfaces only, respectively. This new method was considered to provide an alternative and yet e cient numerical technique for the study of cracked 2-D anisotropic media, and for the simulation of quasi-static crack propagation. Ke et al. [101] presented a systematic procedure for determining fracture toughness of an anisotropic marble using the diametral compression test (Brazilian test) with a central crack on the discs. Their new fracture criterion is based on the examination of mode I, mode II and mixed mode (I{II) fracture toughness for di erent crack angles and anisotropic orientation. The method has also been utilized in modeling hydrofracture phenomenon [109, 29, 67]. Legan et al. [109] modeled the fracturing process taking into account the inhomogeneity of the stress state near the hole in a cylindrical bodies using the boundary elements method (in the variant of the ctitious load method) and the gradient fracture criterion. Cao et al. [29] proposed an improved Boundary Element Method for modeling uid ow through fractured porous medium. In their proposed method they developed a theoretically sound, and practically robust numerical algorithm to accurately capture the ow behavior and dynamics in fractured reservoirs.
Molecular dynamics{MD
Because of the exponential growth of computing power, large-scale atomic simulations are be-ing developed rapidly to study the failure mechanisms of materials [197]. Molecular Dynamics is a time-dependent numerical solution of Newton’s equation of motion for all particles in atomic-scale [143]. The model in MD is composed of a collection of interacting spherical atoms under assumed interaction potential. The interaction are described using potential functions, i.e., Hooke’s law, Lennard{Jones potential, embedded atom method potential and the reactive force eld interatomic potential. Several studies have investigated the di erent aspects of crack initiation and propagation mechanism, such as the plastic deformation process at the crack, CZM parameters and dynamic crack processes using MD [198]. Generally, the MD simulation is a very useful tool for studying the change in the microstructure and therefore, it is a suitable technique for investigating crack nucle-ation and propagation at the micro-scale. However, the small computational system sizes and short time scales are two major limitations of this technique. Additionally, the nano/microstructures of rock materials are too complicated to model due to there being a multi-phase material. Molecular dynamics couldn’t nd its place among popular methods for fracture modeling in anisotropic rocks nor is it popular method among researchers for hydraulic fracturing modeling.
Crack propagation along non-prede ned paths in TI medium
In this section we discuss fracture propagation along non-prede ned paths in transverse isotropic media. The proposed model is an extension of the model proposed by [141] where we integrate the transverse isotropic characteristic of the medium. The crack initiation is controlled by the CZM used by [141], where the turning (bifurcation) angle is determined a posteriori based on the post-processing of the cohesive state.
The procedure for the propagation on non-prede ned paths has originally been suggested by [64] for crack propagation in concrete. It is based on the cohesive zone model depicted previously (Section 1.2). The originality of the procedure proposed by [64] lies in the a posteriori compu-tation of the crack advance based on the computed cohesive state, instead of a most common determination beforehand from the stress state ahead of the front.
In previous studies the discontinuities in the HM-XFEM model have been represented by a nor-mal level set eld only. Along this immutable discontinuity surface, the CZM allows to distinguish the adherent zone from the debonding zones [141]. For the crack propagation on non-prede ned paths, [64] suggested the crack front to be implicitly located with a tangential level-set eld.
Table of contents :
1 Fracture propagation in TI media
1.1 Phenomenology and observation
1.1.1 The inuence of transverse isotropy
1.1.2 Hydromechanical couplings
1.2 Theoretical aspects
1.2.1 Generalized Grith-Irwin theory
1.2.2 Calculating SIF using J-integral
1.2.3 Bifurcation angle
1.2.4 Cohesive Zone Model{CZM
1.3 Numerical modeling
1.3.1 Continuum{based methods
1.3.2 Discontinuum{based methods
1.4 Conclusion
2 A fully coupled 3D HM-XFEM model for transverse isotropic media
2.1 Introduction
2.2 HM-XFEM model representation
2.2.1 Domain denition and hypotheses
2.2.2 Governing equations
2.2.3 Variational formulations of conservation equations
2.2.4 Discretization with XFEM
2.2.5 Crack propagation along non-predened paths in TI medium
2.2.6 Numerical implementation
2.3 Validation of the HM-XFEM model
2.3.1 Properties and geometry of the medium and boundary conditions
2.3.2 Numerical results
2.4 Conclusion
3 Modeling fracture propagation in TI media using the BPM
3.1 A BPM for TI materials
3.1.1 Formulation
3.1.2 Introduction of TI
3.2 Calibration of the model
3.2.1 Reference behaviors
3.3 Fracture propagation
3.3.1 Homogeneous/isotropic material
3.3.2 TI medium
3.4 Conclusion
4 Conclusion
References
