A modeling approach for analysis of coupled multiphase fluid flow, heat transfer, and deformation in fractured porous rock

This paper presents the methodology in which two computer codes—TOUGH2 and FLAC3D—are linked and jointly executed for coupled thermal–hydrologic–mechanical (THM) analysis of multiphase fluid flow, heat transfer, and deformation in fractured and porous rock. TOUGH2 is a well-established code for geohydrological analysis with multiphase, multicomponent fluid flow and heat transport, while FLAC3D is a widely used commercial code that is designed for rock and soil mechanics with thermomechanical and hydromechanical interactions. In this study, the codes are sequentially executed and linked through external coupling modules: one that dictates changes in effective stress as a function of multi-phase pore pressure and thermal expansion, and one that corrects porosity, permeability, and capillary pressure for changes in stress. The capability of a linked TOUGH-FLAC simulator is demonstrated on two complex coupled problems related to injection and storage of carbon dioxide in aquifers and to disposal of nuclear waste in unsaturated fractured porous media.     

Assessment of an equivalent porous medium for coupled stress and fluid flow in fractured rock

The equivalent porous medium (EPM) for coupled stress and fluid flow in fractured rock was assessed. The assessment was focused on the contributions of the seepage-induced force to the system's equilibrium and deformation as well as fluid flux through the fractures, with the special emphasis on the inter-relations among these three aspects. Seepage-induced average stress was calculated for the EPM. And it was proved, with respect to the fluid flow system (FFS) proposed by Long et al. [Long JCS, Remer JC, Wilson CR, Witherspoon PA. Porous medium equivalents for networks of discontinuous fractures. Water Resour Res 1982;18:645–658], that the equivalence of the equilibrium contribution of the seepage-induced force is unconditionally satisfied in the EPM. The equivalence of fluid flux and seepage-induced deformation was evaluated using numerical techniques based on the FFS, in which the seepage-induced deformation was considered as ‘tensile’ deformation under fluid pressures. The numerical simulation results of the two cases suggest that the EPM established in terms of fluid flux based on the FFS guarantees equivalence of coupled stress and fluid flow. Also, the numerical simulation results indicate that the percolation theory can not be applied to a system with fractures of remarkably different size and hydraulic conductivity.     

裂隙岩体饱和-非饱和渗流研究进展

本文从裂隙岩体非饱和渗流特性及参数、裂隙岩体浸润面的确定方法、渗流数学模型、工程应用研究等方面对高坝工程裂隙岩体渗流研究现状进行了系统的分析总结。当前,主要沿用多孔介质渗流的方法研究裂隙渗流,所建立的裂隙岩体非饱和渗流水力参数关系还不能真实反映裂隙渗流的特性;连续介质-非连续介质渗流模型能够较好反映裂隙岩体的空隙结构和渗流特性,要继续研究相应的数值计算方法;库水位升降变化及降雨影响下库岸边坡岩体渗流是饱和-非饱和非稳定渗流,应将饱和区与非饱和区统一起来分析。     

Estimation of REV size and three-dimensional hydraulic conductivity tensor for a fractured rock mass through a single well packer test and discrete fracture fluid flow modeling

