Comparison of cell- and vertex-centered finite-volume schemes for flow in fractured porous media

https://doi.org/10.1016/j.jcp.2021.110715Get rights and content

Highlights

  • The newly presented schemes lead to smaller linear systems on unstructured grids.

  • Additional degrees of freedom are introduced locally to capture discontinuities.

  • Flux mortars at the interfaces allow for strong incorporation of coupling conditions.

Abstract

Flow in fractured porous media is of high relevance in a variety of geotechnical applications, given the fact that they ubiquitously occur in nature and that they can have a substantial impact on the hydraulic properties of rock. As a response to this, an active field of research has developed over the past decades, focusing on the development of mathematical models and numerical methods to describe and simulate the flow processes in fractured rock. In this work, we present a comparison of different finite volume-based numerical schemes for the simulation of flow in fractured porous media by means of numerical experiments. Two novel vertex-centered approaches are presented and compared to well-established numerical schemes in terms of convergence behavior and their performance on benchmark studies taken from the literature. The new schemes show to produce results that are in good agreement with those of established methods while leading to less degrees of freedom on unstructured simplex grids.

Introduction

Fractures are a common feature in many geological formations throughout the planet, and therefore, they are of particular importance in a variety of geotechnical applications ranging from groundwater management over energy production and storage to the disposal of waste [7]. Although fractures are typically characterized by rather small apertures, they can have significant lateral extents and form interconnected networks that have substantial influence on the overall hydraulic and mechanical behavior of rock [29]. Besides naturally occurring fractures, they might also be introduced into a rock mass as a consequence of human activities, for instance, in geothermal energy or unconventional shale gas production techniques. In the context of radioactive waste disposal, fractures can occur in the vicinity of emplacement tunnels as a consequence of their excavation, and thus, their presence has to be taken into account in safety assessment efforts [38].

As a response to the importance of fractures, a number of modeling approaches and simulation tools have been developed over the past decades. Typically, the approaches are classified into continuum fracture models and discrete fracture matrix (dfm) models. In the first class of models, the fracture geometry is not explicitly captured, but the overall medium is described in terms of effective properties in so-called single-continuum models [44], [46] or by describing the fracture network with one or more individual overlapping continua in dual- or multi-continuum models [55], [33], [43], [58]. In dfm models, the geometries of the most important fractures are resolved, while small scale fractures might still be accounted for by means of effective properties or additional continua. This allows for an explicit representation of fractures whose sizes are comparable to that of the domain and which are difficult to upscale, while maintaining the computational efficiency of continuum approaches where applicable. In such models, the fractures are incorporated either by describing them as heterogeneities in the rock, resolving the interior of the fractures with the computational mesh [40], or as lower-dimensional interfaces embedded in the bulk medium. The latter, mixed-dimensional approach has received significant attention in the past and a variety of models and numerical schemes have been developed, which can be further subdivided into conforming and non-conforming models. In the former, the computational grid used in the bulk medium is constrained such that the grid element faces are aligned with the fractures, while in the latter, bulk medium and fractures can be discretized independently.

For single-phase flow, a conforming, cell-centered approach using two-point flux approximation (tpfa) was presented in Karimi-Fard et al. [31], which has been extended to multi-point flux approximation (mpfa) techniques in Sandve et al. [47], Ahmed et al. [3], [4]. Models based on mixed finite elements have been proposed in Martin et al. [39], Boon et al. [10], while Nordbotten et al. [42] present a general mathematical framework, which includes well-established schemes such as linear Lagrange elements, mixed finite elements, virtual elements and cell-centered finite-volume schemes (tpfa and mpfa). Models using the vertex approximate gradient and the hybrid finite-volume schemes can be found in Brenner et al. [13], while Kadeethum et al. [30] present an approach using the enriched Galerkin method. A benefit of non-conforming approaches is that mesh generation is facilitated, allowing for structured grids to be used in the bulk medium. However, they often struggle to capture the discontinuities that arise, in particular, in the presence of low-permeable fractures. Non-conforming methods on the basis of finite-volume schemes can be found in Karvounis and Jenny [32], Nikitin and Yanbarisov [41], while Köppel et al. [36], Schädle et al. [49] and Schwenck et al. [48], Flemisch et al. [19] present models using finite element and extended finite element methods, respectively. Note that both conforming and non-conforming approaches have been extended to more complex physics in recent works. For instance, conforming methods involving two-phase flow are presented in Brenner et al. [12], Aghili et al. [2], Droniou et al. [16], while Bonaldi et al. [9] additionally include mechanical deformations of the bulk medium. Berge et al. [6] consider single-phase flow, but include the transition to a contact mechanical problem on closing fractures, a model which is further extended by Stefansson et al. [50] to additionally account for thermally induced tractions. Examples of non-conforming schemes for two-phase flow can be found in Fumagalli and Scotti [23], Tene et al. [52], and for coupled geomechanics and two-phase flow in Cusini et al. [14]. Beyond that, it should be mentioned that other approaches to modeling fractured porous media have been developed, such as phase-field models [56], [37] or discrete fracture networks in which the bulk porous medium is neglected [22], [11].

