Effective Elastic and Hydraulic Properties of Fractured Rock Masses with High Contrast of Permeability: Numerical Calculation by an Embedded Fracture Continuum Approach

Dang, Hong-Lam; Do, Duc Phi; Hoxha, Dashnor

doi:https://doi.org/10.1155/2019/7560724

Advances in Civil Engineering

On this page

Abstract Introduction Conclusions Data Availability Conflicts of Interest References Copyright Related Articles

Research Article | Open Access

Volume 2019 | Article ID 7560724 | https://doi.org/10.1155/2019/7560724

Effective Elastic and Hydraulic Properties of Fractured Rock Masses with High Contrast of Permeability: Numerical Calculation by an Embedded Fracture Continuum Approach

Hong-Lam Dang,^1,2Duc Phi Do,¹and Dashnor Hoxha¹

Academic Editor: Stefano de Miranda

Received28 Aug 2018

Revised06 Dec 2018

Accepted19 Dec 2018

Published03 Feb 2019

Abstract

In this work, the hydromechanical modeling of the fractured rock masses was conducted based on a new numerical simulation method named as embedded fracture continuum (EFC) approach. As the principal advantage, this approach allows to simplify the meshing procedure by using the simple Cartesian meshes to model the fractures that can be explicitly introduced in the porous medium based on the notion of fracture cells. These last elements represent the grid cells intersected by at least one fracture in the medium. Each fracture cell in the EFC approach present a continuum porous medium whose hydromechanical properties are calculated from ones of the matrix and ones of the intersected fractures, thanks for using the well-known solution of the joint model. The determination of the hydromechanical properties of the fracture cells as presented in this work allows to provide the theoretical base and to complete some simple approximations introduced in the literature. Through different verification tests, the capability of the developed EFC approach to model the hydromechanical behavior of fractured rock was highlighted. An analysis of different parameters notably the influence of the fracture cell size on the precision of the proposed approach was also conducted. This novel approach was then applied to investigate the effective permeability and elastic compliance tensor of a fractured rock masses taken from a real field, the Sellafield site. The comparison of the results calculated from this approach with ones conducted in the literature based on the distinct element code (UDEC) presents a good agreement. However, unlike the previous studies using UDEC, which limits only in the case of fractured rock masses without dead-end fractures, our approach allows accounting for this kind of fractures in the medium. The numerical simulations show that the dead-end fractures could have a considerable contribution on the effective compliance moduli, while their effect can be neglected to calculate the overall permeability of the of fractured rock masses.

1. Introduction

Evaluation of effective mechanical and hydraulic properties of a geological reservoir, represented as porous rock masses with complex fracture networks, is a major issue in many engineering applications such as oil and gas field recovery, CO₂ storage, and geothermal exploitation to just mention a few. The porosity of the matrix constitutes the primary porosity and the essential storage capacity of the reservoir, while the capacity storage of fractures represents only some fraction of the total porosity. Likewise, the hydraulic conductivity of such fractured reservoirs is composed by the conductivity of the rock matrix (pores and their connectivity) and the conductivity of cracks/fractures; the later one is usually some orders of magnitude greater than the former. Thus, the presence of fractures induces in the reservoir a quite complex fluid flow and hence significantly affects the production rate of oil recovery or geothermal energy. The idea of dissociation of flow inside the fracture network from that of the matrix is first proposed in [1, 2]. Since then, the double permeability system is still used as a common nomination for fractured reservoirs.

Albeit a huge number of studies on fractured rock masses in the past, their modeling is yet a challenging and dynamic research topic, with two principal trends grouping continuum and discontinuum approaches [3–10]. In the former approach, the fractured medium is considered as an equivalent continuum medium whose properties are derived from an appropriate upscaling procedure by using a homogenization technique. For this approach, the effect of fractures is implicitly accounted for in the equivalent constitutive model that can be described by the principles of continuum mechanics: the properties of the equivalent continuum medium are homogenous and defined at any point of the domain. This approach has the advantage for solving problems of large scales with the high densities of fractures, but its results may be sensitive to the domain’s size especially when the studied domain is smaller than the representative elementary volume [10]. Otherwise, the homogenization-based continuum model may not take into account the connectivity effect of fractures, the effect of clustering and spatial distribution of fractures, or the individual characteristics of a fracture since these characteristics are, at best, taken into account through some statistical features of all crack sets.

For the discontinuous approach, the medium with an embedded discrete fracture network (DFN) is represented as an assemblage of blocks (discrete elements) bounded by a number of intersecting discontinuities [3–7, 10, 11]. The overall behavior of the fractured medium is mainly affected by the interactions between the matrix blocks and fractures (represented by interfaces between blocks). Although this approach can better describe the behavior of the fractured medium on a small scale, it can be very expensive on computer memory and time simulation particularly for a large scale and/or longtime transient problems. The main issue for this approach is the reproduction by numerical models of the complexity of in field fractures. In fact, because of the complexity and variation of fracture properties at a given site, the statistics of geometric characteristics and the properties of fractures/cracks are issued from limited and potentially biased field measurements (e.g., one-dimensional (1D) borehole imaging or two-dimensional (2D) outcrop mapping). The challenge here is not only how to input a great number of fractures into the model but also the conceptual model to express the fracture network and its connectivity.

To overcome limitations, as well as to take advantage of these two-mentioned approaches, some hybrid approaches borrowing the concept of continuum models incorporating the effect of fractures explicitly were introduced under different names: “fracture continuum approach” [9, 12–17], “fracture-cell model” [18], or “embedded discrete fracture model” [19–24]. The common point and the principal idea of these approaches lies on the concept of “fracture cell” which represents the grid cell of the discretized medium intersected by one or many fractures in the medium. Each fracture cell presents an equivalent porous medium with its own properties calculated from the contributed properties of intact matrix and fractures intersecting the cell. Thus, the fracture cell properties are different from that of the intact matrix and vary from cell to cell in respect with the geometry and properties of intersecting fractures. This approach was successfully applied to the study flow problem in fractured rocks [9, 12–16, 18] or the coupled flow and heat transfer process [17] in the fractured reservoir. Some applications were recently conducted for the coupled hydromechanical problem [21–24].

In this work, a hybrid method, called from now on as the embedded fracture continuum (EFC) approach, is chosen to study the coupled hydromechanical behavior of the fractured porous media. The proposed approach can be applied for all fractured media owing a weak to a high density of fracture. But for the demonstration purpose, in the present work, we focused notably in the case of high density with an application in a real fractured rock of the well-known Sellafield site. The work is limited in the case of single phase flow. For some approaches on multiphase flow, the interested reader could consult [25, 26].

The structure of this paper is organized as follows. Following this introduction, the EFC approach based on the notion of the fracture cell with its equivalent (poroelastic and hydraulic) properties will be detailed. The implementation of this approach in the open source code deal.II [27] is validated through different numerical verification tests. As the first application, this developed approach will be utilized to investigate the effective (hydraulic and hydromechanical) properties of the rock mass taken from the Sellafield site [28, 29]. Finally, the paper will be finished with some conclusions.