A new methodology is presented to determine the representative elementary volume (REV) size and three-dimensional (3-D) hydraulic conductivity tensor for a fractured rock mass. First, a 3-D stochastic fracture network model was built and validated for a gneissic rock mass based on the fracture data mapped from scanline surveys at the site. This validated fracture network model was combined with the fracture data observed on a borehole to generate a stochastic-deterministic fracture network system in a cubic block around each packer test conducted at a different depth region in the same borehole. Each packer test was simulated numerically applying a developed discrete fracture fluid flow model to estimate the influenced region or effective range for the packer test. A cubic block of size 18 m, with the packer test interval of length about 6.5 m located at the centre of this block, was found to be suitable to represent the influenced region. Using this block size, the average flow rate per unit hydraulic gradient (defined as the transmissivity multiplied by mean width of flow paths) field for fractures was calibrated at different depth regions around the borehole by numerically simulating the packer tests conducted at different depth regions. The average flow rate per unit hydraulic gradient of the fractures that intersect the borehole was considered to be quite different to the average flow rate per unit hydraulic gradient of the fractures that do not intersect the borehole. A relation was developed to quantify the ratio between these two parameters. By studying the directional hydraulic conductivity behaviour of different cubic block sizes having the validated stochastic fracture network and calibrated hydraulic parameters, a REV for the hydraulic behaviour of the rock mass was estimated to be a block size of 15 m. The hydraulic conductivity tensor in 3-D computed through regression analysis using the calculated directional hydraulic conductivity values in many directions was found to be significantly anisotropic. The principal directions of the hydraulic conductivity tensor were found to be agreeable with the existing fracture system of the site. Further, the geometric hydraulic conductivity calculated was found to be comparable to the hydraulic conductivity estimated through the radial flow assumption in continuum porous media.     

Stress and water flow relation in fractured rock

Bulletin of Engineering Geology and the Environment - A field test was performed to investigate the relation between rock stresses and the hydraulic properties of a fractured rock mass. In thisin...     

Continuous approach for coupled mechanical and hydraulic behavior of a fractured rock mass during hypothetical shaft sinking at Sellafield, UK

An indicator simulation technique was used to infer the continuous hydraulic conductivity distribution in the BVG for Task 1A of DECOVALEX II. Young's modulus from the wireline geophysical survey was used as soft data. To calibrate the model with measured data in Task 1B, a back analysis with Bayesian estimate was applied to determine the hydraulic conductivity of all the elements of the finite element mesh. To avoid the ill-posed problem, a geostatistically kriged hydraulic conductivity distribution was used as mean a priori information. While there are some remaining problems, the solution was improved in comparison with the ordinary analysis. For the shaft sinking problem defined as Task 1C, the continuous approach to coupled mechanical and hydraulic behavior was used by combining crack tensor theory and the Barton and Bandis model for rock joints. The mechanical and hydraulic properties were derived with resultant anisotropy and heterogeneity during the simulated shaft excavation. It is found that the calculated hydraulic conductivity becomes larger than measured one as a consequence of excavation.     

矿山岩体破坏突水非达西流模型及数值求解

 为模拟突水流体瞬态流动全过程并揭示突水流体流动机制,针对矿山岩体破坏突水过程中流体的高速非达西流特征和流场动力学统一性特征,基于流体质量守恒和压力平衡原理建立耦合Darcy方程、Forchheimer方程和Navier-Stokes方程的矿山岩体破坏突水非达西流模型,将含水层水源、破碎岩体突水通道和采场巷道出口整个水流路径有机联系在统一流动场中;并基于有限元弱形式和有限体积法耦合积分方程,提出有限元和有限体积法相结合的数值计算方法,应用FEPG有限元软件编译FORTRAN源程序,模拟突水瞬态流动全过程。通过算例对比达西和非达西流结果表明,对于矿山岩体破坏突水问题采用非达西流模型计算十分必要,破碎岩体作为导水通道沟通含水层和巷道的水力联系,导水通道内高速水流的惯性作用是含水层充足的补给水量沿破碎带进入巷道形成突水的主要原因。     

A numerical study on differences in using Navier–Stokes and Reynolds equations for modeling the fluid flow and particle transport in single rock fractures with shear