This work considers conforming mixed-dimensional methods for flow in fractured porous media, adopting the model concepts presented in Martin et al. [39], Boon et al. [10]. As mentioned above, several cell-centered schemes, among others, have been developed for the solution of the mixed-dimensional flow model. Cell-centered schemes have the advantage that they allow for a straightforward incorporation of the coupling conditions between the bulk porous medium and the fractures. However, in comparison to vertex-centered schemes they lead to significantly more degrees of freedom on unstructured simplex grids, that is, triangular or tetrahedral grids in two or three dimensions, respectively. Moreover, while the tpfa scheme is inconsistent on unstructured grids, mpfa[1] schemes exhibit comparatively large flux stencils, leading to computationally more expensive flux expressions and a large number of nonzero entries in the Jacobian matrices. Therefore, in this work we want to discuss potentially more efficient discretizations on the basis of the vertex-centered box scheme (see e.g. Hackbusch [27], Helmig [28]), and compare their performance to well-established approaches. In an early contribution, Reichenberger et al. [45] presented an extension of the box scheme to account for discrete lower-dimensional fractures, however, the underlying model involves a number of assumptions and simplifications that do not generally hold in fractured porous media, and which are incompatible with our model concept. Therefore, we introduce novel extensions to the box scheme that allow us to discretize our model equations and to incorporate the coupling conditions.

The comparisons are carried out in a number of test cases for incompressible single-phase flow in two and three dimensions with varying geometrical complexity. The discrete solutions are compared in terms of the primary unknown and the transfer fluxes between the bulk porous medium and the fractures. A good approximation of the fluxes is of crucial importance, for instance, when solving transport problems on the resulting flow fields. Our implementation is capable of performing transient transport simulations (see e.g. our contributions to Berre et al. [8] hosted at git.iws.uni-stuttgart.de/benchmarks/fracture-flow-3d), however, such investigations go beyond the scope of this work.

This paper is structured as follows. In section 2, we present the equi-dimensional mathematical model for flow in fractured porous media, where the fractures are treated as heterogeneities. Subsequently, the mixed-dimensional approximation of the model is presented in section 3, before we introduce our extensions to the box scheme in section 4. Section 5 contains the numerical examples, and we conclude the paper in section 6 with a summary and outlook. All schemes involved in the comparisons have been implemented into the same software framework, namely DuMux [18], [35], and the source code to all examples can be found at git.iws.uni-stuttgart.de/dumux-pub/glaeser2020b.

Section snippets

Equi-dimensional model

This section introduces the equi-dimensional model for flow in fractured porous media, in which we treat fractures, their intersections, and intersections of intersections as regions with different material properties than the surrounding medium. To this end, let us consider a domain ΩR3 with boundary ∂Ω. Furthermore, let Ωb, Ωf, Ωis and Ωj be disjoint partitions of Ω representing the bulk porous medium, the fractures, the intersections of fractures and the intersections of intersections,

Mixed-dimensional model

Following the ideas presented in Boon et al. [10], Nordbotten et al. [42], we consider a mixed-dimensional decomposition of fractured porous media, that is, the fractures and their intersections are described by entities of codimension 1 and 2, respectively. In this mixed-dimensional setting, we will from now on refer to the sub-domains with the subscript d indicating the dimension, i.e. the domain is decomposed into a 3-dimensional bulk domain Ω3, as well as the lower-dimensional domains Ω2, Ω1

Vertex-centered finite-volume schemes

In this section, we want to present two extensions to the box scheme, a vertex-centered finite-volume scheme (see e.g. Hackbusch [27], Helmig [28]), to discretize the mixed-dimensional equations (8a), (8b)-(11a), (11b). Section 4.1 presents the construction of the finite-volume discretizations of the subdomains, which is common to both schemes. Subsequently, the ebox-dfm scheme is presented in section 4.2, which uses additional degrees of freedom around the fractures to capture discontinuities