2. Fracture Modeling by the Embedded Fracture Continuum (EFC) Approach

Shearing the same idea with some other approaches proposed in the literature such as the “fracture continuum approach” [9, 12–17], the “fracture-cell model” [18], or the “embedded discrete fracture model” [19–22], the embedded fracture continuum (EFC) approach is developed in this work to model explicitly the fracture network on the coupled hydromechanical behavior of fractured rocks. The common point of this approach with abovementioned models is the use of the fracture cell notion as the fundamental basic brick. Follow that, each fracture is modeled explicitly in the discretized porous medium by a group of grid cells (fracture cells) intersected by this fracture (Figure 1) [30]. Thus, we distinguish a fracture cell which represents the grid cell of the discretized media intersected at least by one fracture [9, 12–18] with the grid cells that are not intersected by any fracture (called the matrix cell). If the matrix cells represent the continuum porous matrix owing the physical properties of the intact rock, the fracture cell properties represent those of an equivalent continuum porous medium counting for hydraulic and mechanical properties of intact rock and ones of fractures intersecting the cell. Note also that, in the EFC approach, we distinguish two types of fracture cell: type I represents the cell intersected by only one fracture, while type II gathers all fracture cells intersected by multiple fractures (more than one fracture), as indicated in Figure 2.

(a)

(b)

Based on the fracture-cell concept, the meshing procedure of fractured media in the EFC approach will be an easy task because the mesh is not necessary to be conformed to fracture geometry, and fracture could intersect whatever element and hence a simple uniform Cartesian mesh can be used. At this point, the present approach shares the same advantage as the well-known extended finite element method (XFEM) largely developed in the last decade [31–33]. However, if this latter approach uses the lower dimensional interface elements to model the fractures, in our proposed approach, they are described by the equidimensional elements.

As highlighted in some previous works [9, 12–18, 21, 23, 24], it is required that the fracture cell size could be quite small in order to satisfy the computational error tolerance and to bring the geometrical description of fractures by a mesh to an acceptable level. This requirement could lead to a very heavy numerical calculation if the uniform Cartesian mesh (called hereafter as the global refinement to distinguish to the local refinement) is utilized in the whole studied domain. In order to satisfy the exigencies of the tolerated numerical error and minimize the number of degree of freedom and hence the computational time, the local mesh refining strategies may be one of the most appropriated choice. In this work, the local refinement technique based on the hanging node method [27, 28, 34–36] is used. The efficiency of this technique in combination with the EFC approach will be demonstrated in the next section through some verification tests.

2.1. Equivalent Poroelastic Properties of Fracture Cell

2.1.1. Equivalent Poroelastic Based on Micromechanical Solution

As previously mentioned, each fracture cell in the model represents an equivalent continuum porous medium that their hydraulic and hydromechanical properties will be calculated from the properties of the intact rock and ones of the embedded fracture in the cell. In this subsection, the determination procedure of these equivalent properties will be presented.

In the literature, some scholars attempt to determine the equivalent properties of fracture cell by using the micromechanical solution. As for instance, in [10], the authors divided the fractured domain into different subregions, and they used Oda’s crack tensor to calculate the upscaled properties for each element representing each subregion. By comparing the results with the alternative discrete fracture network (DFN) models, a reasonably good agreement was observed by these authors. However, the authors pointed out that some difficulties can be encountered in this approach based on Oda’s crack tensor where the size of element becomes small or/and only one or a few fractures intersect the element. Recently, Figueiredo and his collaborators [23, 24] propose to determine Young’s modulus of the elements intersected by a fracture from the following equation:where E_m, k_n, and h indicate, respectively, Young’s modulus of the intact matrix, the normal stiffness of fracture, and the element size. Note that the same formula was used in [37] to determine Young’s modulus of the elements containing fracture that are defined as weak elements.

In effect, this last formula is extracted from the solution of the joint model as detailed below. Moreover, it is important to note that, in their works, Figueiredo and coauthors [23, 24] assumed that the equivalent properties are isotropic and homogeneous for all fracture cells. However, the determination of other equivalent poroelastic properties like the Biot coefficient and Biot’s modulus was not discussed in the abovementioned references. To validate the proposed solution, Figueiredo and coauthors [23, 24] compared their numerical results with the analytical solution in a simple case of one vertical and completely open fracture. A relative difference smaller than 5% was reported by these authors in terms of the stress concentration result at the tip of fracture.

Despites that the use of a micromechanical approach to determine properties of the fracture cell could not be justifiable in all cases, various authors state its efficiency [10, 23, 24]. The main concern for the use of such approach in the effective properties of the fracture cell is the size of the fracture cell often smaller than that of the representative elementary volume (REV). While the micromechanical approach based on Oda’s crack tensor can be a good choice for the medium with high density of fractures embedded in a quite large fracture cell, the joint model seems to be more appropriate for small cell with only one or few fractures intersecting the element. This observation is coherent with the discussion of [38] who studied the poroelasticbehavior of jointed rock. Following these authors, the homogenization techniques based on the well-known Eshelby’s theory [39] can be used for short joints as a limiting situation of embedded ellipsoidal inclusions (known as the cracked medium) while for the long joints cross-cutting the REV, the joint model is suggested.

In this work, our goal is to model the rock masses weakened by a quite long fracture network, and thus the joint model could be appropriate to determine the equivalent properties of fracture cell. The detail about this type of model was largely discussed in the literature ([38] and different references cited therein), and hence in what it follows, only some important results of this model will be captured. Furthermore, in our study, we assume that all the joints have the same Biot coefficient (which is equal to 1), while the Biot modulus of fracture is very large and hence, the overall poroelastic properties of the jointed medium can be written as follows [38]:where represent, respectively, the compliance tensor of the overall medium, intact matrix, and i^th joint family. and designate the normal and shear stiffness, while n_i, t_i, and indicate the unit normal and tangential vectors, and the parameter d_i is the spacing of the i^th joint set.

Consider now the simplest case of a fracture cell of type I that is intersected by one horizontal fracture (the normal vector of fracture coincides with the vertical direction n = e₃). For an isotropic intact matrix characterized by elastic properties (E_m, ν_m), one can deduce straightforwardly from Equation (2) the equivalent transversely isotropic poroelastic properties of the fracture cell:where represent the five elastic modulus of the fracture cell; are, respectively, the Biot coefficients and Biot modulus of the fracture cell, while the parameter G_m = E_m/(2(1 + ν_m)) indicates the shear modulus of intact isotropic rock.

For the general case of a fracture cell of type I intersected by an arbitrary orientated fracture (Figure 3(a)), the equivalent properties calculated in the global coordinates system (e₁, e₂, and e₃) can be determined from ones calculated in the local coordinate system associated with the fracture (n, t, and ) through a rotation with an angle as shown in Figure 3 around the axis e₁. In the local reference (n, t, and ), the equivalent poroelastic properties of the fracture cell are merely those previously obtained in Equation (3) with n coinciding with e₃ and d in that case is calculated as follows:

For that case, the general form of the equivalent poroelastic properties (compliance tensor and Biot tensor) of the fracture cell in the global reference (e₁, e₂, and e₃) after rotation becomes as follows:where

For the fracture cell intersected by multiple inclined fractures (fracture cell type II), its equivalent poroelastic properties by using the solution of the joint model (Equation (2)) have the same matrix form as presented in Equation (5) and their components depend also on the elastic modulus (E_m and ν_m) of the intact matrix, the mechanical properties (), and the corresponding inclination angle θ_i of each fracture i that intersects the cell. The expression of the diagonal components S₁₁, S₂₂, and S₃₃ written in this case are as follows:

The other nonzero components of the compliance tensor and Biot tensor are not presented here just to simplify the presentation.

2.1.2. Isotropic Approximation of the Equivalent Poroelastic Properties of Fracture Cell

As previously mentioned, in their work, Figueiredo and coauthors [23, 24] considered that the elastic properties of all fracture cells generated in the whole medium are isotropic and are characterized by the same Young’s modulus of the elements calculated from the fracture cell intersected by a horizontal fracture. Otherwise, Poisson’s ratio is equal to one of the intact rock:where is defined in Equation (3).

More precisely, in these last works, one does not distinguish fracture cells intersected by one or many fractures. It means also that this approximation is only applicable for fractures owing the same mechanical properties.

In this work, another isotropic approximation of the equivalent poroelastic properties of the fracture cell is proposed. Following this, the isotropic compliance tensor of this proposed approximation has the in-plane compliances coefficients equal to the weakest ones calculated from the joint model. Thus, in difference with the previous approximation (Equation (8)), the isotropic approximation introduced here allows to model the fractures having different mechanical properties by accounting for the contribution of each fracture in the isotropic properties of the fracture cell. More precisely, for the fracture cell of type I intersected by an inclined fracture, their approximated isotropic moduli are proposed as follows:

Similarly, for the fracture cells of type II, intersected by more than one fracture, the effective properties are obtained using the hypothesis of noninteraction leading to the summation of all fracture contributions:

To distinguish with the isotropic approximation proposed in this work, we call hereafter the approximation proposed in Equation (8) as Figueiredo isotropic approximation.

The corresponding isotropic Biot coefficient and Biot modulus of the fracture cell can be deduced as function of the isotropic elastic modulus of the intact matrix and ones of fracture cell (Equation (2)) as follows:

Thanks for its simplicity, the isotropic approximation of the equivalent poroelastic properties of the fracture cell will be used in our proposed EFC approach. The accuracy of the two approximations as detailed here will be investigated in the next section of the paper.

2.2. Equivalent Permeability of Fracture Cell

The determination of the equivalent permeability of the fractured rock masses has been discussed in a lot of work since several decades [9, 15, 16, 18, 40, 41]. It is largely accepted that, in fractured rocks where the permeability of the intact matrix is negligible as compared to the permeability of fracture, the flow is principally conducted through the fracture inside the domain. The permeability of the fracture K_f is calculated from the fracture’s hydraulic aperture that, in the 2D case, assuming smooth fractures coincides with the width of the fracture. Combined with laminar flow hypothesis, the permeability K_f of the fracture is finally obtained as follows [42, 43]:

The effective permeability of the fracture cell of type I is obtained by considering the general case of an arbitrary fracture orientated to an angle θ to the direction of pressure gradient as shown in Figure 3(b). In Figure 3(b), two constant heads p₁ and p₂ are applied, respectively, on the left and on the right of a fracture cell, while the top and bottom sides of the cell are considered to be impermeable. Considering that the flow in the fracture cell is exclusively through the fracture, permeability of the matrix is small with respect to the permeability of fracture, the equivalent permeability of the fracture cell in the x direction [40, 41]:

Similar result is obtained for the permeability of the fracture cell in the vertical direction as follows:

These equations show that the effective permeability of the fracture cell equals that of the fracture intersecting the cell, corrected for its orientation and multiplied by the ratio between the cell size and hydraulic aperture. Clearly, the orientation of the anisotropic permeability of the fracture cell depends on the orientation of the fracture with the highest value in the direction of the plane of the fracture. The permeability tensor in the global coordinate could be obtained again by calculating firstly the permeability in the local fracture coordinate system and then making a rotation to bring results on the global coordinate system.

Similar to the mechanical case, for simplicity, we propose to use for a fracture cell (type I) an isotropic hydraulic permeability that is defined as the maximum between K_x and K_y:

Correspondingly, for the fracture cell of type II, the effective permeability is obtained by adding the contribution of all cracks intersecting the fracture cell:

The isotropic approximation of the permeability in this last case is obtained by generalizing Equation (15):

3. Verification Tests

In this section, through some numerical verification tests, the accuracy of the developed EFC approach will be verified. In the first stage, two examples which represent correspondingly the purely mechanical and purely hydraulic problems are considered with the focus on the convergence solution of EFC with respect to the fracture-cell size and on the correctness of the isotropic approximations of the equivalents (hydraulic and hydromechanical) properties of the fracture cell. For the first purpose, the size of fracture cell will be gradually decreased by using both the global refinement with uniform meshes and local refinement technique utilizing the hanging nodes. The efficiency of this local refinement technique in combination with the EFC approach will be highlighted. Finally, the same methodology will also be applied in the third example to verify the validity of the proposed EFC approach in the coupled hydromechanical modeling of fractured rocks.

3.1. Purely Mechanical Test

As the first test case, we consider an elastic medium with dimension L = 1.0 m of width and H = 1.25 m in height which is intersected by a fracture with 1 mm of width as shown in Figure 4(a). The fracture is inclined with respect to the horizontal axis by an angle ϕ. Due to the symmetry of the considered problem, we consider only the inclination angle ϕ ranging from 0° to 45° (with interval 5° for each case of study). The elastic properties of the isotropic elastic medium consist of elements of Young’s modulus E_m = 84.6 (GPa) and Poisson ratio ν_m = 0.24, while the normal and shear stiffness of the fracture are, respectively, k_n = 434 (GPa/m) and k_t = 86.8 (GPa/m). These are typical values published for a granitic fractured rock at the Sellafield site in UK [28, 29]. On the top boundary of the considered model, a normal compressive stress 24.46 MPa is applied, while on the bottom boundary, all the displacements are fixed.

(a)

(b)

