Numerical Simulation of 3D Hydraulic Fracturing Based on an Improved Flow-Stress-Damage Model and a Parallel FEM Technique

The failure mechanism of hydraulic fractures in heterogeneous geological materials is an important topic in mining and petroleum engineering. A three-dimensional (3D) finite element model that considers the coupled effects of seepage, damage, and the stress field is introduced. This model is based on a previously developed two-dimensional (2D) version of the model (RFPA2D-Rock Failure Process Analysis). The RFPA3D-Parallel model is developed using a parallel finite element method with a message-passing interface library. The constitutive law of this model considers strength and stiffness degradation, stress-dependent permeability for the pre-peak stage, and deformation-dependent permeability for the post-peak stage. Using this model, 3D modelling of progressive failure and associated fluid flow in rock are conducted and used to investigate the hydro-mechanical response of rock samples at laboratory scale. The responses investigated are the axial stress–axial strain together with permeability evolution and fracture patterns at various stages of loading. Then, the hydraulic fracturing process inside a rock specimen is numerically simulated. Three coupled processes are considered: (1) mechanical deformation of the solid medium induced by the fluid pressure acting on the fracture surfaces and the rock skeleton, (2) fluid flow within the fracture, and (3) propagation of the fracture. The numerically simulated results show that the fractures from a vertical wellbore propagate in the maximum principal stress direction without branching, turning, and twisting in the case of a large difference in the magnitude of the far-field stresses. Otherwise, the fracture initiates in a non-preferred direction and plane then turns and twists during propagation to become aligned with the preferred direction and plane. This pattern of fracturing is common when the rock formation contains multiple layers with different material properties. In addition, local heterogeneity of the rock matrix and macro-scale stress fluctuations due to the variability of material properties can cause the branching, turning, and twisting of fractures.     

天然地应力作用下不同产状岩体裂隙渗流特性研究

为了研究天然地应力作用下裂隙产状等因素对深部岩体裂隙渗流特性的影响,基于单裂隙面渗透性服从负指数变化规律,建立了三维应力作用下不同产状裂隙的渗透系数计算公式,利用Lagrange乘子法分析裂隙面产状变化对其渗透性的影响,并分析了岩体裂隙有、无充填物对其渗透性的影响及敏感性;然后,以我国大陆地区地应力统计规律为例,分析了地表以下5 000 m范围内在天然地应力作用下裂隙渗透性随深度、产状的变化规律。结果表明：裂隙产状的变化对其渗流特性有明显影响,对于浅层岩体,在大主应力大致呈水平方向分布时,随着裂隙面倾角的增加,裂隙渗透系数逐渐降低;但随着深度的增加,在裂隙深度超过约200 m和裂隙面走向与大或中主应力方向大致一致时,裂隙渗透性反而会随着裂隙面倾角的增加逐渐增加,在裂隙面走向与小主应力方向垂直时增加最为明显;对于深部岩体,裂隙的渗透性很小,裂隙面产状的变化对其渗透性影响很弱;对于有充填物裂隙,岩块与充填物的弹模比和充填物泊松比的变化对裂隙渗透性的影响很小。研究结果可为深入研究我国深部岩体渗透特性变化规律提供借鉴意义。     

Numerical simulation of fracture flow with a mixed-hybrid FEM stochastic discrete fracture network model

We introduce a discrete fracture network model of stationary Darcy flow in fractured rocks. We approximate the fractures by a network of planar circle disks, which is generated on the basis of statistical data obtained from field measurements. We then discretize this network into a mesh consisting of triangular elements placed in three-dimensional space. We use geometrical approximations in fracture planes, which allow for a significant simplification of the final triangular meshes. We consider two-dimensional Darcy flow in each fracture. In order to accurately simulate the channeling effect, we assign to each triangle an aperture defining its hydraulic permeability. For the discretization we use the lowest order Raviart-Thomas mixed finite element method. This method gives quite an accurate velocity field, which is computed directly and which satisfies the mass balance on each triangular element. We demonstrate the use of this method on a model problem with a known analytical solution and describe the generation and triangulation of the fracture network and the computation of fracture flow for a particular real situation.     

Numerical analysis of thermal fracturing in subsurface cold water injection by finite element methods