Numerical experiments

Ideally, a numerical method should be both efficient and accurate in order to be useful, for instance, in the design process of a geotechnical engineering application. Therefore, the test cases presented in the sequel aim at comparing the vertex-centered schemes outlined in the previous section to the well-established mpfa-dfm and tpfa-dfm schemes with respect to these two aspects. In addition, we include the vertex-centered box-dfm scheme [45] in the comparison. For an overview over the

Summary and outlook

In this work we have presented a comparison, by means of numerical experiments, of different finite-volume schemes for the simulation of flow in fractured porous media. The underlying mathematical model is based on a mixed-dimensional formulation in which the bulk medium, the fractures and intersections of fractures are described by geometries of decreasing dimensionality. On each of these subdomains, an individual system of partial differential equations is solved and the interaction between

CRediT authorship contribution statement

Dennis Gläser: Investigation, Methodology, Software, Visualization, Writing – original draft. Martin Schneider: Investigation, Writing – review & editing. Bernd Flemisch: Supervision, Writing – review & editing. Rainer Helmig: Supervision.

Declaration of Competing Interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

Acknowledgements

We thank the Deutsche Forschungsgemeinschaft (DFG, German Research Foundation) for supporting this work by funding SFB 1313, Project Number 327154368.

References (58)

  • K.D. Nikitin et al.

    Monotone embedded discrete fractures method for flows in porous media

    J. Comput. Appl. Math.

    (2020)
  • V. Reichenberger et al.

    A mixed-dimensional finite volume method for two-phase flow in fractured porous media

    Adv. Water Resour.

    (2006)
  • T. Sandve et al.

    An efficient multi-point flux approximation method for discrete fracture–matrix simulations

    J. Comput. Phys.

    (2012)
  • P. Schädle et al.

    3d non-conforming mesh model for flow in fractured porous media using Lagrange multipliers

    Comput. Geosci.

    (2019)
  • M. Tene et al.

    Projection-based embedded discrete fracture model (PEDFM)

    Adv. Water Resour.

    (2017)
  • M.F. Wheeler et al.

    Ipacs: integrated phase-field advanced crack propagation simulator. An adaptive, parallel, physics-based-discretization phase-field framework for fracture propagation in porous media

    Comput. Methods Appl. Mech. Eng.

    (2020)
  • I. Aavatsmark

    An introduction to multipoint flux approximations for quadrilateral grids

    Comput. Geosci.

    (2002)
  • J. Aghili et al.

    Two-phase Discrete Fracture Matrix models with linear and nonlinear transmission conditions

    GEM Int. J. Geomath.

    (2019)
  • T. Arbogast et al.

    Mixed finite element methods on nonmatching multiblock grids

    SIAM J. Numer. Anal.

    (2000)
  • R.L. Berge et al.

    Finite volume discretization for poroelastic media with fractures modeled by contact mechanics

    Int. J. Numer. Methods Eng.

    (2020)
  • F. Bonaldi et al.

    Gradient discretization of two-phase poro-mechanical models with discontinuous pressures at matrix fracture interfaces

  • W.M. Boon et al.

    Robust discretization of flow in fractured porous media

    SIAM J. Numer. Anal.

    (2018)
  • A. Borio et al.

    Comparison of the response to geometrical complexity of methods for unstationary simulations in discrete fracture networks with conforming, polygonal, and non-matching grids

    Comput. Geosci.

    (2021)
  • K. Brenner et al.

    Hybrid-dimensional modelling of two-phase flow through fractured porous media with enhanced matrix fracture transmission conditions

    J. Comput. Phys.

    (2018)
  • K. Brenner et al.

    Gradient discretization of hybrid-dimensional Darcy flow in fractured porous media with discontinuous pressures at matrix–fracture interfaces

    IMA J. Numer. Anal.

    (2016)
  • M. Cusini et al.

    Simulation of coupled multiphase flow and geomechanics in porous media with embedded discrete fractures

    Int. J. Numer. Anal. Methods Geomech.

    (2021)
  • T.A. Davis

    Algorithm 832: Umfpack v4.3—an unsymmetric-pattern multifrontal method

    ACM Trans. Math. Softw.

    (2004)
  • J. Droniou et al.

    Numerical analysis of a two-phase flow discrete fracture matrix model

    Numer. Math.

    (2019)
  • B. Flemisch et al.

    A review of the XFEM-based approximation of flow in fractured porous media

  • Cited by (3)

    View full text