Upscaling methods for effective permeability estimation of large digital rock models

In recent decades, Digital Rock Physics (DRP) has gained significant attention as an alternative to traditional experimental core analysis. Within the DRP workflow, the absolute permeability of a digital rock sample, reconstructed from micro-computed tomography data, is determined through direct numerical simulation (DNS) by solving the Navier-Stokes equations within the pore space. However, performing full-scale DNS on large digital rock models can demand prohibitive computational resources [1]. A common strategy to mitigate this is to simulate fluid dynamics in smaller subdomains, assuming they constitute a Representative Elementary Volume (REV). But, this is challenging for tight sandstones, an important type of unconventional reservoirs characterized by complex pore structures, narrow and tortuous flow channels, and pronounced heterogeneity. For such digital models, the adequate REV size often becomes so large that its direct simulation exceeds available computational capacities. This limitation necessitates the use of robust numerical upscaling methods to bridge the gap between detailed pore-scale simulations and the macro-scale flow properties of the entire rock.

This work implements and compares two upscaling approaches to predict the effective permeability of full-scale digital rocks: hierarchical homogenization and an analytical method based on an analogy between porous media and electrical resistor networks. Both methods are based on the domain decomposition into smaller subdomains and subsequently calculating local permeability via high-resolution pore-scale DNS. The pore-scale flow simulations are performed using an efficient GPU-accelerated finite difference solver based on the matrix-free relaxation method [2, 3]. The key distinction between the approaches lies in the integration step for obtaining the effective permeability: hierarchical homogenization estimates the overall permeability by considering the rock as a composite medium of homogenized cells, governed by upscaled parameters, while the resistor-network method employs analytical summation based on resistor-network summation rules. 

This work is supported by the Russian Science Foundation under grant 24-77-10022.

1. Yakovlev, M., & Konovalov, D. (2023). Multiscale geomechanical modeling under finite strains using finite element method. Continuum Mechanics and Thermodynamics, 35(4), 1223–1234.

2. Räss, L., Utkin, I., Duretz, T., Omlin, S., & Podladchikov, Y. Y. (2022). Assessing the robustness and scalability of the accelerated pseudo-transient method towards exascale computing. Geoscientific Model Development, 15(14), 5757–

3. Alkhimenkov, Y., & Podladchikov, Y. Y. (2025). Accelerated pseudo-transient method for elastic, viscoelastic, and coupled hydro-mechanical problems  with applications. Geoscientific Model Development, 18(2), 563–583.