Integration of poromechanics and fracture mechanics plays an important role in understanding a series of thermal fracturing phenomena in subsurface porous media such as cold water flooding for enhanced oil recovery, produced‐water reinjection for waste disposal, cold water injection for geothermal energy extraction, and CO₂ injection for geosequestration. Thermal fracturing modeling is important to prevent the potential risks when fractures propagate into undesired zones, and it involves the coupling of heat transfer, mass transport, and stress change as well as the fracture propagation. Analytical method, finite element method, and finite difference method as well as boundary element method have been used to perform the thermal fracturing modeling considering different degrees and combinations of coupling. In this paper, extended finite element method is employed for the thermal fracturing modeling in a fully coupled fashion with remeshing avoided, and the stabilized finite element method is employed to account for the convection‐dominated heat transfer in the fracturing process with numerical oscillation circumvented. With the thermal fracturing model, a hypothetical numerical experiment on cold water injection into a deep warm aquifer is conducted. Results show that parameters such as injection rate, injection temperature, aquifer stiffness, and permeability can affect the fracture development in different ways and extended finite element method and stabilized finite element method provide effective tools for thermal fracturing simulation. Copyright © 2012 John Wiley & Sons, Ltd.     

Discrete fracture model for simulating waterflooding processes under fracturing conditions

In our study, we develop a model for simulating fracturing processes in a poroelastic medium. The proposed approach combines the discrete fracture model enriched with contact plane mechanics. The model captures mechanical interactions of fractures and a deformable medium, fluid, and heat transfer in fractures and in a porous medium. Both effects of poroelasticity and thermoelasticity are accounted in our model. Mass and heat conservation equations are approximated by the finite volume method, and mechanical equilibrium equations are discretized by means of the Galerkin finite element approach. Two‐dimensional grid facets between 3‐dimensional finite elements are considered as possible fracture surfaces. Most of these facets are inactive from the beginning and are activated throughout the simulation. A fracture propagation criterion, based on Irwin's approach, is verified on each nonlinear iteration. When the criterion is satisfied, additional contact elements are added into finite element and discrete fracture model formulations respectively. The proposed approach allows modeling of existing natural and artificially created fractures within one framework. The model is tested on single‐ and multiple‐phase fluid flow examples for both isothermal and thermal conditions and verified against existing semianalytical solutions. The applicability of the approach is demonstrated on an example of practical interests where a sector model of an oil reservoir is simulated with different injection and production regimes.     

A fracture mapping and extended finite element scheme for coupled deformation and fluid flow in fractured porous media

This paper presents a fracture mapping (FM) approach combined with the extended finite element method (XFEM) to simulate coupled deformation and fluid flow in fractured porous media. Specifically, the method accurately represents the impact of discrete fractures on flow and deformation, although the individual fractures are not part of the finite element mesh. A key feature of FM‐XFEM is its ability to model discontinuities in the domain independently of the computational mesh. The proposed FM approach is a continuum‐based approach that is used to model the flow interaction between the porous matrix and existing fractures via a transfer function. Fracture geometry is defined using the level set method. Therefore, in contrast to the discrete fracture flow model, the fracture representation is not meshed along with the computational domain. Consequently, the method is able to determine the influence of fractures on fluid flow within a fractured domain without the complexity of meshing the fractures within the domain. The XFEM component of the scheme addresses the discontinuous displacement field within elements that are intersected by existing fractures. In XFEM, enrichment functions are added to the standard finite element approximation to adequately resolve discontinuous fields within the simulation domain. Numerical tests illustrate the ability of the method to adequately describe the displacement and fluid pressure fields within a fractured domain at significantly less computational expense than explicitly resolving the fracture within the finite element mesh. Copyright © 2013 John Wiley & Sons, Ltd.     

A coupled compressible flow and geomechanics model for dynamic fracture aperture during carbon sequestration in coal

This paper presents the development of a discrete fracture model of fully coupled compressible fluid flow, adsorption and geomechanics to investigate the dynamic behaviour of fractures in coal. The model is applied in the study of geological carbon dioxide sequestration and differs from the dual porosity model developed in our previous work, with fractures now represented explicitly using lower-dimensional interface elements. The model consists of the fracture-matrix fluid transport model, the matrix deformation model and the stress-strain model for fracture deformation. A sequential implicit numerical method based on Galerkin finite element is employed to numerically solve the coupled governing equations, and verification is completed using published solutions as benchmarks. To explore the dynamic behaviour of fractures for understanding the process of carbon sequestration in coal, the model is used to investigate the effects of gas injection pressure and composition, adsorption and matrix permeability on the dynamic behaviour of fractures. The numerical results indicate that injecting nonadsorbing gas causes a monotonic increase in fracture aperture; however, the evolution of fracture aperture due to gas adsorption is complex due to the swelling-induced transition from local swelling to macro swelling. The change of fracture aperture is mainly controlled by the normal stress acting on the fracture surface. The fracture aperture initially increases for smaller matrix permeability and then declines after reaching a maximum value. When the local swelling becomes global, fracture aperture starts to rebound. However, when the matrix permeability is larger, the fracture aperture decreases before recovering to a higher value and remaining constant. Gas mixtures containing more carbon dioxide lead to larger closure of fracture aperture compared with those containing more nitrogen.     