Different numerical simulations can be conducted to study this problem. In the first step, the commercial code Flac3D (Itasca Consulting Group Inc.) is used in which the fracture is explicitly modeled as the elastic interface elements (Figure 5(a)). The results obtained from this simulation will serve as the referent results for the comparison purpose of the results calculated from the EFC approach. In addition, for this last approach, two meshing configurations of the fractured medium are considered. The first configuration uses a uniform mesh, resulting from the global refinement technique as shown in Figure 5(b). For this configuration, the commercial code Flac3D is also chosen in which we implemented the EFC approach. This simulation demonstrates the easy use of the EFC approach presented in this work in an existence commercial code when only a simple algorithmic to search the mesh cells intersected by the fracture is required. The other configuration corresponding to the case uses the nonuniform meshing based on the local refinement technique (Figure 5(d)) to describe the fracture in the EFC approach. The numerical simulation of this second configuration is performed using open source library deal.II. The reasons for the choice of this code are not only on its capacity to deal with different problems, as demonstrated in numerous works ([44–47] and different references cited therein) but also on its capacity for local refinement of meshes using the hanging nodes technique [27, 30, 34–36]. Finally, as shown in Figure 5(c), the fracture can also be modeled with its real geometry by using the conformed mesh model which is also simulated in deal.II. The grid cells in this case coincide with the boundaries of the fracture, i.e., the size of elements is exactly equal to the fracture aperture (h = = 1 mm) and mechanical properties of elements are equal to the elastic moduli of the fractures (E_f, ν_f, and G_f). These latter parameters can be determined from the normal and shear stiffness of fracture (k_n and k_t) through the following relationships [48]:

(a)

(b)

(c)

(d)

In Figure 6 are depicted the results of the vertical displacement following a vertical profile at the middle of the domain. These results are calculated with different models of fracture inclined at 45° to the horizontal direction. As the first statement, a superposition of the results determined from the models based on the interface elements and the conformed mesh can be noted as expected which demonstrated the correctness of the numerical simulation by deal.II and justifies Equation (18) deriving the relationships of the mechanical properties of fracture. Our attention is then focused to the results evaluated from the EFC approach, which were conducted on the two other models using the uniform and nonuniform meshes and the isotropic approximations of the fracture cell properties. For the clarity purpose, we detail here that, in the uniform mesh model 128 × 156, elements are used which is similar to the interface element model. The fracture cell size of the uniform mesh model (which is about h = 7.8125 mm) is similar as one chosen in the simulation conducted on the nonuniform mesh model. However, in this last model, thanks to the local refinement technique, the matrix cell size can be taken much more important (about 8 times greater) with respect to the one of the uniform mesh model. Consequently, both models use the same size of the fracture cell, but only 920 elements are necessary in the nonuniform mesh model as compared to 19968 elements for the uniform mesh model. As illustrated in Figure 6, the results achieved from these two models are similar which demonstrated the efficiency of the local refinement based on the hanging node technique in combination with the EFC approach to model the fractured medium.

Note that in all these simulations the isotropic approximations of effective properties of the fracture cell are used. The results exhibit a quite significant difference between the two isotropic approximations of the equivalent mechanical properties of the fracture cell. In comparison with the referent results provided from the interface element model, a maximum difference of about 2.86% is noted for the isotropic approximation proposed in this work, while the difference of up to 12.5% is evaluated for the approximation presented by Figuereido et al. [23, 24]. The gap between the two isotropic approximations can be explained by the fact that, in the Figuereido isotropic approximation, Young’s modulus is assumed constant and independent on the inclination angle of fracture. In addition, as shown in Figure 7, the oriented fracture (represented by an inclination angle θ) can present an important effect on the equivalent elastic properties of the fracture cell, particularly in case of large fracture cell size (with respect to the fracture aperture) and high inclination angle (around 45°). For a chosen fracture cell size h, the elastic Young modulus of the fracture cell is smallest in case of θ = 45° and highest in case of θ = 0°(or θ = 90°). Thus, the isotropic elastic Young’s modulus proposed by Figuereido et al. [23, 24] is overestimated in case of inclined fracture. As consequence, the approximation of these last authors seems appropriated only in cases of quasi-horizontal or quasi-vertical fractures.

(a)

(b)

As indicated by different scholars [9, 12–18, 21, 22], the correctness of the approaches based on the fracture cell concept could depend on the size of this element. This discussion can be verified now when we investigate the influence of this parameter on the convergence of the EFC approach. Note that from the results presented in the previous paragraph, in the last part of this paper, we discuss only on the EFC approach using our isotropic approximation and the local refinement technique just to simplify presentation. In deal.II this could be done easily by changing the local refinement level around the fracture while keeping a constant prechosen matrix cell size. In addition, our investigation will be carried on with different inclination angles of fracture ranging from 0° to 45°.

In Figure 8(a), we present the vertical displacement determined at the controlled point “A” located on the middle of the top boundary. These results are evaluated from different values of fracture cell size and inclination angle of fracture ranging from 0° to 45°. In these simulations, the matrix cell size is fixed at 62.5 mm, while the fracture cell size is gradually decreased (from h = 62.5 mm to h = 1.9531 mm) by changing the local refinement level around the fracture in deal.II. In general, we can state that the vertical displacement at this controlled point decreases by decreasing the size of the fracture cell. This displacement tends to an asymptotic value which approaches to the results calculated with the conformed mesh (case h = ). However, from a fracture cell of size h = 7.81 mm, it can be noted that the relative difference of results is smaller than 1% for the whole considered range of inclination angle of fracture.

(a)

(b)

Figure 8(b) highlights the displacement of the controlled point calculated from the EFC approach in comparison with the referent results based on the elastic interface element model (which is exactly the results conducted with the conformed mesh model). It shows that the EFC approach could underestimate the result for the inclination angle of fracture ranging from 0° to about 30°, while it could overestimate for the case of fracture owning the inclination angle from 30° to 45°. However, for all considered cases of inclined fracture, the relative difference is smaller than 3%. This difference is small enough to validate the EFC approach using the isotropic approximation proposed in this work in modeling the mechanical behavior of the fractured medium.

3.2. Purely Hydraulic Test

In this second test, the same geometric model as the previous example (Figure 4(b)) is considered to simulate the hydraulic behavior in the fractured medium. The validation of the EFC approach is also done by comparing with the referent results provided from the calculation conducted on the conformed mesh model. More precisely, for the referent case, the fracture is modeled explicitly with his real aperture ( = 1 mm), and the permeability determined from Equation (12) is equal to K_f = /12 = 8.33 × 10⁻⁸ (m²). The permeability of the porous intact matrix is taken equal to 2.4 × 10⁻¹⁵ (m²), a typical value for a Sellafield granitic rock (UK) as reported in [49]. Thus, the permeability of fracture is seven times higher than the matrix permeability. The viscosity and compressibility of fluid are 0.001(Pa.s) and 5.0 × 10⁻¹⁰ (Pa⁻¹), respectively. Note here that the permeability of the fracture cell is supposed to be isotropic and calculated from the permeability of fracture following Equation (15). As shown in Figure 4(b), the top and bottom boundaries of the model are considered to be impermeable, while a constant gradient pressure is imposed on two lateral boundaries by applying the pressure p₁ = 10 (MPa) on the left and p₂ = 0 (MPa) on the right of the model.