The study on fluid flow and transport processes of rock fractures in most practical applications involves two fundamental issues: the validity of Reynolds equation for fluid flow (as most often assumed) and the effects of shear displacements on the magnitudes and anisotropy of the fluid flow velocity field. The reason for such concerns is that the impact of the surface roughness of rock fractures is still an unresolved challenging issue. The later has been systematically investigated with results showing that shear displacement plays a dominant role on evolutions of fluid velocity fields, for both magnitudes and anisotropy, but the former has not received examinations in details due to the numerical complexities involving solution of the Navier–Stokes (NS) equations and the representations of fracture geometry during shear. The objective of this paper aims to solve this problem through a FEM modeling effort.Applying the COMSOL Multiphysics code (FEM) and assuming a 2D problem, we consider the coupled hydromechanical effect of fracture geometry change due to shear on fluid flow (velocity patterns) and particle transport (streamline/velocity dispersion), using measured topographical data of natural rock fracture surfaces. The fluid flow in the vertical 2D cross-sections of single rock fractures was simulated by solving both the Navier–Stokes and the Reynolds equation, and the particle transport was predicted by the streamline particle tracking method with calculated flow velocity fields (vectors) from the flow simulations, obtaining results such as flow velocity profiles, total flow rates, particle travel time, breakthrough curves and the Péclet number, Pe, respectively.The results obtained using NS and Reynolds equations were compared to illustrate the degree of the validity of the Reynolds equation for general applications in practice since the later is mush more computationally efficient for solving large-scale problems.The flow simulation results show that both the total flow rate and the flow velocity fields in a rough rock fracture predicted by the NS equation were quite different from those predicted by the Reynolds equation. The results show that a roughly 5–10% overestimation on the flow rate is produced when the Reynolds equation is used, and the ideal parabolic velocity profiles defined by the local cubic law, when Reynolds equation is used, is no longer valid, especially when the roughness feature of the fracture surfaces changes with shear. These deviations of flow rate and flow velocity profiles across the fracture aperture have a significant impact on the particle transport behavior and the associated properties, such as the travel time and Péclet number. The deviations increase with increasing flow velocity and become more significant when fracture aperture geometry changes with shear.     

Numerical determination of the equivalent elastic compliance tensor for fractured rock masses using the distinct element method

The purpose of this paper is to establish a methodology to determine the equivalent elastic properties of fractured rock masses by explicit representations of stochastic fracture systems, and to investigate the conditions for the application of the equivalent continuum approach for representing mechanical behavior of the fractured rock masses. A series of numerical simulations of mechanical deformation of fractured rock masses at different scales were conducted with a large number of realizations of discrete fracture networks (DFN), based on realistic fracture system information and using the two-dimensional distinct element program, UDEC. General theory of anisotropic elasticity was used for describing the macroscopic mechanical behavior of fractured rock masses as equivalent elastic continua. Verification of the methodology for determining the elastic compliance tensor was conducted against closed-form solutions for regularly fractured rock mass, leading to very good agreements. The main advantage of the developed methodology using the distinct element method is that it can consider complex fracture system geometry and various constitutive relations of fractures and rock matrix, and their interactions. Two criteria for the applicability of equivalent continuum approach were adopted from the investigations: (i) the existence of a properly defined REV (representative elementary volume) and (ii) the existence of an elastic compliance tensor. For the problems with in situ conditions studied in this paper, the results show that a REV can be defined and the elastic properties of the fractured rock mass can be represented approximately by the elastic compliance tensor through numerical simulations.     

Stochastic analysis of flow response in a three-dimensional fractured rock mass block

The pressure variations observed during an interference test performed in a fractured volcanic tuff are modeled using a stochastic continuum approach. Two sets of spatially heterogeneous conductivity fields are generated. The first one is only conditioned to hard and soft data on conductivities, whereas the second one is also conditioned to transient pressure information. Conditioning to the interference tests introduces an important modification in the conductivity fields, reducing the presumed contrast between the faulted blocks and those that represent the background fracturing. It was possible to reproduce most of the pressure responses simply using a heterogeneous distribution of conductivities without the need to resort to a hydro-mechanical process coupling.     

Numerical modeling of closure effect of natural fracture surfaces of rock on behavior of fluid flow

Bulletin of Engineering Geology and the Environment - It is essential to understand the hydraulic behavior of rock apertures as the fluid flow path in the rock mass. The geometric parameters of...     