An extended finite volume method and fixed-stress approach for modeling fluid injection–induced tensile opening in fractured reservoirs

Fluid injection–induced tensile opening is modeled using an extended finite volume method (XFVM). An embedded fracture strategy is used for the flow problem, that is, the fractures are discretized using finite volume segments without resolving the grid around them. Further, the discontinuities across fractures are modeled using special basis functions. The fracture openings due to enhanced fluid pressure and the associated shear slip due to traction free boundary condition on the fracture segments are both modeled using these special discontinuity basis functions. Mass transfer between fractures and matrix is modeled using the pressure difference. The enhancement of fracture storativity due to tensile opening leads to stronger coupling between flow and mechanics. An iterative scheme relying on the fixed-stress approach for fractures, which conserves the stress dependent terms over each iteration of the flow problem, has been introduced. Tensile opening has been simulated for single fractures embedded in two- and three-dimensional matrices. The convergence criterion for sequentially implicit fixed-stress scheme for fractures embedded in elastic media is established and has been validated numerically. Further, for 2D simulations, the effect of the matrix permeability for fracture propagation due to tensile opening has been studied.     

Modeling fluid injection in fractures with a reservoir simulator coupled to a boundary element method

We describe an algorithm for modeling saturated fractures in a poroelastic domain in which the reservoir simulator is coupled with a boundary element method. A fixed stress splitting is used on the underlying fractured Biot system to iteratively couple fluid and solid mechanics systems. The fluid system consists of Darcy’s law in the reservoir and is computed with a multipoint flux mixed finite element method, and a Reynolds’ lubrication equation in the fracture solved with a mimetic finite difference method. The mechanics system consists of linear elasticity in the reservoir and is computed with a continuous Galerkin method, and linear elasticity in the fracture is solved with a weakly singular symmetric Galerkin boundary element method. This algorithm is able to compute both unknown fracture width and unknown fluid leakage rate. An interesting numerical example is presented with an injection well inside of a circular fracture.     

Simulation of non‐planar three‐dimensional hydraulic fracture propagation

Hydraulic fracturing is the method of choice to enhance reservoir permeability and well efficiency for extraction of shale gas. Multi‐stranded non‐planar hydraulic fractures are often observed in stimulation sites. Non‐planar fractures propagating from wellbores inclined from the direction of maximum horizontal stress have also been reported. The pressure required to propagate non‐planar fractures is in general higher than in the case of planar fractures. Current computational methods for the simulation of hydraulic fractures generally assume single, symmetric, and planar crack geometries. In order to better understand hydraulic fracturing in complex‐layered naturally fractured reservoirs, fully 3D models need to be developed. In this paper, we present simulations of 3D non‐planar fracture propagation using an adaptive generalized FEM. This method greatly facilitates the discretization of complex 3D fractures, as finite element faces are not required to fit the crack surfaces. A solution strategy for fully automatic propagation of arbitrary 3D cracks is presented. The fracture surface on which pressure is applied is also automatically updated at each step. An efficient technique to numerically integrate boundary conditions on crack surfaces is also proposed and implemented. Strongly graded localized refinement and analytical asymptotic expansions are used as enrichment functions in the neighborhood of fracture fronts to increase the computational accuracy and efficiency of the method. Stress intensity factors with pressure on crack faces are extracted using the contour integral method. Various non‐planar crack geometries are investigated to demonstrate the robustness and flexibility of the proposed simulation methodology. Copyright © 2014 John Wiley & Sons, Ltd.     

沁水盆地柿庄北区块3号煤层裂缝预测

柿庄北地区位于沁水盆地东南部,是我国重要的煤层气产区。煤层裂缝作为主要的渗流通道,对煤层气高产富集起到重要控制作用。该区山西组3号煤层中褶皱轴部走向为近南北向,野外露头区发育一组共轭剪裂缝,优势方位为NW与近EW向,其次是近SN向。根据构造形迹确定燕山期应力场以SEE105°挤压为主。利用水力加砂压裂法对现今井点目的层地应力进行计算,通过古地磁定向与声速各向异性确定现今最大水平主应力方向为NE55°。将煤层构造与地表起伏形态作为主要影响因素,采用有限元法对燕山期及现今地应力场进行数值模拟,利用弹性力学理论对裂缝参数进行计算,计算结果认为柿庄北3号裂缝孔隙度、渗透率普遍较小,整体上中部背斜区 > 中部平缓区 > 东部单斜区 > 中部向斜区 > 西部地区。模拟的裂隙特征与实测数据所反映的裂缝特征相近,单井产气量高的地区,裂隙较发育。      