Similar to the previous mechanical test, in the first stage, we will investigate the convergence of the EFC approach by modeling the fracture with different local refinement levels corresponding to different values of fracture-cell size ranging from h = 62.5 mm to h = 1.95 mm. The results depicted in Figure 9(a) present the horizontal fluid flux across the fractured medium, which tends to decrease with respect to the fracture cell size. However, the obtained results seem to converge at a constant value from a local refinement of mesh corresponding to the fracture cell size about h = 15.625 mm with relative differences of about 1% for all considered inclination angles of fracture.

(a)

(b)

In the second stage, the results calculated with the EFC approach (with the fracture cell size h = 7.8125 mm) are compared with the referent results conducted with conformed mesh configuration as illustrated in Figure 9(b). A quite similar tendency as in the purely mechanical problem can be observed. An overestimation of the fluid flux for the inclination angle between 0° and 40° is noted, while for the inclined fracture of about 45°, it seems that the EFC approach proposed here may slightly underestimate the result. In general, a maximum difference of about 5% (case of fracture with inclination angle 25°) between the two models is observed which seems acceptable to validate the capability of the EFC approach to solve the hydraulic problem of fractured media.

3.3. Coupled Hydromechanical (HM) Test

As the last verification test, in this part, we use the EFC approach to model the coupled hydromechanical behavior of the fractured medium. To this end, we retain the outer geometry (L = 1.0 m and H = 1.25 m) and the mechanical properties and hydraulic properties of the porous matrix as the two previous examples. The two other necessary coupling parameters of the porous matrix, namely, the Biot coefficient and Biot modulus are taken, respectively, equal to b = 0.8 and M = 20.8 (GPa). Two configurations are considered: the first one consists of a short fracture (l = 0.6 m) placed in the center of the domain and orients at 45° (Figure 10(a)) in respect with the horizontal direction, and the second case accounts for two fractures (with the same length l = 0.6 m) with two (arbitrary) inclination angles of 25° and 90° (vertical fracture), respectively, as illustrated in Figure 10(b). All the parameters of the fractures like the aperture, mechanical, and hydraulic properties are taken as in the previous examples, i.e., = 1 mm, k_n = 434 (GPa/m), k_t = 86.8 (GPa/m), and K_f = 8.33 × 10⁻⁸ (m²). The boundary conditions for this test are also highlighted in Figure 11 where the left, right, and bottom boundaries of the saturated medium are fixed and impermeable, while on the top boundary, it is permeable, and a compressive normal stress σ₁ = 24.46 (MPa) is applied. Again, the validation of the EFC approach will be conducted by comparing with the referent results calculated from the conformed mesh model (Figure 11) in which the real geometry and poroelastic properties of fractures are used. For this purpose, two controlled points “A” and “B” located, respectively, at (0.5 m and 1.25 m) and (1.0 m and 0.0 m) are selected to observe the results during the transient calculations (Figure 10). By using the same size, the hydraulic and mechanical properties of the fracture cells (case of the nonuniform mesh model) are similar as ones of the previous tests while their coupling HM parameters are determined from Equation (11).

(a)

(b)

(a)

(b)

The comparison of the results determined from the two models is presented in Figure 12, where the vertical displacement at the controlled point A and excess pore pressure at the controlled point B are traced against the elapsed time. As expected, during the consolidation, the decrease of the excess pore pressure induces an increase of effective stress which explains the increase of the vertical displacement in the medium versus elapsed time. A very good agreement of the results was calculated from two models (nonuniform mesh and conformed mesh) for both configurations of embedded fractures in the porous medium. The relative error seems more significant in the case of two fractures which is about 4.5% for the displacement and 1.5% for the result of pore pressure. Again, these moderate differences allow confirming the validity of our proposed EFC approach to model the HM behavior of the fractured media.

(a)

(b)

(c)

4. Application of the EFC Approach to Investigate the Effective Permeability and Elastic Compliance Tensor of Fractured Rock Masses in Real Field

The capacity of the EFC approach will be also demonstrated in this section through a numerical investigation of the effective permeability and elastic compliance tensor of the real rock mass taken from the Sellafield site. The abundancy of information about this site, thanks to the investigation program performed by United Kingdom Nirex Limited [28, 29] as part of the second Bench Mark Test (BMT2) of the international collaborative program DECOVALEX III, is the main reason of this choice. Additionally, the information of the site is also well-synthetized in various works [4–7, 49–51]. The necessary information about the spatial distribution of fractures is summarized in Table 1 and is sufficient to generate without difficulty the fracture network in the rock mass. These parameters consist of input data about the fractures network such as the distribution of fracture’s length, orientation, location, and fractures’ aperture. Following this, the fracture trace lengths can be characterized by the power-law distribution [4–7, 49]:where N_F is the number of fractures (with fracture length greater than the length L) per unit area, C is the constant density, and D is the fractal dimension.

Concerning the orientations of fractures, the Fisher distribution [4–7, 28, 29, 50] is largely adopted. It was also pointed out that, in this site, there are four principal sets of fractures (Table 1) and the probability of a fracture of angle deviating from the mean orientation angle of a fracture set is given by Baghbanan and Jing [51]:where Fisher’s constant K is assigned for each fracture set [28, 29].

Besides, in our investigation, the fractures with uniform aperture = 65 μm as indicated in the works of Min and his collaborators [4, 5] will be chosen, while the spatial location of fractures’ center is represented by the Poisson distribution as largely accepted in the literature [4–7, 51].

In the upscaling procedure to determine the effective properties of the fractured rock masses, the investigation of REV is an important step. The REV size is defined as the minimum size beyond which the upscaled values are stable [52]. For the fractured rock masses in the Sellafield site, Min and his collaborators [4, 5] by using the UDEC code demonstrated the existence of the REV whose size can be chosen from 2 m to 6 m with the coefficient of variation varying from 10% to 5%, respectively [4], for the purely mechanical problem. A quite similar REV size was also proposed for the hydraulic problem [5], which can range from 5 m to 8 m with the corresponding 20% and 10% acceptable variations, respectively. Based on these discussions, in this work, the REV size of 5 m is chosen to calculate the effective permeability and mechanical properties of this fractured rock masses using the EFC approach.

4.1. Methodology for the Determination of the Effective Elastic Properties

The effective elastic properties of fractured rock masses in the Sellafield site will be calculated in this work based on the 2D plane strain hypothesis. The same methodology to determine the effective compliance matrix in Equation (21) as proposed in [4, 53] will be adopted here. For plane strain hypothesis, the linear elastic behavior of a general anisotropic media written in the matrix form is as follows [4, 53]:

Note here that the coordinate Oxyz presents the global coordinate (z = 1, x = 2, y = 3) as shown in Figure 3(a).

Assuming that the Young modulus and Poisson ratios in the z direction are equal to the ones of the intact matrix (), it follows . The study in [4] also revealed that the values of are very small in comparison with the other components, while are equal to zero due to the fact that shear stress does not affect the z direction deformation. Thus, the compliance tensor in Eq. (21) can be rewritten in an explicit form as follows:

Hence, calculating the effective elastic properties of the fractured rock masses reduces to the determination of only five components which can be conducted through three numerical experiments (independent loads), as illustrated in Figure 13 (also see the details in [52, 53]).

(a)

(b)

(c)

The first loading set (Figure 13(a)) represents a uniaxial compressive test with the normal stress applied in the x direction. Two components and can be calculated from this test taking into account the condition of (hence, ):

Similarly, the two other components and can be calculated from the second loading set (Figure 13(b)):

The last component can be determined from the third loading set (Figure 13(c)):

Note that , and designate the averaged normal and lateral strain under uniaxial stress in x and y directions correspondingly, while indicates the average shear strain under pure stress loading.

From these components of the effective compliance tensor, we can deduce without difficulty the in-plane effective Young modulus, Poisson ratio, and shear modulus of the equivalent medium:

4.2. Methodology for the Determination of the Effective Permeability

The evaluation of the effective permeability of the fractured rock masses in each direction x or y is calculated with a constant hydraulic pressure gradient imposed on the boundaries of the corresponding direction (x or y) and no-flow boundaries of the other direction (Figure 14). With the assumption that the permeability tensor is diagonal in the global coordinates system and by adopting Darcy’s law for the upscaled fractured media [5], their effective permeability and following the corresponding the x and y directions are determined from the following equation:where and are the steady flow rates at the boundaries of the corresponding x and y directions, A is the cross-sectional area, and μ is the dynamic viscosity of fluid.

(a)

(b)

In the numerical simulations, all terms on the right-hand side of Equation (25), except for and , are specified, while the steady flow rate and are calculated numerically. Thus, and can be back calculated.

4.3. Monte Carlo Simulation of the Effective Hydraulic and Mechanical Properties of Sellafield Fractured Rock Masses

To obtain the representative results of the effective hydraulic and mechanical properties, the Monte Carlo simulations were carried out to obtain the effective mechanical and hydraulic properties. For that, ten realizations of the REV of Sellafield fractured rock mass were generated using the synthetized parameters in Table 1. As in the previous verification tests, an investigation of the fracture cell size on the results of the effective properties was also carried out. However, to not complicate the presentation, this study will not be detailed here. The results of this investigation indicate that a fracture cell size of about h = 9.77 mm can be used. This fracture cell size presents an acceptable convergence of the effective properties of fractures rock masses when a maximum relative difference of about 6% can be stated in comparison with the case of the fracture cell having two times smaller size (h = 4.89 mm).

It is worth to note that, in the literature, the effective permeability and elastic compliance tensor of fractured rock masses are usually numerically determined by utilizing the universal distinct element code (UDEC) for only sample in which the dead ends and isolated fractures were eliminated [4–7, 54–57]. In this work, the EFC approach as proposed in work allows to simulate both two types of samples with and without dead-end fractures. The ten realizations of the REV of Sellafield fractured rock masses accounting for the dead-end and without dead-end fractures are depicted in Figure 15.

In Table 2, the mean value, standard deviation, and coefficient of variation of the effective permeability and elastic properties of the Sellafield fractured rock masses calculated with these ten realizations of REV are summarized. The obtained results show that this fractured rock mass can be considered as isotropic. More precisely, the mean values of Young’s moduli and permeability following the two directions x and y calculated with the dead-end fractures model are correspondingly equal to , meaning that the anisotropic mechanical and hydraulic ratios are and . Averaging the results of two directions we have, Young’s modulus, Poisson’s ratio, and permeability of the fractured rock masses are . The corresponding results of this rock calculated in the model without dead-end fractures are . These last results seem consistent with the ones presented in [4, 5, 51]. More precisely, by using the distinct element method (UDEC), these authors showed that the normalized elastic moduli in x and y directions are about 31% and 34% than that in intact rock meaning that the averaged Young’s modulus of this fractured rock masses is about . Concerning Poisson’s ratio, a mean value presented by these authors seems higher than our result. With a 95% confidence level, the permeability in x and y directions determined by [5] are, respectively, which present an averaged value of

A comparison of the mean results calculated from the dead-end and without dead-end fractures models highlighted a quite significant difference of about 17% of the effective Young and shear moduli, while for the effective permeability, a more moderate difference of about 4% is observed. This difference confirmed the discussion of Yang et al. [53] who concluded that the results calculated from the model without dead-end fractures overestimate the deformation modulus of the fractured rock masses.

5. Conclusions

We introduce in this work a new numerical simulation method to model the coupled hydromechanical behavior of fractured porous media. This method called the embedded fracture continuum (EFC) approach allows to model explicitly the embedded fractures network in the continuum porous media by using the fracture-cell concept. Each fracture cell in this approach presents an equivalent continuum porous medium that has its own properties calculated from the contributed properties of intact matrix and fractures. The details of the determination procedure of equivalent poroelastic properties of the fracture cell are derived in this paper using the well-known joint model. Through different verification tests, the capability of the EFC approach is demonstrated, while the isotropic approximation of both the hydraulic and poroelastic properties of the fracture cell seems an appropriate choice. This EFC approach which can be easily implemented in a commercial or open-source code introduces an efficient numerical tool to model the hydromechanical behavior of the fractured media. To overcome the required small size of fracture cell in the EFC approach, the local refinement technique based on the hanging nodes is considered and presents its efficiency. As the first application, in this paper, we use this approach to investigate the effective permeability and elastic properties of the Sellafield fractured rock mass. The Monte Carlo simulation conducted with different generations of the fractures network in this rock presents a good consistency with the available results in the literature calculated from the distinct element method.

Data Availability

The numerical data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The authors declare that there are no conflicts of interest regarding the publication of this paper.

References

