Modeling the orientation- and stress-dependent permeability of anisotropic rock with particle-based discrete element method

Abstract

A novel physics-based numerical model is proposed to simulate the orientation and effective confining pressure dependent permeability of anisotropic rock. In the two-dimensional discrete element model, presence of anisotropy is explicitly represented by inserting a set of smooth joints. Based on the experimentally obtained effective stress law at sample scale, a physics-based relation is deduced to describe the reduction of pipe aperture upon the normal contact force at grain scale. Darcy test is conducted to validate the model by comparing the flow rate and pressure distributions with analytical solutions. Different parameters are assigned to represent the difference in flow capacity of rock matrix and beddings. Fluid flow tests performed on the isotropic model and anisotropic models with horizontal and vertical beddings reveal that the macro permeability decreases with increasing effective stress, following the same effective stress law. The initial aperture dominates the intrinsic permeability while the reduction of permeability is due to the closure of pipe aperture. Permeability anisotropy is caused by the different apertures assigned to the rock matrix and the bedding while the force sensitivity factors determine the stress-dependence of the permeability anisotropy. Simulations of the stress concentration and fluid dissipation around borehole confirm the capacity of the model in capturing the hydro-mechanical coupled responses of anisotropic rock formation. This study provides a fluid flow model for the exploration of mechanisms underlying the orientation and stress dependent permeability of anisotropic rocks and for the simulation of their engineering responses subjected to hydro-mechanical coupling.

Introduction

Anisotropy is one of the most prominent features characterizing geomaterials. Many rocks encountered by human activities show clear geological structures such as bedding, stratification, layering, foliation or fissuring.¹ Such rocks are said to be anisotropic as their physical, mechanical and hydraulic properties vary with direction. These structures may act as weak planes and preferential flow paths in many engineering applications, and lead to the preferred stress concentration and initiation of fractures.

Rock permeability is a key physical parameter describing the ability of fluid to move within a rock mass. Permeability of anisotropic rock is of great importance in many geo-mechanical and geological processes, including the borehole instability, long-term production from gas shale reservoirs, hydraulic fracturing treatment, CO₂ storage, fluid percolation in the crust, and fluid flow in a mantle wedge.2, 3, 4 For instance, depletion of fluid pressure during hydrocarbon exploration may raise the effective stress and thus decline the bulk permeability of shale reservoir. This process has significant implications on the analysis and forecasting of production data.⁵ Therefore, accurate estimation and representation of the permeability of anisotropic rocks are of great importance.

Fluid transmitting in anisotropic rocks is a complicated process which might be affected by many factors including the effective stress, the flow direction with respect to beddings, and the preferred orientation of minerals.⁶ So far, numerous fluid flow experiments have been conducted on shales,5, 6, 7 sandstones,^3,8 limestones⁹ and siliciclastic caprock,¹⁰ to elucidate the variation of permeability subjected to various stress and flow conditions. In these tests, the transport of gases (e.g., argon gas, CO2 and nitrogen gas),^5,7,10,11 water^4,9 and NaCl solution^6,12 are measured by either the steady state flow method^3,9 or the transient pulse decay method.^4,7,10 Although the magnitude of measured permeability differs by several orders, in general, the permeability decreases irreversibly with externally applied stress for a specific rock.13, 14, 15 Various relations have been proposed on an empirical basis to describe the reductions in permeability with increasing effective pressure, known as the effective stress law, among which, the cubic law,^6,12 the exponential law,^4,10 and the power law^3,9 are the commonly used.

Fluid flow experiments also revealed a significant permeability anisotropy in anisotropic rocks, that is, the magnitude of permeability varies with flow direction relative to bedding.¹⁰ In general, the permeability parallel to bedding is initially higher than that normal to bedding.¹⁶ Kwon et al.⁶ found that the permeability of illite‐bearing shale becomes increasing isotropic, showing little directional dependence at the effective pressure of 10–12 MPa. However, Bahandari et al.⁵ reported that both the horizontal and vertical permeability of Barnett shale are equally sensitive to stress, that is, the permeability anisotropy ratio does not vary with effective stress.⁵ Experiments conducted on the sheared serpentinite reveals that the ratio between permeability parallel to foliation and perpendicular to foliation is around one under low confining pressure, but increases to several orders under a confining pressure of 50 MPa.⁴ Although possible mechanisms were provided to explain these phenomena,⁴ it is difficult to find a definite answer through laboratory experiments due to the limited test methods and sample heterogeneity. The details of fluid flow processes remain poorly understood and the mechanisms contributing to these discrepancies are still unclear.

Numerical simulations provide an essential way to study the micro-mechanism underlying many geomechanical processes.17, 18, 19, 20 In particular, the discrete element method (DEM) can explicitly represent the micro-structure and physical processes of rock by combining simple mechanisms of granular deformation, cracking and pore channel fluid flow. DEM has been extensively applied in the modeling of rock mechanical properties and has been proved to be able to capture many of the dominant physical processes with sufficient completeness.21, 22, 23, 24, 25, 26 Regarding the simulation of anisotropic rocks, efforts have been made to represent the stratified structure at multi-scale, by inserting individual smooth joints,^27,28 a set of parallel beddings,^29,30 and by combining the Distinct Fracture Network.^31,32

DEM has also been coupled with lattice Boltzmann method (LBM) and the computational fluid dynamics (CFD) to simulate thermally induced damage of rocks,³³ and fluid flow in fractures.³⁴ Regarding the fluid flow within rock matrix, simulation is commonly implemented through the pipe-network algorithm,35, 36, 37 which has been applied in the simulation of stress-induced permeability change,²⁴ hydraulic fracture propagation^38,39 and gas invasion in sediments.⁴⁰ Nevertheless, there is a challenge to develop a physics-based numerical approach to explicitly represent the existence of weak layers and to explore the fundamental process contributing to permeability alternation in anisotropic rocks.33, 34, 35, 36

In this study, a novel DEM model is developed to simulate the orientation and stress dependent permeability of anisotropic rock. The presence of layered structures is explicitly represented by inserting a set of smooth joints into the bonded particle model. Based on the experimentally obtained effective stress law at macro scale, a new physics-based relation is developed to describe the reduction of pipe aperture versus normal contact force at the grain scale. Validation of the model is examined by comparing the Darcy test results with analytical solutions. Effective confining pressure dependence of the permeability and the permeability anisotropy are investigated through a set of parametric studies. The model is also applied in the investigation of stress concentration and fluid transmission in the borehole system. The fundamental processes underlying permeability alternation are explored at the grain scale.