Identification of fracture development period and stress field analysis based on fracture fabrics in tectonic superposition areas

As a direct consequence of multiple periods of stress applied on areas with tectonic superposition, the multiple-periods fractures have complex abutting relationships, and the field study of fractures is usually restricted by outcrop conditions, such as section direction. Therefore, previous studies of superposed stress fields based on fractures have been generally performed in areas with proper observation conditions and clear abutting relationships. In contrast, in many other areas, the identification of fracture development period based on field observation is often infeasible. Compared to abutting relationships, fracture fabrics obtained from field measurement are not affected by the restriction of outcrops and consequently are more representative of the fractures. According to the analysis of fracture fabrics and fracture features, this paper has separated and extracted the superposed fracture sets and identified the fracture development period in the area without available abutting relationships. Taking the southern segment of the Longmen Mountain thrust belt as an example, fractures of two development periods are identified and timed in the tectonic superposition area between two adjacent fold belts. The analysis of stress direction in each period suggests that the structural boundaries, consisting of such pre-existing structures as faults and anticlines, could have induced directional rotation in the subsequent stress. An equivalent result was achieved using a finite element simulation of the stress field. Based on the stress analysis of the field sites and the stress field simulation, the stress variation in the tectonic superposition area is well modeled.     

渤中25-1低渗透油田地应力分布特征及对开发的影响

利用有限元方法, 对渤中25-1低渗透油田现今地应力状态进行数值模拟研究, 并结合低渗透储层地质特征, 分析地应力对油田开发的影响。数值模拟结果表明, 渤中25-1低渗透油田沙二和沙三段的现今地应力为压应力, 以北东东—南西西向为水平最大主应力的主要优势方位, 地应力在平面上受断层、构造起伏和沉积微相控制明显。根据地应力状态, 分析研究区的人工压裂裂缝为垂直裂缝。在远离断层的断块中部, 人工压裂缝的展布方向主要为北东东—南西西方向; 在断层附近, 由于水平最大主应力方向发生偏转, 人工压裂缝的延伸方向也会发生一定偏转, 部署开发井网时需要根据断层附近地应力的实际分布规律进行相应调整。      

Theoretical and numerical study of the steady‐state flow through finite fractured porous media

This paper investigates the two‐dimensional flow problem through an anisotropic porous medium containing several intersecting curved fractures. First, the governing equations of steady‐state fluid flow in a fractured porous body are summarized. The flow follows Darcy's law in matrix and Poiseuille's law in fractures. An infinite transversal permeability is considered for the fractures. A multi‐region boundary element method is used to derive a general pressure solution as a function of discharge through the fractures and the pressure and the normal flux on the domain boundary. The obtained solution fully accounts for the interaction and the intersection between fractures. A numerical procedure based on collocation method is presented to compute the unknowns on the boundaries and on the fractures. The numerical solution is validated by comparing with finite element solution or the results obtained for an infinite matrix. Pressure fields in the matrix are illustrated for domains containing several interconnected fractures, and mass balance at the intersection points is also checked. Copyright © 2013 John Wiley & Sons, Ltd.     

Numerical modeling of stress-permeability coupling in rough fractures

A numerical model is described for coupled flow and mechanical deformation in fractured rock. The mechanical response of rock joints to changes in hydraulic pressure is strongly influenced by the geometric characteristics of the joint surfaces. The concept of this work is to combine straightforward finite element solutions with complex and realistic fracture surface geometry in order to reproduce the non-linear stress-deformation-permeability coupling that is commonly observed in fractures. Building on the numerous studies that have expanded the understanding of the key parameters needed to describe natural rough-walled fractures, new methods have been developed to generate a finite element mesh representing discrete fractures with realistic rough surface geometries embedded in a rock matrix. The finite element code GeoSys/Rockflow was then used to simulate the coupled effects of hydraulic stress, mechanical stress, and surface geometry on the evolving permeability of a single discrete fracture. The modeling concept was experimentally verified against examples from the literature. Modeling results were also compared to a simple interpenetration model.     

复杂碳酸盐岩储层裂缝预测方法对比性研究