Numerical manifold method modeling of coupled processes in fractured geological media at multiple scales

The greatest challenges of rigorously modeling coupled hydro-mechanical(HM) processes in fractured geological media at different scales are associated with computational geometry.These challenges include dynamic shearing and opening of intersecting fractures at discrete fracture scales as a result of coupled processes,and contact alteration along rough fracture surfaces that triggers structural and physical changes of fractures at micro-asperity scale.In this paper,these challenges are tackled by developing a comprehensive modeling approach for coupled processes in fractured geological media based on numerical manifold method(NMM) at multiple scales.Based on their distinct geometric features,fractures are categorized into three different scales:dominant fracture,discrete fracture,and discontinuum asperity scales.Here the scale is relative,that of the fracture relative to that of the research interest or domain.Different geometric representations of fractures at different scales are used,and different governing equations and constitutive relationships are applied.For dominant fractures,a finite thickness zone model is developed to treat a fracture as a porous nonlinear domain.Nonlinear fracture mechanical behavior is accurately modeled with an implicit approach based on strain energy.For discrete fractures,a zero-dimensional model was developed for analyzing fluid flow and mechanics in fractures that are geometrically treated as boundaries of the rock matrix.With the zero-dimensional model,these fractures can be modeled with arbitrary orientations and intersections.They can be fluid conduits or seals,and can be open,bonded or sliding.For the discontinuum asperity scale,the geometry of rough fracture surfaces is explicitly represented and contacts involving dynamic alteration of contacts among asperities are rigorously calculated.Using this approach,fracture alteration caused by deformation,re-arrangement and sliding of rough surfaces can be captured.Our comprehensive model is able to handle the computational challenges with accurate representation of intersections and shearing of fractures at the discrete fracture scale and rigorously treats contacts along rough fracture surfaces at the discontinuum asperity scale.With future development of three-dimensional(3 D) geometric representation of discrete fracture networks in porous rock and contacts among multi-body systems,this model is promising as a basis of 3 D fully coupled analysis of fractures at multiple scales,for advancing understanding and optimizing energy recovery and storage in fractured geological media.     

裂隙岩体渗流场与应力场耦合研究进展与展望

在分析岩体多介质渗透特性的基础上，提出了岩体多介质渗透类型，综述了裂隙岩体渗流场与应力场耦合模型的研究进展，并提出耦合模型进一步研究的方向。     

裂隙岩体渗流溶质运移耦合离散裂隙模型数值计算方法

研究裂隙岩体渗流溶质运移问题对于岩土工程地下水污染物预测控制具有重要意义。基于离散裂隙网络模型,采用实体单元模拟基质岩块、无厚度单元模拟复杂裂隙网络,提出了裂隙岩体渗流溶质运移耦合的三维数值计算方法。针对无反应项和含反应项两种情况,通过算例分析了单裂隙中溶质迁移行为,并与精细模拟方法、解析方法的结果进行对比验证;进一步将该法应用于预测大规模裂隙岩体溶质浓度分布规律及发展趋势,并评价了主要影响因素。结果表明,该法可有效模拟裂隙网络、基质岩块中水分溶质传输行为;由于贯通裂隙网络的优势流影响,溶质羽主要受控于裂隙水的对流作用,出现了高度非均匀分布现象;通过参数敏感性分析发现,相较于岩块基质的扩散作用,裂隙开度产生的对流作用是影响浓度场分布的主控因素。在保证精度的前提下,该法可大幅减小计算量和计算时长,对于解决含复杂裂隙网络岩体渗流传质的三维数值模拟问题具有明显优势。     

Coupled fluid-flow and geomechanics for triple-porosity/dual-permeability modeling of coalbed methane recovery

