Direct numerical simulations of immiscible two-phase flow in rough fractures: Impact of wetting film resolution

The immiscible displacement of a wetting fluid by a non-wetting fluid in rough fractures is crucial in many subsurface applications. Hydrodynamic-scale modeling of such drainage flows is challenging due to the complex interaction between the forces at play, the intricate geometry, and the required modeling of moving contact lines. In addition, a remaining critical open question is to what extent not resolving the films of wetting fluid deposited on fracture walls degrades numerical predictions. We address this question by solving the Navier–Stokes equations, employing the volume-of-fluid method to capture fluid–fluid interfaces and considering numerical meshes that result in either resolved films (RF) or unresolved films (UF) in the simulations. The numerical model, implemented in OpenFOAM, is validated in the classical Saffman–Taylor (ST) viscous instability configuration using the original ST experimental measurements; at capillary numbers (Ca) larger than 10−3, UF simulations overpredict ST finger widths. We then address two synthetic fracture geometries: one with sinusoidally varying apertures and one with stochastic geometric properties typical of geological fractures. Predictions of RF and UF simulations are compared quantitatively for Ca ranging between 10−5 and 10−3. Wetting film thicknesses follow a power law of Ca similar to Bretherton's law. RF and UF approaches both predict similar invasion patterns, but the latter underestimates interfacial lengths and macroscopic pressure drops, as compared to RF simulations, while overpredicting invading fluid saturations and breakthrough times. These discrepancies increase with Ca, whereas the disordered nature of the geological fracture tends to limit them. For Ca10−5, the discrepancies are negligible.