Three‐dimensional simulations of fracture dissolution

Abstract

Numerical studies of fracture dissolution are frequently based on two‐dimensional models, where the fracture geometry is represented by an aperture field h(x,y). However, it is known that such models can break down when the spatial variations in aperture are rapid or large in amplitude; for example, in a rough fracture or when instabilities in the dissolution front develop into pronounced channels (or wormholes). Here we report a finite‐volume implementation of a three‐dimensional reactive transport model using the OpenFOAM® toolkit. Extensions to the OpenFOAM source code have been developed which displace and then relax the mesh in response to variations in the surface concentration; up to 100‐fold increases in fracture aperture are possible without remeshing. Our code has simulated field‐scale fractures with physical dimensions of about 10 m. We report simulations of smooth fractures, with small, well‐controlled perturbations in fracture aperture introduced at the inlet. This allows for systematic convergence studies and for detailed comparisons with results from a two‐dimensional model. Initially, the fracture aperture develops similarly in both models, but as local inhomogeneities develop the results start to diverge. We investigate numerically the onset of instabilities in the dissolution of fractures with small random variations in the initial aperture field. Our results show that elliptical cross sections, which are characteristic of karstic conduits, can develop very rapidly, on time scales of 10–20 years in calcite rocks.

The Development of Wormholes in Laboratory‐Scale Fractures: Perspectives From Three‐Dimensional Simulations

“…The phase diagram is bounded from below by a straight line with slope Pe D a eff ∞ = n S h L /2 h 0 that depends only on the fracture aspect ratio; here D a eff ∞ is the transport limit of D a eff . In Figure the largest experimentally accessible aspect ratios (2,000–4,000) are from Garcia‐Rios et al () and thus can reach higher D a eff ; simulations can access even larger aspect ratios (up to 10 5 ; Starchenko et al, ).…”

“…The simulation data from Figure have been replotted (open circles) using colors to represent the different dissolution patterns: face dissolution (turquoise), compact wormholes (red), necked wormholes (blue), and uniform dissolution (green). Additional simulation data with a rough initial aperture field are shown as open symbols: Elkhoury et al () —triangles and Starchenko et al ()—squares. Experimental data (solid symbols) are shown using the values for P e and D a eff determined in : (Detwiler et al, )—diamonds, (Elkhoury et al, )—triangles, (Garcia‐Rios et al, )—squares, and (Osselin et al, )—circles.…”

“…It has the effect of creating a single well‐defined channel running along the centerline of the fracture. We used a graded mesh to discretize the fluid volume with 260 × 300 × 8 cells; in previous work we have shown that this is sufficient for resolution‐independent dissolution (Starchenko et al, ).…”

“…Numerical simulations are developing an increasing capability to model dissolution, including more complex chemistry (Alt‐Epping et al, ; Steefel et al, ; Steefel & Maher, ) and more realistic descriptions of the three‐dimensional pore structure or fracture surfaces (Molins et al, , ; Starchenko et al, ; Soulaine et al, ; Soulaine & Tchelepi, ; Trebotich et al, ; Trebotich & Graves, ). Starchenko et al () described a fully three‐dimensional code based on the OpenFOAM library, which uses a mesh that conforms to the pore or fracture surfaces. Computationally, this is more complex than the representation of a mineral as a volume field (Molins et al, , ; Soulaine et al, ; Soulaine & Tchelepi, ; Trebotich & Graves, ), but it allows for a precise and unambiguous representation of the rock‐fluid interfaces.…”

Mineral Precipitation in Fractures: Using the Level‐Set Method to Quantify the Role of Mineral Heterogeneity on Transport Properties

“…For uniform surface mineralogy, it is relatively straightforward to convert local reaction rates into changes in fracture aperture by applying a 1-D alteration of the local domain (e.g., Andre & Rajaram, 2005;Cheung & Rajaram, 2002). However, when surface roughness increases, it is necessary to develop more rigorous approaches for tracking the fracture surface (Starchenko et al, 2016;Starchenko & Ladd, 2018;Yu & Ladd, 2010) and account for evolving reactive surface area (Molins, 2015).…”