A coupled fluid-flow and geomechanics model for simulating coalbed methane (CBM) recovery is presented. In the model, the fluid-flow process is simulated with a triple-porosity/dual-permeability model, and the coupling effects of effective stress and micro-pore swelling/shrinkage are modeled with the coupled fluid-flow and geomechanical deformation approach. The mathematical model is implemented with a finite volume method. First, a case without considering coupling between fluid-flow and geomechanics is simulated and compared with an existing simulator. The effects of coupled fluid-flow and geomechanics are then studied in detail with two illustrative examples. The first one is designed for testing the effective stress effect without micro-pore swelling/shrinkage effect, and the other for testing the coupling effects of the effective stress and micro-pore swelling/shrinkage on the methane production. The numerical results indicate that both the effective stress and the micro-pore shrinkage make a significant contribution to fluid-flow in CBM reservoir and to methane production. The methane production sensitivity to Young’s modulus and Langmuir sorption strain are investigated as well. Finally, we make a dynamic analysis of the coupling effects of fluid-flow process and geomechanics.     

A fully coupled thermo-hydro-mechanical model for simulating multiphase flow,deformation and heat transfer in buffer material and rock masses

This paper presents a numerical method for modeling coupled thermo-hydro-mechanical processes of geomaterials with multiphase fluid flow. A FEM code has been developed and validated for modeling the behavior of porous geological media, and is equally applicable for modeling coupled THM processes in rocks. The governing equations are based on the theory of mixtures applied to the multiphysics of porous media, considering solid phase deformation, multiphase fluid flow, and heat transport. New numerical techniques have been developed for more efficient FEM formulation and equation solution for modeling saturated or partially saturated water flow, gas flow and heat transfer in deformable porous media, as are commonly encountered in performance and safety assessment of underground radioactive repositories. The code has been validated against an experimental benchmark test, which involves bentonite under laboratory conditions, with good results. Several critical outstanding issues for modeling coupled processes of geomaterials are discussed in depth.     

Thermal response testing of a fractured hard rock aquifer with and without induced groundwater flow

This study quantifies experimentally the influence of groundwater on the thermal conductivity measurements via thermal response tests (TRT) in a fractured hard rock with low matrix permeability. An artificial groundwater flow towards a nearby well was induced by groundwater extraction and a TRT performed simultaneously. The results were compared with a TRT performed 24?days later in the same well without groundwater extraction. The effect of the groundwater flow is shown indirectly by the enhanced effective thermal conductivity and directly through the comparison of temperature profiles before and 4?h after both TRTs. Simulations in FEFLOW show that a groundwater flow velocity of 130–1,300?m?d^?1 through one open horizontal fracture of small volume increases the effective thermal conductivity by 11?% in the studied system.     

A thermo-hydro-mechano-chemical formulation for modeling water transport around a ventilated tunnel in an argillaceous rock

Argillaceous rocks are being considered as potential host rocks for deep geological disposal. For the research work in DECOVALEX-2011, 5 participant research teams performed simulations of a laboratory drying test and a ventilation experiment for Mont Terri underground laboratory built in argillaceous rock formation. Our study starts with establishing a coupled thermo-hydro-mechano-chemical (THMC) processes model to simulate water transport in rock around the ventilated tunnel. Especially in this THMC formulation, a three-phase and two-constituent hydraulic model is introduced to simulate the processes which occur during tunnel ventilation, including desaturation/resaturation in the rock, phase change and air/rock interface, and to explore the Opalinus clay parameter set. It can be found that water content evolution is very sensitive to intrinsic permeability, relative permeability and capillary pressure in clay rock. Water loss from surrounding rock is sensitive to the change of permeability in clay which is induced by excavation damaged zone. Chemical solute transport in the rock near ventilation experiment tunnel is simulated based on the coupled THMC formulation. It can be estimated that chemical osmotic flow has little significance on water flow modeling. Comparisons between simulation results from 5 teams and experimental observations show good agreement. It increases the confidence in modeling and indicates that it is a good start for fully THMC understanding of the moisture transportation and mechanical behavior in argillaceous rock.