围绕着复杂碳酸盐岩储层裂缝的预测方法,以普光地区嘉陵江组为例,在统计研究区裂缝发育规律的基础上,分别利用构造曲率法、纵波宽方位各向异性法和有限元法对裂缝进行预测,比较各个方法预测结果和适用条件。对比结果表明:曲率法基于几何学理论,可以定性的分析裂缝发育强度,无法预测裂缝方位,该方法适用于基础资料较少、构造简单、褶皱发育的地区;纵波宽方位各向异性法是地球物理法预测裂缝的首选方法,可以直接反映有效裂缝的方位和密度,但是该方法适用于倾斜裂缝和高角度裂缝发育区,对地震资料质量和处理要求较高;有限元数值模拟法从正演角度再现了裂缝形成过程,计算分析裂缝的发育强度和方位,该方法需要对研究区裂缝形成机理认识充分,预测结果的精度取决于模型的精细程度。      

A shear‐dilation‐based model for evaluation of hydraulically stimulated naturally fractured reservoirs

The role of shear dilation as a mechanism of enhancing fluid flow permeability in naturally fractured reservoirs was mainly recognized in the context of hot dry rock (HDR) geothermal reservoir stimulation. Simplified models based on shear slippage only were developed and their applications to evaluate HDR geothermal reservoir stimulation were reported. Research attention is recently focused to adjust this stimulation mechanism for naturally fractured oil and gas reservoirs which reserve vast resources worldwide. This paper develops the overall framework and basic formulations of this stimulation model for oil and gas reservoirs. Major computational modules include: natural fracture simulation, response analysis of stimulated fractures, average permeability estimation for the stimulated reservoir and prediction of an average flow direction. Natural fractures are simulated stochastically by implementing ‘fractal dimension’ concept. Natural fracture propagation and shear displacements are formulated by following computationally efficient approximate approaches interrelating in situ stresses, natural fracture parameters and stimulation pressure developed by fluid injection inside fractures. The average permeability of the stimulated reservoir is formulated as a function of discretized gridblock permeabilities by applying cubic law of fluid flow. The average reservoir elongation, or the flow direction, is expressed as a function of reservoir aspect ratio induced by directional permeability contributions. The natural fracture simulation module is verified by comparing its results with observed microseismic clouds in actual naturally fractured reservoirs. Permeability enhancement and reservoir growth are characterized with respect to stimulation pressure, in situ stresses and natural fracture density applying the model to two example reservoirs. Copyright © 2002 John Wiley & Sons, Ltd.     

Modelling fluid flow through a single fracture with different contacts using cellular automata

The representation of the surface of a rock fracture and a numerical method to simulate fluid flow in single fractures are the keys to understanding the hydraulic behaviour of rock fractures. In this paper, a cellular automaton (CA) approach is used to generate the single fracture structure, which is assumed to be composed of contacts and voids. We develop a CA evolution rule to produce a contact area, and randomly model a single rock fracture with different contact ratios to reflect natural fracture properties such as dead voids, islands and tortuous flow path. Then, based on the localisation theory of a CA, a numerical method to simulate fluid flow in single fractures with contacts is developed. In this method, the fracture is discretised into a system composed of cell elements. Different apertures, i.e., zero for contacts and non-zero for voids, are assigned to each cell element. Therefore, the contribution of the cell elements in a contact on a cell’s transmissivity can be ignored completely. The local transmissivity is assumed to conform to the cubic law. The fluid flow in a fracture with different contact situations is then modelled using the method established in this paper. The fluid flow path, flow velocity and fluid head distributions as well as the channel flow in the fracture are well-modelled. The flow behaviour of the fracture strongly depends on the effective fluid flow path.     

Modeling hydraulic fracture propagation in permeable media with an embedded strong discontinuity approach

The paper presents an embedded strong discontinuity approach to simulate single hydraulic fracture propagation in the poroelastic medium under plane-strain conditions. The method enriches the strain field with the discontinuous deformation mode and allows the fracture to be modeled inside elements. The Mode-I fracture initiation and propagation are described by the trilinear cohesive law, which is implemented by the penalty method. The enhanced permeability inside the fractured elements is dependent on the fracture aperture. Hydraulic fracture propagation is driven by the high pressure gradient near the fracture. Fluid transfer between the fracture and bulk rock is automatically captured within the poroelastic framework. The numerical framework is verified by the comparisons with the asymptotic analytical solutions for single hydraulic fracture propagation.     

塔河油田12区块奥陶系裂缝分布规律研究

通过对露头、铸体薄片、岩性分析和成像测井资料观察和分析发现,塔河油田12区块奥陶系储层内缝洞发育,并且裂缝在一间房和鹰山组的储层中起了至关重要的作用,它既是地层流体的存储空间,又是连通孔洞缝的主要通道。裂缝的发育也增加了储层的非均质性。所以,精细研究裂缝的发育和分布规律,并定量地计算裂缝的相关参数对于识别奥陶系储层和对...