Three-dimensional simulations of fracture dissolution

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.