Development of models and dedicated numerical methods for dynamics in fractured rocks is an active research field, with research moving towards increasingly advanced process couplings and complex fracture networks. The inclusion of coupled processes in simulation models is challenged by the high aspect ratio of the fractures, the complex geometry of fracture networks, and the crucial impact of processes that completely change characteristics on the fracture-rock interface. This paper provides a general discussion of design principles for introducing fractures in simulators, and defines a framework for integrated modeling, discretization, and computer implementation. The framework is implemented in the open-source simulation software PorePy, which can serve as a flexible prototyping tool for multiphysics problems in fractured rocks. Based on a representation of the fractures and their intersections as lower-dimensional objects, we discuss data structures for mixed-dimensional grids, formulation of multiphysics problems, and discretizations that utilize existing software. We further present a Python implementation of these concepts in the PorePy open-source software tool, which is aimed at coupled simulation of flow and transport in three-dimensional fractured reservoirs as well as deformation of fractures and the reservoir in general. We present validation by benchmarks for flow, poroelasticity, and fracture deformation in porous media. The flexibility of the framework is then illustrated by simulations of non-linearly coupled flow and transport and of injection-driven deformation of fractures. All results can be reproduced by openly available simulation scripts.

No abstract

Mixed-dimensional elliptic equations exhibiting a hierarchical structure are commonly used to model problems with high aspect ratio inclusions, such as flow in fractured porous media. We derive general abstract estimates based on the theory of functional a posteriori error estimates, for which guaranteed upper bounds for the primal and dual variables and two-sided bounds for the primal-dual pair are obtained. We improve on the abstract results obtained with the functional approach by proposing four different ways of estimating the residual errors based on the extent the approximate solution has conservation properties, i.e.: (1) no conservation, (2) subdomain conservation, (3) grid-level conservation, and (4) exact conservation. This treatment results in sharper and fully computable estimates when mass is conserved either at the grid level or exactly, with a comparable structure to those obtained from grid-based a posteriori techniques. We demonstrate the practical effectiveness of our theoretical results through numerical experiments using four different discretization methods for synthetic problems and applications based on benchmarks of flow in fractured porous media.

In this paper, we derive a posteriori error estimates for mixed-dimensional elliptic equations exhibiting a hierarchical structure. Exploiting the exterior calculus perspective of such equations, we introduce mixed-dimensional variables and operators, which, together with careful construction of the functional spaces, allow us to recast the set of partial differential equations as a regular linear elliptic problem structure-wise. We therefrom apply the well-established theory of functional a posteriori error estimates to our model to derive guaranteed abstract as well as fully computable upper bounds. Our estimators are tested using three different families of locally-mass conservative methods on synthetic problems and verification benchmarks of flow in fractured porous media. The numerical results support our theoretical findings while showcasing satisfactory effectivity indices.

In this paper, we present a model for saturated-unsaturated flow in fractured soils using an explicit representation of the fractures. The model is based on Richards' equation in the matrix and hydrostatic equilibrium in the fractures. While the first modeling choice is standard, the latter is motivated by the difference in flow regimes between matrix and fractures, i.e., the water velocity inside the fractures is considerably larger than in the soil even under saturated conditions. On matrix/fracture interfaces, the model permits water exchange between matrix and fractures only when the capillary barrier offered by the presence of air inside the fractures is overcome. Thus, depending on the wetting conditions, fractures can either act as impervious barriers or as paths for rapid water flow. Since in numerical simulations each fracture face in the computational grid is a potential seepage face, solving the resulting system of non-linear equations is a non-trivial task. Here, we propose a general framework based on a discretefracture matrix approach, a finite volume discretization of the equations, and a practical iterative technique to solve the conditional flow at the interfaces. Numerical examples support the mathematical validity and the physical applicability of the model.