G. I. Barenblatt and I. P. Zheltov, “Fundamental equations of filtration of homogeneous liquids in fissured rocks,” Soviet physics, Doklady, vol. 13, 1960.
View at: Google Scholar
J. E. Warren and P. Root, “The behavior of naturally fractured reservoirs,” Transactions of the Society of Mining Engineers of AIME, vol. 228, 1963.
View at: Google Scholar
M. Sahimi, Flow and Transport in Porous Media and Fractured Rocks, Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim, Germany, 2011.
K. B. Min and L. Jing, “Numerical determination of the equivalent elastic compliance tensor for fractured rock masses using the distinct element method,” International Journal of Rock Mechanics and Mining Sciences, vol. 40, no. 6, pp. 795–816, 2003.
View at: Publisher Site | Google Scholar
K. B. Min, L. Jing, and O. Stephansson, “Determining the equivalent permeability tensor for fractured rock masses using a stochastic REV approach: method and application to the field data from Sellafield, UK,” Hydrogeology Journal, vol. 12, no. 5, pp. 497–510, 2004.
View at: Publisher Site | Google Scholar
K. B. Min, J. Rutqvist, C. F. Tsang, and L. Jing, “Stress-dependent permeability of fractured rock masses: a numerical study,” International Journal of Rock Mechanics and Mining Sciences, vol. 41, no. 7, pp. 1191–1210, 2004.
View at: Publisher Site | Google Scholar
K. B. Min, J. Rutqvist, C. F. Tsang, and L. Jing, “Thermally induced mechanical and permeability changes around a nuclear waste repository—a far-field study based on equivalent properties determined by a discrete approach,” International Journal of Rock Mechanics and Mining Sciences, vol. 42, no. 5-6, 2005.
View at: Publisher Site | Google Scholar
J. C. Jaeger, N. G. W. Cook, and R. W. Zimmerman, Fundamentals of Rock Mechanics, Blackwell Publishing, Hoboken, New Jersey, USA, 4th edition, 2007.
F. E. Botros, A. E. Hassan, D. M. Reeves, and G. Pohll, “On mapping fracture networks onto continuum,” Water Resources Research, vol. 44, no. 8, 2008.
View at: Publisher Site | Google Scholar
J. Rutqvist, C. Leung, A. Hoch, Y. Wang, and Z. Wang, “Linked multicontinuum and crack tensor approach for modeling of coupled geomechanics, fluid flow and transport in fractured rock,” Journal of Rock Mechanics and Geotechnical Engineering, vol. 5, no. 1, pp. 18–31, 2013.
View at: Publisher Site | Google Scholar
J. Jiang and R. M. Younis, “Numerical study of complex fracture geometries for unconventional gas reservoirs using a discrete fracture-matrix model,” Journal of Natural Gas Science and Engineering, vol. 26, pp. 1174–1186, 2015.
View at: Publisher Site | Google Scholar
U. Svensson, “A continuum representation of fracture networks. Part I: method and basic test cases,” Journal of Hydrology, vol. 250, no. 1–4, 2001.
View at: Publisher Site | Google Scholar
U. Svensson, “A continuum representation of fracture networks. Part II: application to the Äspö Hard Rock Laboratory,” Journal of Hydrology, vol. 250, no. 1–4, pp. 187–205, 2001.
View at: Publisher Site | Google Scholar
S. A. McKenna and P. C. Reeves, “Fractured continuum approach to stochastic permeability modeling,” Computer Applications in Geology, vol. 5, 2006.
View at: Google Scholar
D. M. Reeves, D. A. Benson, and M. M. Meerschaert, “Transport of conservative solutes in simulated fracture networks: 1. Synthetic data generation,” Water Resources Research, vol. 44, no. 5, 2008.
View at: Publisher Site | Google Scholar
D. M. Reeves, D. A. Benson, M. M. Meerschaert, and H. P. Scheffler, “Transport of conservative solutes in simulated fracture networks: 2. Ensemble solute transport and the correspondence to operator-stable limit distributions,” Water Resources Research, vol. 44, no. 5, 2008.
View at: Publisher Site | Google Scholar
Y. Hao, P. Fu, S. M. Johnson, and C. R. Carrigan, “Numerical studies of coupled flow and heat transfer processes in hydraulically fractured geothermal reservoirs,” GRC Applications, vol. 36, 2012.
View at: Google Scholar
A. Sakhaee-Pour and M. F. Wheeler, “Effective flow properties for cells containing fractures of arbitrary geometry,” SPE Journal, vol. 21, no. 3, pp. 965–980, 2016.
View at: Publisher Site | Google Scholar
C.-C. Lee, C.-H. Lee, H.-F. Yeh, and H. I. Lin, “Modeling spatial fracture intensity as a control on flow in fractured rock,” Environmental Earth Sciences, vol. 63, no. 6, pp. 1199–1211, 2011.
View at: Publisher Site | Google Scholar
L. Li and S. H. Lee, “Efficient field-scale simulation of black oil in a naturally fractured reservoir through discrete fracture networks and homogenized media,” SPE Reservoir Evaluation & Engineering, vol. 11, no. 4, pp. 750–758, 2008.
View at: Publisher Site | Google Scholar
A. Moinfar, K. Sepehrnoori, R. T. Johns, and A. Varavei, “Coupled geomechanics and flow simulation for an embedded discrete fracture model,” in Proceedings of SPE Reservoir Evaluation & Engineering, no. 1998, The Woodlands, TX, USA, 2013.
View at: Google Scholar
X. Yan, Z. Huang, J. Yao, Y. Li, and D. Fan, “An efficient embedded discrete fracture model based on mimetic finite difference method,” Journal of Petroleum Science and Engineering, vol. 145, pp. 11–21, 2016.
View at: Publisher Site | Google Scholar
B. Figueiredo, C. F. Tsang, J. Rutqvist, and A. Niemi, “A study of changes in deep fractured rock permeability due to coupled hydro-mechanical effects,” International Journal of Rock Mechanics and Mining Sciences, vol. 79, pp. 70–85, 2015.
View at: Publisher Site | Google Scholar
B. Figueiredo, C. F. Tsang, J. Rutqvist, and A. Niemi, “Study of hydraulic fracturing processes in shale formations with complex geological settings,” Journal of Petroleum Science and Engineering, vol. 152, pp. 361–374, 2017.
View at: Publisher Site | Google Scholar
D. P. Do and D. Hoxha, “Effective transfer properties of partially saturated geomaterials with interfaces using the immersed interface method,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 37, no. 18, pp. 3237–3257, 2013.
View at: Publisher Site | Google Scholar
H. Hoang, D. Hoxha, N. Belayachi, and D.-P. Do, “Modelling of two-phase flow in capillary porous medium by a microscopic discrete approach,” European Journal of Environmental and Civil Engineering, vol. 17, no. 6, pp. 444–452, 2013.
View at: Publisher Site | Google Scholar
W. Bangerth, R. Hartmann, and G. Kanschat, “Deal. II—a general purpose object oriented finite element library,” ACM Transactions on Mathematical Software, vol. 33, no. 4–25, 2007.
View at: Publisher Site | Google Scholar
Nirex, “Evaluation of heterogeneity and scaling of fractures in the borrowdale volcanic group in the Sellafield area,” Nirex, Harwell, UK, 1997, Nirex Rep. SA/97/028.
View at: Google Scholar
Nirex, “Data summary sheets in support of gross geotechnical predictions,” Nirex, Harwell, UK, 1997, Nirex Rep. SA/97/052.
View at: Google Scholar
H. L. Dang, “A hydro-mechanical modeling of double porosity and double permeability fractured reservoirs,” University of Orleans, Orléans, France, 2018, Ph.D. thesis.
View at: Google Scholar
N. Watanabe, W. Wang, C. I. McDermott, T. Taniguchi, and O. Kolditz, “Uncertainty analysis of thermo-hydro-mechanical coupled processes in heterogeneous porous media,” Computational Mechanics, vol. 45, no. 4, pp. 263–280, 2010.
View at: Publisher Site | Google Scholar
N. Watanabe, W. Wang, J. Taron, J. U. Görke, and O. Kolditz, “Lower-dimensional interface elements with local enrichment: application to coupled hydro-mechanical problems in discretely fractured porous media,” International Journal for Numerical Methods in Engineering, vol. 90, 2012.
View at: Google Scholar
A. R. Mohammadnejad and T. Khoei, “Hydro-mechanical modeling of cohesive crack propagation in multiphase porous media using the extended finite element method,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 37, no. 10, pp. 1247–1279, 2013.
View at: Publisher Site | Google Scholar
W. C. Rheinboldt and C. K. Mesztenyi, “On a data structure for adaptive finite element mesh refinements,” ACM Transactions on Mathematical Software, vol. 6, no. 2, pp. 166–187, 1980.
View at: Publisher Site | Google Scholar
L. Demkowicz, J. T. Oden, W. Rachowicz, and O. Hardy, “Toward a universal h-p adaptive finite element strategy, part 1. Constrained approximation and data structure,” Computer Methods in Applied Mechanics and Engineering, vol. 77, no. 1-2, pp. 79–112, 1989.
View at: Publisher Site | Google Scholar
F. Bornemann, B. Erdmann, and R. Kornhuber, “Adaptive multilevel methods in three space dimensions,” International Journal for Numerical Methods in Engineering, vol. 36, no. 18, pp. 3187–3203, 1993.
View at: Publisher Site | Google Scholar
P. Z. Pan, X. T. Feng, X. H. Huang, Q. Cui, and H. Zhou, “Coupled THM processes in EDZ of crystalline rocks using an elasto-plastic cellular automaton,” Environmental Geology, vol. 57, 1009.
View at: Google Scholar
S. Maghous, L. Dormieux, D. Kondo, and J. F. Shao, “Micromechanics approach to poroelasticbehavior of a jointed rock,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 37, no. 2, pp. 111–129, 2011.
View at: Publisher Site | Google Scholar
J. D. Eshelby, “The determination of the elastic field of an ellipsoidal inclusion, and related problems,” in Proceedings of the Royal Society A, vol. 241, no. 1226, 1957.
View at: Google Scholar
R. W. Parsons, “Permeability of idealized fractured rock,” Society of Petroleum Engineers Journal, vol. 6, no. 2, pp. 126–136, 1966.
View at: Publisher Site | Google Scholar
R. A. Nelson, Geologic Analysis of Naturally Fractured Reservoirs, Gulf Professional Publishing, Houston, TX, USA, 2nd edition, 2001.
R. L. Kranz, A. D. Frankel, T. Engelder, and C. H. Scholz, “The permeability of whole and jointed Barre Granite,” International Journal of Rock Mechanics and Mining Sciences & Geomechanics Abstracts, vol. 16, no. 4, pp. 225–234, 1979.
View at: Publisher Site | Google Scholar
P. A. Witherspoon, J. S. Y. Wang, K. Iwai, J. E. Gale, and J. E., “Validity of Cubic Law for fluid flow in a deformable rock fracture,” Water Resources Research, vol. 16, no. 6, pp. 1016–1024, 1980.
View at: Publisher Site | Google Scholar
T. Wick, G. Singh, and M. F. Wheeler, “Fluid-filled fracture propagation with a phase-field approach and coupling to a reservoir simulator,” SPE Journal, vol. 21, no. 3, pp. 981–999, 2016.
View at: Publisher Site | Google Scholar
T. Heister, M. F. Wheeler, and T. Wick, “A primal-dual active set method and predictor-corrector mesh adaptivity for computing fracture propagation using a phase-field approach,” Computer Methods in Applied Mechanics and Engineering, vol. 290, pp. 466–495, 2015.
View at: Publisher Site | Google Scholar
A. Mikelić, M. F. Wheeler, and T. Wick, “Phase-field modeling of a fluid-driven fracture in a poroelastic medium,” Computers & Geosciences, vol. 19, no. 6, 2015.
View at: Google Scholar
S. Lee, M. F. Wheeler, and T. Wick, “Pressure and fluid-driven fracture propagation in porous media using an adaptive finite element phase field model,” Computer Methods in Applied Mechanics and Engineering, vol. 305, pp. 111–132, 2016.
View at: Publisher Site | Google Scholar
M. FatehiMarji, A. Pashapour, and J. Gholamnejad, “Relationship between fracture dip angle, aperture and fluid flow in the fractured rock masses,” Journal of Mining and Environment, vol. 2, 2011.
View at: Google Scholar
S. Namdari, A. Baghbanan, and M. J. Habibi, “Effects of matrix permeability and fracture density on flow pattern in dual porous rock masses,” Rock Mechanics and Rock Engineering, vol. 1, 2016.
View at: Google Scholar
P. Blum, R. Mackay, M. S. Riley, and J. L. Knight, “Performance assessment of a nuclear waste repository: upscaling coupled hydro-mechanical properties for far-field transport analysis,” International Journal of Rock Mechanics and Mining Sciences, vol. 42, no. 5‐6, pp. 781–792, 2005.
View at: Publisher Site | Google Scholar
BaghbananA and L. Jing, “Stress effects on permeability in a fractured rock mass with correlated fracture length and aperture,” International Journal of Rock Mechanics and Mining Sciences, vol. 45, no. 8, 2008.
View at: Publisher Site | Google Scholar
A. Pouya and M. Ghoreychi, “Determination of rock mass strength properties by homogenization,” International Journal for Numerical and Analytical Methods in Geomechanics, vol. 25, 2011.
View at: Google Scholar
J. P. Yang, W. Z. Chen, Y. Dai, and H. D. Yu, “Numerical determination of elastic compliance tensor of fractured rock masses by finite element modeling,” International Journal of Rock Mechanics and Mining Sciences, vol. 70, pp. 474–482, 2014.
View at: Publisher Site | Google Scholar
Y. M. Alshkane, A. M. Marshall, and L. R. Stace, “Prediction of strength and deformability of an interlocked blocky rock mass using UDEC,” Journal of Rock Mechanics and Geotechnical Engineering, vol. 9, no. 3, pp. 531–542, 2017.
View at: Publisher Site | Google Scholar
N. BidgoliM and L. Jing, “Water pressure effects on strength and deformability of fractured rocks under low confining pressures,” Rock Mechanics and Rock Engineering, vol. 48, no. 3, pp. 971–985, 2015.
View at: Google Scholar
N. BidgoliM and L. Jing, “Anisotropy of strength and deformability of fractured rocks,” Journal of Rock Mechanics and Geotechnical Engineering, vol. 6, no. 2, 2014.
View at: Google Scholar
A. Khani, A. Baghbanan, S. Norouzi, and H. Hashemolhosseini, “Numerical investigation of the effect of fracture intensity on deformability and REV of fractured rock masses,” International Journal of Rock Mechanics and Mining Sciences, vol. 60, pp. 345–352, 2013.
View at: Publisher Site | Google Scholar

Copyright

Copyright © 2019 Hong-Lam Dang et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

PDF Download Citation

Download other formats

Order printed copies

Views

910

Downloads

1187

Citations