Source-linked AI summary
Benchmarks for single-phase flow in fractured porous media
Bernd Flemisch, Inga Berre, Wietse Boon, Alessio Fumagalli, Nicolas Schwenck, Anna Scotti, Ivar Stefansson, Alexandru Tatomir
TL;DR
The paper addresses the need for systematic comparisons of numerical schemes for single-phase flow in fractured porous media. It proposes computational benchmarks spanning increasingly difficult fracture geometries and physical parameters, and compares several DFM-based methods. Across applicable cases, the methods agree well for permeable fractures, while blocking fractures expose method-specific failures and substantial differences in system condition numbers.
Problem
Benchmarking studies for fractured porous-media flow remain limited, despite challenges from discontinuous properties, complex fracture networks, and differing numerical models.
Method
The paper systematically compares several DFM-based numerical schemes using four computational benchmark cases and publicly shares grids and results.
Results
Methods agree quite well within their applicability ranges, accurately describing more-permeable fracture networks, whereas some fail to model blocking fractures correctly.
Takeaways & Limitations
The benchmarks provide a basis for distinguishing methods and inviting further comparisons by the scientific community.
Takeaways & Limitations
The benchmarks are restricted to simple physics and two-dimensional computational domains.
Abstract
from arXiv · showhide
This paper presents several test cases intended to be benchmarks for numerical schemes for single-phase fluid flow in fractured porous media. A number of solution strategies are compared, including a vertex and a cell-centered finite volume method, a non-conforming embedded discrete fracture model, a primal and a dual extended finite element formulation, and a mortar discrete fracture model. The proposed benchmarks test the schemes by increasing the difficulties in terms of network geometry, e.g. intersecting fractures, and physical parameters, e.g. low and high fracture-matrix permeability ratio as well as heterogeneous fracture permeabilities. For each problem, the results presented by the participants are the number of unknowns, the approximation errors in the porous matrix and in the fractures with respect to a reference solution, and the sparsity and condition number of the discretized linear system. All data and meshes used in this study are publicly available for further comparisons.
1 Introduction
Fractured porous media contain geometrically thin features with strongly discontinuous coefficients that can act as conduits or barriers, making flow prediction challenging. This work systematically compares DFM-based numerical schemes using computational benchmarks designed to distinguish methods and support further evaluation.
- Motivation: Fractures are thin anisotropic inclusions with material coefficients spanning several orders of magnitude.Their small normal size relative to tangential dimensions motivates representing them as fractures.
- Motivation: Fractures can act as conduits or barriers and severely affect flow patterns across applications including groundwater, energy, petroleum, waste, mining, and geothermal systems.
- Motivation: Flow prediction remains challenging despite decades of modeling development and many existing conceptual and numerical approaches.The paper distinguishes discrete fracture-matrix models from continuum fracture models.
- Related approaches: DFM approaches represent flow in both fractures and surrounding rock, using hybrid-dimensional or non-conforming formulations alongside classical matching-grid models.Examples include EDFM, XFEM-based approaches, and mortar-type methods.
- Benchmarking gap: Benchmarking verifies, tests, and compares modeling tools, but only a limited number of such studies exist in geoscience.
- Contribution: The paper systematically compares existing and new computational benchmarks and provides grid and result files for further method evaluation.The paper is accompanied by a Git repository containing the benchmark materials.
2 The model problem
The paper formulates incompressible single-phase flow through a porous medium using Darcy’s law in a domain containing a thin fracture network. Small fracture apertures permit reducing the full-dimensional fractures to a lower-dimensional network, leading to hybrid-dimensional formulations.
- Governing model: The model considers incompressible single-phase flow through a porous medium governed by Darcy’s law.
- Boundary conditions: The domain is bounded and subject to Dirichlet and Neumann boundary conditions on disjoint boundary portions.The boundary is partitioned as ∂D = ∂D_D ∪ ∂D_N with empty intersection.
- Variables: The formulation uses u for macroscopic fluid velocity and K and p for absolute permeability and pressure, respectively.
- Fracture reduction: When fracture aperture is small relative to other characteristic dimensions, the full-dimensional fracture domain is reduced to an (N−1)-dimensional fracture network.Different modeling choices in this reduction produce hybrid-dimensional problem formulations.
3 Participating discretization methods
The paper compares established discretization methods for hybrid-dimensional flow in fractured porous media, emphasizing their grid requirements and fracture–matrix coupling choices. The methods span conforming and non-conforming meshes, vertex- and cell-centered formulations, and several treatments of intersections and transmissibilities.
- Hybrid-dimensional models reduce thin, full-dimensional fractures to lower-dimensional networks embedded in the surrounding porous matrix.This avoids resolving the fracture interior width directly and supports lower-dimensional fracture elements.
- Box-DFM: Box-DFM combines finite-element and finite-volume grids in a vertex-centered, locally conservative scheme on unstructured meshes.Its control volumes are constructed from a dual finite-element/finite-volume mesh around mesh nodes.
- CC-DFM: CC-DFM uses cell-centered TPFA with fractures aligned to matrix-cell faces and harmonic-average face transmissibilities.The method assigns matrix-cell pressure unknowns and uses permeability, face geometry, and cell-center distances in transmissibility calculations.
- CC-DFM: CC-DFM removes intermediate fracture-intersection cells and uses star-delta transmissibilities to directly couple neighboring fracture cells.The elimination reduces the condition number and avoids small-cell time-step restrictions, but disregards intersection permeability information.
- EDFM: EDFM overlays independently generated fracture meshes on a background rock grid and computes approximate matrix–fracture and fracture–fracture transmissibilities.Fracture intersections do not determine the underlying fracture and rock grids, while transmissibilities use fracture-cell geometry, permeability, aperture, and intersection measures.
- Mortar-DFM, P-XFEM, and D-XFEM: Mortar-DFM couples matrix and fracture flows through a flux mortar variable and supports non-matching grids, including unified treatment across arbitrary dimensions.P-XFEM permits fracture vertices independent of the matrix grid but requires matching branch grids at intersections; D-XFEM uses a saddle-point velocity–pressure formulation for intersecting fractures with different permeabilities.
4 Benchmark Problems
Four benchmark cases compare the participating flow methods across increasingly complex fracture geometries and material settings. Results report discretization quality and linear-system characteristics, with line plots added for selected cases.
- Four benchmark cases compare the methods, progressing from intersecting fractures and a non-straight surface to regular, immersed, intersecting, and application-derived networks.
- Each case specifies the domain, boundary conditions, fracture geometry, aperture, and permeability, followed by method-specific grids.
- Reference-based cases evaluate errors separately in the porous matrix and fracture network.
- Matrix and fracture errors are complemented by linear-system density, condition number, and selected line-based solution comparisons.
4.1 Benchmark 1: Hydrocoin
The Hydrocoin benchmark models groundwater flow through a domain with two crossing fractures and prescribed head and no-flow boundaries. Most methods agree well with the reference, while fracture errors and especially condition numbers vary substantially.
- The modeled domain contains two crossing, highly permeable fractures and a non-straight top surface.
- The top boundary prescribes Dirichlet piezometric head, while the other three boundaries impose Neumann no-flow conditions.
- 10−6 m/s fracture conductivity contrasts with 10−8 m/s rock-matrix conductivity.
- Most participating methods show good agreement with the reference head profile at 200 m depth, although EDFM is somewhat offset.
- Matrix errors remain within narrow bounds, whereas fracture errors vary more across methods.
- Condition numbers are about 10^4 for Box-DFM, CC-DFM, and EDFM, but are five to seven orders of magnitude larger for P-XFEM, mortar-DFM, and D-XFEM.The larger values are attributed to the saddle-point nature of the latter methods.
4.2 Benchmark 2: Regular Fracture Network
Benchmark 2 compares methods on a regular fracture network under highly conductive and blocking fracture conditions. Conductive cases yield comparable results, whereas blocking fractures expose differences in discontinuity capture, errors, and conditioning.
- Conductive fractures: For conductive fractures, all methods produce pressure results relatively close to the reference solution along horizontal and vertical comparison lines.P-XFEM produces a more diffuse pressure profile in the vertical fracture.
- Conductive fractures: Conductive-fracture methods have comparable matrix and fracture errors, with similar sparsity; degrees of freedom and condition numbers provide the main distinctions.Mortar-DFM and D-XFEM are outliers with many degrees of freedom and high condition numbers, while P-XFEM has the highest condition number but significantly fewer degrees of freedom.
- Blocking fractures: Blocking-fracture results are compared along the diagonal line from (0.0, 0.1) to (0.9, 1.0), with errors, sparsity densities, and condition numbers reported in Table 6.The pressure field exhibits discontinuities associated with the low fracture permeability.
- Blocking fractures: For blocking fractures, Box-DFM and EDFM cannot capture pressure discontinuities and consequently show large matrix and fracture errors.Methods that handle discontinuities have similar errors, while mortar-DFM and P-XFEM show significantly improved condition numbers.
4.3 Benchmark 3: Complex Fracture Network
Benchmark 3 tests ten intersecting fractures with both conductive and blocking behavior under predominantly vertical and horizontal pressure gradients. The cases reveal stronger method differences when the flow configuration makes blocking effects more consequential.
- Setup: The network contains ten immersed fractures, with fractures 4 and 5 blocking and the others conductive; two pressure-gradient orientations define cases a and b.All fractures have aperture ε = 10^-4, with permeabilities 10^4 for conductive and 10^-4 for blocking fractures.
- Setup: P-XFEM could not participate because its implementation allows at most one fracture branch to cut each matrix-element face.Constructing a compliant matrix grid would conflict with the method’s independent fracture and matrix grids.
- Flow from top to bottom: In the top-to-bottom case, small intersection angles and mixed permeable-blocking intersections create geometric and coupling difficulties.Methods accounting explicitly for intersection permeability use a harmonic average for permeable-blocking crossings.
- Flow from top to bottom: For top-to-bottom flow, pressure-continuity methods show slightly higher matrix errors, although the difference is not particularly sharp.The average pressure gradient is almost parallel to the blocking fractures in this subcase.
- Flow from left to right: For left-to-right flow, the gap between continuous and discontinuous methods increases, although all errors remain of the same order of magnitude.Removing fracture-intersection cells makes CC-DFM ill-suited to mixed-permeability crossings; CC-DFM* substantially reduces error but loses the condition-number benefit of elimination.
4.4 Benchmark 4: a Realistic Case
Benchmark 4 applies the methods to a realistic fracture set interpreted from an outcrop on Sotra island. Despite mesh-generation difficulties, the computed solutions are generally similar, with EDFM showing greater robustness to avoiding unphysical connections.
- Case setup: The realistic case contains 64 fractures grouped into 13 connected networks in a 700 m × 600 m domain.The matrix permeability is 10^-14 m^2, fracture permeability is 10^-8 m^2, and fracture aperture is 10^-2 m.
- Case setup: Because some fractures were represented by multiple segments, the simulations replaced them with single segments, changing system connectivity.This modification was made to accommodate methods requiring one geometrical object per fracture.
- Results: The methods produce similar and comparable pressure solutions except in the top-right portion of the domain for Box-DFM.The solutions range from 0 to 101325 Pa in the matrix-pressure representation shown for the realistic case.
- Results: Mesh generation is the main concern because close fractures can create unphysical connections unless the meshes are carefully tuned.EDFM is the only method described as more robust with respect to this constraint.
- Results: Pressure profiles along y = 500 m and x = 625 m behave similarly across methods, while Box-DFM slightly overestimates some peaks.Observed oscillations are attributed to mesh effects.
5 Summary and Outlook
The benchmarks show generally good agreement among methods within their intended applicability ranges, while blocking fractures and implementation constraints expose important differences. The study remains limited to simple physics and two-dimensional domains, with future extensions proposed for coupled processes.
- Four benchmark cases generally produce good agreement when methods are used within their intended applicability ranges.
- Fracture networks more permeable than the surrounding matrix are accurately described by all methods.
- Blocking fractures reveal that some methods fail to predict the correct flow patterns.
- Condition numbers differ substantially across system matrices, but their effects on linear solvers are difficult to compare because solver requirements vary.
- The benchmarks are restricted to simple physics and two-dimensional computational domains.
- Future benchmark extensions may add transport, deformation, and reaction processes through broader international workshops.
B Domain modifications for Benchmark 1
Benchmark 1 modifies the original domain geometry by removing corner plateaus, adjusting fracture endpoints, and recalculating selected node positions.
- The modified domain omits the plateaus near the upper-left and upper-right corners.
- The upper ends of the two fractures are modified through changes to nodes 2–4 and 6–8.
- Nodes 16–19 are repositioned with higher precision in the modified setup.
C Fracture coordinates for Benchmark 3
The section directs readers to Table 12 for the coordinates of the modeled fracture points.
- Table 12 lists the coordinates for the modeled points.
- The coordinate data are presented as a tabulated reference for the section’s modeled geometry.
- Readers should use Table 12 to obtain the point coordinates.