Two-branch retention behavior in unsaturated fractured rock driven by fracture-matrix flow partitioning

Abstract

Upscaling unsaturated flow in fractured rock remains challenging because fractures and matrix often exhibit sharply contrasting hydraulic behaviors across saturation states. Here, we demonstrate that unsaturated flow undergoes a transition between matrix- and fracture-dominated regimes. Three-dimensional direct numerical simulations reveal that both relative permeability and capillary pressure curves display a robust two-branch structure. We analytically derive a generalized retention formulation that identifies a critical saturation marking the transition between the two distinct retention regimes and reproduces the two-branch behavior captured in the numerical simulations. An analytical expression for the critical pressure head is further derived to represent the limiting case of fully connected fracture networks, providing a physical explanation for the retention regime shift and showing good agreement with the numerical results for systems above the percolation threshold. Our results provide a mechanistic framework for understanding and upscaling unsaturated flow in fractured rock, with broad implications for hydrology and geophysics.

Figures (4)

• Figure 1: Generated three-dimensional fracture networks associated with different power-law exponents of fracture radii $a$ and fracture intensities $\gamma$. Ten realizations are generated for each case, with only one representative realization presented here for illustration. The percolation parameter $\chi$ is also computed to quantify fracture network connectivity.
• Figure 2: Numerical simulation results of upscaled (a) relative permeability $k_r$ and (b) normalized capillary pressure $p_c/\rho gL$ as functions of effective saturation $S_e$ in a fractured rock with the power-law exponent $a = 3$ and fracture intensity $\gamma = 20/L$. Matrix and fracture contributions are also separately shown to illustrate their relative roles. The black cross indicates the critical saturation $S_c$, marking the transition between matrix-dominated and fracture-dominated flow regimes.
• Figure 3: Comparison of numerically simulated and analytically predicted relative permeability and capillary pressure curves. The scatter points represent simulation results of relative permeability $k_r$ and normalized capillary pressure $p_c/\rho gL$ as functions of effective saturation $S_e$, while the solid lines indicate the prediction by the two-branch analytical formulation.
• Figure 4: Numerical results (red dots) and analytical solution (black dashed line; Equation 16) for the dimensionless critical pressure head plotted against the percolation parameter $\chi$. The shaded gray area indicates the range of the percolation threshold $\chi_c \approx 0.7-2.8$.