# Micro-Continuum Modeling: An Hybrid-Scale Approach for Solving Coupled Processes in Porous Media

Cyprien Soulaine

► To cite this version:

Cyprien Soulaine. Micro-Continuum Modeling: An Hybrid-Scale Approach for Solving Coupled Processes in Porous Media. Water Resources Research, 2024, 60 (2), pp.106976. 10.1029/2023WR035908. [https://doi.org/10.1029/2023WR035908](https://doi.org/10.1029/2023WR035908)

HAL Id: hal-04491044

https://hal.science/hal-04491044v1

Submitted on 5 Mar 2024

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons CC BY 4.0 - Attribution - International License
# Water Resources Research®

## REVIEW ARTICLE

10.1029/2023WR035908

**Special Section:**

Impacts on Water Resources of Coupled Hydrological, Chemical, and Mechanical Processes in the Fractured Subsurface

**Key Points:**

* Micro-continuum models are hybrid-scale approaches for solving flow and transport in porous media

* State of the art micro-continuum models handle full Temperature-Hydrodynamics-Mechanics-Chemistry coupling in unsaturated environments

* Applications include Digital Rock Physics with sub-voxel porosity, pore-scale reactive transport, and poromechanics

**Correspondence to:**

C. Soulaine,
cyprien.soulaine@cnrs.fr

**Citation:**

Soulaine, C. (2024). Micro-continuum modeling: An hybrid-scale approach for solving coupled processes in porous media. *Water Resources Research, 60*, e2023WR035908. https://doi.org/10.1029/2023WR035908

Received 28 JULY 2023

Accepted 23 JAN 2024

© 2024. The Authors.

This is an open access article under the terms of the Creative Commons Attribution License, which permits use, distribution and reproduction in any medium, provided the original work is properly cited.

<figure><img src="image_2.png" /></figure>

## Micro-Continuum Modeling: An Hybrid-Scale Approach for Solving Coupled Processes in Porous Media

Cyprien Soulaine¹

¹Institut des Sciences de la Terre d'Orléans, CNRS, Université Orléans, BRGM, Orléans, France

**Abstract** Micro-continuum models are versatile and powerful approaches for simulating coupled processes in two-scale porous systems. Initially oriented for modeling static single-phase flow in microtomography images with sub-voxel porosity, the concept has been extended over the years to multi-phase flow, reactive transport, and poromechanics. This paper introduces an integrated micro-continuum framework to model coupled processes in porous media. It reviews state-of-the-art models and discusses applications in geosciences including Digital Rock Physics with sub-voxel porosity, moving fluid-solid interface at the pore-scale due to geochemical reactions, fracture-matrix interactions, and solid deformation. Finally, the paper discusses future developments in micro-continuum models.

### 1. Introduction

Critical zone and subsurface environments are synonyms of multi-scale and multiphysics processes. These two aspects constitute one of the biggest challenges in developing robust and accurate predictive modeling capabilities. For example, in fractured porous media, there are often orders of magnitude differences between the fracture aperture and the typical pore-size in the rock matrix. The multi-scale nature of geological porous materials that span a wide range of spatial scales often leads to the introduction of a cut-off length below which the microstructure exists but is unresolved. Even in high-resolution pore-scale imaging, some geometrical features are not resolved but still need to be considered in the modeling as they may strongly impact larger scales (Bultreys et al., 2015; Perez et al., 2022). Hybrid-scale models have been proposed to describe systems that include multiple characteristic length scales, for which some regions are resolved and described using pore-scale modeling while others are unresolved and modeled with continuum Darcy-based equations.

Hybrid-scale computational approaches proposed in the literature can be classified into two categories: the two-domain and the single-domain approaches. In the first, separate domains linked together through appropriate boundary conditions are used to solve the physics in resolved and unresolved regions (Molins et al., 2019). For example, fluid flow is described by Navier-Stokes equations in the resolved pore-scale domain, and by Darcy's law in the continuum unresolved region, and the two regions are coupled by Beavers-Joseph boundary conditions (Beavers & Joseph, 1967). Some hybrid-scale models adopt a discrete description of the pore-scale regions using Pore-Network Models (PNM) (Balhoff et al., 2007; Weishaupt et al., 2019, 2022). Although they are efficient, two-domain approaches become limited when the resolved-unresolved interface is moving because they require advanced dynamic mesh strategies.

Single-domain approaches, on the other hand, use a single set of partial differential equations throughout the computational domain regardless of the cell content (Soulaine & Tchelepi, 2016; Steefel et al., 2015). The conditions at the interface between the two domains are included in the partial differential equations and automatically satisfied. It means that interfacial conditions between the two regions are described as body source terms within the governing equations. Micro-continuum models belong to this hybrid-scale category. They arise from the volume integration of pore-scale equations. In addition to hybrid-scale modeling, the same formulation can also be used to solve full pore-scale and full continuum-scale problems, respectively.

Darcy-Brinkman-Stokes (DBS) equation (Brinkman, 1947) is the cornerstone of micro-continuum models for describing fluid flow both in resolved and unresolved regions using a single partial differential equation. Asymptotically, this equation tends toward Darcy's law in unresolved domains and Stokes in resolved domains, and captures automatically the Beavers-Joseph coupling conditions between the two regions (Neale & Nader, 1974). The underlying physics in the unresolved regions is described through appropriate sub-grid models. For example, the sub-grid geometry is described by concepts such as cell porosity, permeability, tortuosity, and
specific surface area. Although the approach was initially dedicated to single-phase flow, state-of-the-art micro-continuum models can now handle full Temperature-Hydrodynamics-Mechanics-Chemistry coupling in unsaturated environments using single-field equations and sub-grid models to describe multi-phase flow (Carrillo et al., 2020; Soulaine et al., 2019), solute and heat transport (Soulaine & Tchelepi, 2016), as well as solid deformation (Carrillo & Bourg, 2019).

Unlike two-domain approaches, micro-continuum models are also well-suited to capture the dynamic displacement of the interface between the porous and solid-free regions without involving complex re-meshing strategies. For example, micro-continuum models have been used successfully to simulate the formation and growth of wormholes in acidic environments (Golfier et al., 2002; Ormond & Ortoleva, 2000). Applications include Digital Rock Physics with sub-voxel porosity (Apourvari & Arns, 2014; Arns et al., 2005; Scheibe et al., 2015; Soulaine et al., 2016), biofilm growth in porous media (Kurz et al., 2022), simulation of gas flow in coal (Lanet et al., 2021), coke combustion (Xu, Dai, et al., 2022), and reactive transport modeling at the pore-scale (Soulaine et al., 2017) to name a few.

The objective of this paper is to introduce the main concepts of multiphysics micro-continuum modeling and to review applications of current interest in geosciences. The paper is organized as follows. First (Section 2), we introduce the philosophy of the micro-continuum concept and its mathematical background including governing equations, asymptotic behaviors, sub-grid models, and interfacial conditions between resolved and unresolved regions. Second (Section 3), we focus on the resolution of multiphysics coupled systems using micro-continuum models including multi-phase flow, heat transfer, reactive transport, and solid deformation. Then (Section 4), we present applications of micro-continuum models related to current research interests including Digital Rock Physics with sub-voxel porosity, coupled hydro-geochemistry modeling, and poromechanics. Finally (Section 5), we discuss open-challenges, and we close with a summary.

## 2. The Micro-Continuum Approach

### 2.1. A Description That Depends on Spatial Resolution

The concept of micro-continuum is intimately related to the spatial resolution of the porous domain considered. The base element is the distribution of the fluid volume fraction, $\phi$, in a cell of the computational grid. This volume fraction is therefore a cell porosity field. We have,

$$ \phi = \begin{cases} 1, & \text{if the cell is occupied by fluids only,} \\ ]0, 1[, & \text{if the cell contains a fluid-solid aggregate,} \\ 0, & \text{if the cell is occupied by solids only.} \end{cases} \qquad (1) $$

In the micro-continuum framework, the governing equations are defined in every cells regarding their content. Therefore, a small amount of fluid is always considered in the cells occupied by solid only. In practice, we use $\phi = 0 \equiv 10^{-3}$.

Three cases of interest are illustrated in Figure 1. On the one hand, the medium porosity is fully resolved on the grid and $\phi$ has only two values, namely 0 and 1 (see Figure 1a). This binary case corresponds to a pore-scale description with a sharp delineation between the fluid and the solid. The physics of flow is usually described by Navier-Stokes equations. On the other hand of the resolution spectrum (see Figure 1c), every cells of the grid contain a fluid-solid aggregate—that is, an unresolved porous medium—and $\phi$ is always in the range [0, 1]. There is always a tiny amount of solid—even a minuscule amount—in the grid blocks, and there is no cell containing fluid only. The physics of flow and transport is described by Darcy-based models. In between these two limit cases, there are intermediate situations for which some regions of the domain are fully resolved while others are unresolved (see Figure 1b).

### 2.2. The Darcy-Brinkman-Stokes Equation

Micro-continuum equations are formed by cell-averaging the fundamental equations of continuum mechanics along with fluid-solid interfacial areas within the cell. The integration of Navier-Stokes momentum equation over the cell volume, $V_{cell}$, gives (Whitaker, 1986a, 1986b),
<figure><img src="image_3.png" /><figcaption>**Figure 1.** Schematic representations of a porous medium with two characteristic pore sizes depending on the scale of resolution: (a) full pore scale (Navier-Stokes), (b) intermediate or hybrid scale, and (c) full continuum scale (Darcy). The micro-continuum framework describes the physics of flow and transport using on a single set of equations resolved throughout the entire system regardless the spatial resolution.</figcaption></figure>

$$ \frac{1}{\phi} \left( \frac{\partial \rho \bar{v}}{\partial t} + \nabla \cdot \left( \frac{\rho \bar{v}}{\phi} \right) \right) = - \nabla \bar{p} + \rho g + \nabla \cdot \frac{\mu}{\phi} (\nabla \bar{v} + \nabla^t \bar{v}) - \mu k^{-1} (\bar{v} - \bar{v}_s) + F, \quad (2) $$

where $\bar{v} = \frac{1}{V_{cell}} \int_{V_j} \nu dV$ and $\bar{p} = \frac{1}{\bar{v}_j} \int_{\bar{v}_j} \rho dV$ are the cell-averaged velocity and pressure, respectively ($V_f$ being the volume occupied by fluids within $V_{cell}$). The left-hand side corresponds to inertia effects that are often neglected in subsurface flows. This momentum equation is an extension of the Darcy-Brinkman-Stokes (DBS) equation proposed by Brinkman (1947). The first and second terms of the right-hand side are the hydrodynamic driving forces, namely the pressure gradient and the gravitational acceleration. The third term describes the viscous dissipations due to the mutual friction of the fluid elements with each other. The fourth term is a drag force that describes the friction of the unresolved microstructure with the fluid. It can be seen as a subgrid model. This flow resistance is characterized by the skin factor—also referred to as permeability—$k^{-1}$ that is related to the geometry of the unresolved microstructure. $\bar{v}_s$ is the velocity of the solid structure, for example, because of poromechanical effects. If the structure is motionless, $\bar{v}_s = 0$. Finally, the last term of the right-hand side is a body source term (e.g., surface tension forces, magnetic forces).

The extended DBS equation is generic and valid throughout the computational grid, regardless of the cell content. Nonetheless, the form of the drag force and of the body source term differs whether they describe the flow in a resolved or an unresolved region. If there is no solid in resolved regions there is no drag force. Hence, the skin factor writes,

$$ \mu k^{-1} = \begin{cases} 0, & \text{if the cell is occupied by fluids only,} \\ \mu k^{-1}, & \text{if the cell contains a fluid-solid aggregate.} \end{cases} \quad (3) $$

We will see in Section 2.5 that it exists constitutive laws such as Kozeny-Carman permeability-porosity relationship (Carman, 1937; Kozeny, 1927) for which the inverse of the permeability asymptotically tends toward zero if $\phi = 1$.

The body source force has a different formulation in the resolved and unresolved regions. We simply write,

$$ F = \begin{cases} F_{\text{resolved}}, & \text{if the cell is occupied by fluids only,} \\ F_{\text{unresolved}}, & \text{if the cell contains a fluid-solid aggregate.} \end{cases} \quad (4) $$
For example, for immiscible two-phase flow, the surface tension forces are modeled by computing the interface of the curvature in the resolved region and the capillary pressure gradient in the unresolved region (Carrillo et al., 2020; Soulaine et al., 2019).

In practice, the DBS momentum equation is supplemented by, at least, the fluid mass balance,

$$ \frac{\partial \phi \rho}{\partial t} + \nabla \cdot (\rho \bar{v}) = 0. \qquad (5) $$

Section 3 reviews the additional equations that need to be considered for modeling coupled processes including two-phase flow, heat transfer, reactive transport, and solid deformation.

## 2.3. Asymptotic Behavior

### 2.3.1. Asymptotic Behavior of the DBS Equation

The most significant aspect of the DBS equation is its asymptotic behavior across regions. On the one hand, in clear fluid regions, $\phi = 1$ and the drag force vanishes, and Equation 2 becomes the Navier-Stokes momentum balance,

$$ \frac{\partial \rho \bar{v}}{\partial t} + \nabla \cdot (\rho \bar{v} \bar{v}) = -\nabla \bar{p} + \rho g + \nabla \cdot \mu (\nabla \bar{v} + \nabla^t \bar{v}) + F_{\text{resolved}}, \quad \text{if } \phi = 1. \qquad (6) $$

On the other hand, in porous regions, the drag force is dominant. Because inertia effects and the viscous forces are negligible with respect to the drag force (Auriault, 2009; Tam, 1969; Whitaker, 1986a), Equation 2 tends asymptotically toward

$$ 0 = -\nabla \bar{p} + \rho g - \mu k^{-1} (\bar{v} - \bar{v}_s) + F_{\text{unresolved}}, \quad \text{if } \phi \in ]0, 1[, \qquad (7) $$

which is nothing else than Darcy's law.

Hence, thanks to these two limiting behaviors, the DBS equation, Equation 2, is a unique partial differential equation that solves for flow in computational grids regardless of the spatial resolution (and therefore the content of a grid block). This feature is the cornerstone of the micro-continuum approach, as it offers the possibilities for modeling processes in hybrid-scale domains that contain both resolved and unresolved regions.

### 2.3.2. Pore-Scale Modeling With Immersed Boundaries

If the porous matrix has sufficiently low permeability, fluid velocities in the porous domain are near zero and a no-slip condition is recovered at the interface between solid-free and porous regions (Angot et al., 1999; Khadra et al., 2000). This enables the use of micro-continuum simulations at the pore scale using a penalized approach, that is, the solid phase is described as a low-permeability porous medium (Soulaine & Tchelepi, 2016). Solving flow and transport through complex geometries using the penalized approach does not require complex grids—nor advanced software for meshing complex topology—as it relies on the mapping of the phase indicator function, $\phi$, on regular grids only.

Various techniques have been proposed to describe the immersed boundary conditions at the fluid-solid interface for multiphase flows (Carrillo et al., 2020; Horgue et al., 2014; Sharaborin et al., 2021; Soulaine et al., 2018) and chemical reactions at the solid surface (Maes et al., 2022; Molins et al., 2020; Soulaine et al., 2017). They will be discussed in the following.

Using immersed boundaries on regular grids is a great advantage for pore-scale simulations that involve the displacement of a fluid-solid interface due to geochemical processes such as dissolution and precipitation (Soulaine et al., 2017; Yang, Stack, et al., 2021) or the mechanical deformation of the solid structure (Carrillo & Bourg, 2019).
<figure><img src="image_4.png" /><figcaption>**Figure 2.** Interfacial conditions for single-phase flow between resolved and unresolved regions after Neale and Nader (1974). (Left) Actual velocity profile for coupled parallel flow within a channel and a bounding porous medium (Right) Velocity profile according to the slip-flow hypothesis of Beavers and Joseph (1967).</figcaption></figure>

## 2.4. Interfacial Conditions Between Resolved/Unresolved Regions

In hybrid-scale simulations, both clear-fluid (resolved) and porous (unresolved) regions exist concomitantly in the computational grid. The coupling conditions at the interface between the two regions is not straightforward because the characteristic length of the unresolved region may be several orders of magnitude lower than the one of the resolved porosity. As a consequence, there is no necessary continuity of the physical variables across the interface separating the two regions. In micro-continuum models, the transition between the two regions is achieved by the spatial variation of porosity and permeability (see Figure 1b), and the coupling conditions are included directly within the partial differential equations.

For single-phase flow, the coupling conditions at the interface between a porous and a clear-fluid regions consist of: (a) the continuity of the normal velocities, (b) and a partial slip for the tangential velocities postulated empirically by Beavers and Joseph (1967). This is illustrated in Figure 2, left. The Beavers-Joseph coupling condition for two-region approaches reads, $\frac{\partial v_{res.}}{\partial y} = \frac{a_{BJ}}{\sqrt{k}}(v_{res.} - v_{unres.})$ where $a_{BJ}$ is a slip length. Being single-domain approaches, micro-continuum models include the interfacial conditions directly in the governing equations. Single-phase micro-continuum models satisfy these conditions because, (a) the continuity equation imposes a continuity of the normal component of the velocity field, (b) Neale and Nader (1974) have demonstrated that, for flow parallel to the fluid-porous interface, the DBS equation captures the slip length induced by the continuity of the viscous stresses between the two regions as postulated by Beavers and Joseph (1967). This is illustrated in Figure 2, right. The validity of DBS if the free flow is not oriented parallel to the interface is still an open question (Yang et al., 2019).

For two-phase flow, the derivation of a comprehensive interfacial condition between the two regions still constitutes a major challenge to solve. Indeed, in addition to the momentum exchange discussed for single-phase flow, the discontinuity in porosity also leads to a change in the form of the surface tension force. In the pure fluid regions, the surface tension force is associated with the curvature of the interface (e.g., leading to Laplace pressure in droplets). In the porous media, the surface tension effects are described by capillary pressure curves. The discontinuity in surface tension force at the resolved-unresolved interface is treated by assuming that the fluid-fluid interface of a droplet on a porous substrate forms a contact angle $\theta$ with the solid surface (see Figure 3). The contact angle is an upscaled parameter that depends on various sub-grid scale properties including interfacial
<figure><img src="image_5.png" /><figcaption>**Figure 3.** Surface tension forces acting on a two-fluid system. In the resolved region, surface tension forces relate to the curvature of the interface and can be described numerically using Brackbill's Continuum Surface Force. In the unresolved region, these forces are modeled using the concept of capillary pressure curves. At the fluid-porous interface, the interface curvature is modified locally to enforce a contact angle with the solid substrate.</figcaption></figure>

energies, surface roughness, and the presence of thin precursor films (Cassie & Baxter, 1944; Meakin & Tar-
takovsky, 2009; Wenzel, 1936). In current numerical models, the contact angle is imposed by locally modifying
the orientation of the fluid-fluid interface relative to the solid surface (Horgue et al., 2014; Soulaine et al., 2018).
Further work is still needed to improve this interfacial condition. The recent advances in coupling two-fluid Pore-
Network Models with Navier-Stokes-based solvers (Veyskarami et al., 2023) could inspire the development of an
interfacial condition accounting for the accurate force balance for arbitrary flow directions.

To model scalar transport (species concentration or heat) across the interface between resolved and unresolved regions, two coupling conditions must be satisfied: (a) a continuity of the mass/heat flux, *F*, in the direction normal to the interface, (b) a thermodynamic condition. The second can be an equilibrium condition, for example, equality of temperature and concentration, $C_{\text{res.}} = C_{\text{unres.}}$, or a non-equilibrium condition, $F = \alpha (C_{\text{res.}} - C_{\text{unres.}})$. The kinetic parameter $\alpha$ (in m/s) is a physical property of the unresolved region that controls the delay to reach the equilibrium. The non-equilibrium condition is more generic and has two asymptotic regimes whether $\alpha$ is large or small, that is, the time-scale to reach the equilibrium is very high or very low. On the one hand, if $\alpha \to \infty$, the time delay to reach equilibrium is so small that the equilibrium condition, $C_{\text{res.}} = C_{\text{unres.}}$, is reached instantaneously. On the other hand, if $\alpha \to 0$, then $F \to 0$ and the transport in the two regions is uncoupled and the interface behaves as an inert impervious barrier. In micro-continuum models including scalar transport, the continuity of fluxes across the interface, and the equilibrium conditions are always satisfied. On-going works focus on the development of micro-continuum scalar transport equations that automatically handle non-equilibrium conditions at the interface. It includes the modeling of surface reactions using kinetic reaction rates (Maes et al., 2022).

2.5. Sub-Grid Models

In the micro-continuum framework, the physics of flow and transport in the unresolved regions is filtered and modeled through the concept of sub-grid models. Sub-grid models describe the processes occurring below the grid-size. Because they are not resolved, these processes have to be modeled through constitutive laws that are provided to the micro-continuum model. In practice, sub-grid models rely on the classic concept used for modeling flow and transport in porous media including permeability, specific surface area, tortuosity, and dispersion tensor.

In the momentum equation, Equation 2, for example, the parameter *k*, is a sub-grid model that characterizes the flow resistance due to the sub-grid microstructure. It is described as cell permeability. There is no unique formulation for the sub-grid permeability but a large variety of models that depend on the phenomena involved. The non-exhaustive permeability models shown in Table 1 illustrate the large variety of sub-grid models. It can be constant for low flow rates or varying with the velocity field for fast flow rates due to the Forchheimer correction
<table><caption>**Table 1**
Some Sub-Grid Models Used to Describe the Cell Permeability in Micro-Continuum Simulations</caption><thead><tr><th>Name</th><th>Expression</th><th>Comments</th></tr></thead><tbody><tr><td>Constant</td><td>k = k<sub>0</sub></td><td>k<sub>0</sub> = cell permeability</td></tr><tr><td>Uniform capillary bundle</td><td>k = φ<sup>d/2</sup></td><td>d is the pore-throat diameter</td></tr><tr><td>Kozeny-Carman</td><td>k = <sup>d<sup>2</sup></sup>&frasl;<sub>180</sub> (1-φ)<sup>n/2</sup></td><td>Carman (1937) and Kozeny (1927)</td></tr><tr><td>Power-law</td><td>k = k<sub>0</sub>(φ/φ<sub>0</sub>)<sup>n</sup></td><td>n = user defined variable</td></tr><tr><td>Kozeny-Carman variation</td><td>k = k<sub>0</sub>(φ/φ<sub>0</sub>)<sup>d/2</sup></td><td>used for moving fluid-porous boundary</td></tr><tr><td>Two-phase</td><td>k = k<sub>0</sub>(<sup>1</sup>&frasl;<sub>k<sub>w</sub></sub> + <sup>1</sup>&frasl;<sub>k<sub>ro</sub></sub>)<sup>-1</sup></td><td>k<sub>r<sub>w</sub></sub>(S<sub>w</sub>) and k<sub>r<sub>ro</sub></sub>(S<sub>w</sub>) are relative permeabilities</td></tr><tr><td>Klinkenberg</td><td>k = k<sub>0</sub>(1 + <sup>b<sub>K</sub></sup>&frasl;<sub>p</sub>)</td><td>b<sub>K</sub> is a model parameter</td></tr><tr><td>Forchheimer</td><td>k = <sup>k<sub>0</sub></sup>&frasl;<sub>1 + α<sub>forch</sub>Re<sup>d-1</sup></sub></td><td>α<sub>forch</sub> and p are model parameters. Re = Reynolds number</td></tr></tbody></table>

to Darcy's law (Forchheimer, 1901; Soulaine & Quintard, 2014). If the cell contains several fluid phases, then the permeability depends on the fluid saturation using the relative permeabilities (Soulaine et al., 2019). For sub-micron pore-throats, nanoscale effects are no longer negligible and the permeability is modeled using Klinkenberg's law (Klinkenberg, 1941). This sub-grid model allows the simulation of shale gas recovery from nanoporous kerogen (Soulaine et al., 2019).

The body source term, *F*, has often a different formulation in the resolved and unresolved regions. For example, the surface tension force is described in the resolved domain by computing the curvature of the fluid-fluid interface (Brackbill et al., 1992) and in the unresolved region using the concept of capillary pressure curves including Brooks and Corey (1964) and van Genuchten (1980) models.

For solute transport in the unresolved regions, the concepts of tortuosity and hydrodynamic dispersion apply. The specific surface area in a grid-block is a key sub-grid model in the case of geochemical reactions (see Table 2).

### 3. Micro-Continuum Models for Coupled Processes

Micro-continuum models have the ability to model coupled processes. This is achieved by combining the DBS momentum equation with additional physics. In the following, we discuss coupled micro-continuum models including two-phase flow, heat transfer, geochemistry, and poromechanics.

#### 3.1. Two-Phase Flow Micro-Continuum Models

The transport of two immiscible fluids in porous media (e.g., air and water, oil and water) involves moving interfaces whose description varies significantly in the resolved (solid-free) and unresolved regions. On the one hand, in solid-free regions, the fluid-fluid interface is fully resolved and described by a sharp delineation between the two fluids. Fluid flow is then described by Navier-Stokes equations, and the curvature of the interface leads to discontinuity in the pressure field. Contact angles apply at the junction between the solid substrate and the fluid-

<table><caption>**Table 2**
Some Sub-Grid Models Used in Micro-Continuum Models to Describe the Cell Specific Surface Area</caption><thead><tr><th>Name</th><th>Expression (in m<sup>-1</sup>)</th><th>Comments</th></tr></thead><tbody><tr><td>Constant</td><td>A<sub>e</sub> = A<sub>0</sub></td><td>—</td></tr><tr><td>Power-law</td><td>A<sub>e</sub> = A<sub>0</sub>(Y<sub>s,i</sub>/Y<sub>0</sub>)<sup>n</sup></td><td>Usually, n = -2/3 and n = 2/3 for consolidated and unconsolidated media, respectively.</td></tr><tr><td>Sugar-lump</td><td>A<sub>e</sub> = (A<sub>0</sub> + A<sub>m</sub>(1 - (Y<sub>s,i</sub>/Y<sub>0</sub>)<sup>n<sub>1</sub></sup>)<sup>n<sub>2</sub></sup>) (Y<sub>s,i</sub>/Y<sub>0</sub>)<sup>n<sub>3</sub></sup></td><td>Surface area evolution of a dissolving aggregate (Noiriel et al., 2009). A<sub>m</sub> is the sum of the surface area of all individual particles. n<sub>1</sub>, n<sub>2</sub>, and n<sub>3</sub> are user-defined parameters.</td></tr><tr><td>Volume Of Solid</td><td>A<sub>e</sub> = |∇Y<sub>s,i</sub>|ψ</td><td>Compute the local surface area based on the mineral mapping, Y<sub>s,i</sub>, for pore-scale simulations (Soulaine et al., 2017).</td></tr></tbody></table>
fluid interface. On the other hand, in unresolved regions, the water saturation is governed by two-phase Darcy's law with relative permeabilities and capillary pressure (Muskat, 1949). Micro-continuum models for two-phase flow are based on a combination of the Volume of Fluid (VOF) approach to simulate Navier-Stokes-based two-phase flow (Hirt & Nichols, 1981) and the DBS equation.

The wetting phase saturation, $S_w$, is an additional unknown of the system. In resolved regions, it corresponds to the water phase indicator function ($S_w = 1$ in the wetting phase, $S_w = 0$ in the non-wetting phase, $S_w \in [0; 1]$ in cells containing the fluid-fluid interface) while in unresolved regions, it describes the water content ($S_w \in [0; 1]$). In the micro-continuum momentum equation, Equation 2, the velocity, $\bar{v}$, is, therefore, a mixture velocity, and the fluid viscosity and density are described as single-field variables (e.g., $\rho = S_w\rho_w$ + $(1 - S_w)\rho_{nw}$). The partial differential equation for the wetting phase saturation reads (Carrillo et al., 2020; Soulaine et al., 2018, 2019),

$$ \frac{\partial \phi S_w}{\partial t} + \nabla \cdot (\bar{v} S_w) + \nabla \cdot (\phi S_w (1 - S_w) \bar{v}_r) = 0, \quad (8) $$

where $\bar{v}_r = (\bar{v}_w - \bar{v}_{nw})$ is a multi-scale parameter that describes the average relative velocity between the two fluid phases in a grid block. This term has different forms in the resolved and the unresolved regions. In resolved regions, the relative velocity is described as a velocity normal to the fluid-fluid interface that sharpens the fluid-fluid interface, $\bar{v}_r \propto \max(\bar{v})\nabla S_w$, in agreement with the algebraic Volume of Fluid approach (Rusche, 2003). In unresolved regions, Carrillo et al. (2020) derived the appropriate expression of the relative velocity if the phase-averaged velocities, $\bar{v}_w$ and $\bar{v}_{nw}$ are described by Darcy's laws with relative permeabilities. It means that in unresolved regions, Equation 8 tends toward the classic saturation equation used in porous media.

In two-fluid systems, the sub-grid models also depend on saturation. On the one hand, the permeability, $k$, in Equation 2 relies on relative permeabilities (Soulaine et al., 2019) as shown in Table 1. On the other hand, the surface tension effects between the two immiscible fluids are modeled as a body source term $F$ in the DBS momentum equation. In the resolved region, $F$ corresponds to the Continuum Surface Force computed from the curvature of the interface (Brackbill et al., 1992). In the unresolved region, $F$ depends on the gradient of the capillary pressure curves (Carrillo et al., 2020). Relative permeabilities and capillary pressure curves are sub-grid models classically described using van Genuchten (1980) and Brooks and Corey (1964) formula.

The interfacial condition between resolved and unresolved regions is discussed in Section 2.4 and illustrated in Figure 3.

Carrillo et al. (2020) did intensive verification of the two-phase micro-continuum model. They verified that the model tends asymptotically toward Navier-Stokes solutions in solid-free regions and toward two-phase Darcy's solution in unresolved regions in the presence of viscous, gravity, and capillary effects. They also verify that interfacial conditions were well-captured if the porous region represents a fictitious low-porosity low-permeability medium.

If most of the current models rely on Volume of Fluid approaches to track the fluid-fluid interface in the resolved regions, recent work (Lisitsa et al., 2023) proposes a two-fluid micro-continuum approach based on a Cahn-Hilliard phase-field framework. This approach is promising to simulate systems that include phase change.

## 3.2. Heat Transfer With Micro-Continuum Approaches

Micro-continuum models for heat transfer rely on a thermal equilibrium model in the unresolved regions, that is, a single-field temperature, $T$,—corresponding to a mixture of the fluid and solid temperatures—describes the heat evolution in the system (Soulaine & Tchelepi, 2016). The governing equation writes,

$$ (\phi(\rho_f C p_f) + (1 - \phi)(\rho_s C p_s)) \frac{\partial T}{\partial t} + (\rho_f C p_f) \bar{v} \cdot \nabla T = \nabla \cdot (\lambda \nabla T) + Q_T, \quad (9) $$

where $(\rho_f C p_f)$ and $(\rho_s C p_s)$ are the fluid and solid heat capacities, $\lambda$ is the single-field thermal conductivity, and $Q_T$ is a heat source/sink term (e.g., enthalpy of reaction). The thermal conductivity, $\lambda$, is a mixture of the fluid and
solid conductivities, $\lambda_f$ and $\lambda_s$. It can be described using linear (Alazmi & Vafai, 2001) or harmonic averages (Maes & Menke, 2022).

In the resolved regions, Equation 9 tends to the fluid temperature equation mainly driven by advection. In low-porosity low-permeability unresolved regions, Equation 9 tends to the solid temperature equation driven by thermal conductivity. Other formulations of the micro-continuum heat transfer model are found in the literature. For example, the enthalpy-porosity models proposed by Voller (2009) and Voller et al. (1987) to model the melting/solidification of a solid including mushy zones belong to the family of micro-continuum models.

Notes that the thermal equilibrium base hypothesis of the heat transfer micro-continuum model is limiting in the case of advection-dominated regimes in unresolved regions.

### 3.3. Modeling Hydro-Bio-Geochemical Processes

Micro-continuum models for hydro-geochemical processes describe the evolution of the volume fraction for every $N_s$ minerals, $Y_{s,i}(x, y, z, t)$, in the computational grid. This is achieved using

$$ \rho_{s,i} \frac{\partial Y_{s,i}}{\partial t} = -\dot{m}_{s,i} \quad \text{for } i = 1, \dots, N_s, \qquad (10) $$

where $\rho_{s,i}$ is the mineral density and $\dot{m}_{s,i}$ is the rate of solid change. The latter is calculated from the aqueous species transport using the mass balance equations,

$$ \frac{\partial \phi C_j}{\partial t} + \nabla \cdot (\bar{\nu} C_j) - \nabla \cdot (\phi D_j^* \nabla C_j) = \sum_{i=1}^{N_s} \dot{m}_{f_{j,i}} \qquad (11) $$

where $C_j$ is the concentration of species $j$, $D_j^*$ is the dispersion tensor including molecular diffusion, tortuosity and hydrodynamic dispersion, and $\dot{m}_{f_{j,i}}$ describes chemical reactions at the fluid-rock interface. Here, the species diffusion is described by Fick's law. More sophisticated models such as Maxwell-Stefan equations can be used to model the multicomponent diffusion (Ahmadi et al., 2021).

Surface reactions are often modeled as kinetic laws proportional to the mineral surface area per unit of volume, $A_e$. In micro-continuum models, kinetic law and $A_e$ are two distinct elements. The kinetic laws can be provided by any geochemical package including PHREEQC (Parkhurst & Wissmeier, 2015), Reaktoro (Leal, 2015) or CrunchFlow (Steefel, 2009) regardless of the resolved or unresolved nature of the grid block (Steefel et al., 2015). Therefore, if one knows how to compute surface area per unit of volume in unresolved and resolved regions, one can perform reactive transport modeling at various scales of interest (pore-scale, hybrid-scale, Darcy-scale) (Soulaine et al., 2021). In practice, it is assessed with different sub-grid models whether the surface reaction occurs within the unresolved region or at the interface between the resolved and unresolved domains. In the first case, sub-grid models estimate the accessible reactive surface area as a function of the solid volume fraction and local flow rate (Noiriel et al., 2009). In the other case, the surface area is computed directly from the delineation between the two regions using the so-called Volume of Solid approach, that is, by computing the gradient of the solid volume fraction (Soulaine et al., 2017). A summary of common surface area models is given in Table 2.

The porosity field is updated using:

$$ \phi = 1 - \sum_{i=1}^{N_s} Y_{s,i}. \qquad (12) $$

The variation of the porosity field informs the displacement of the resolved-unresolved interface by updating the local permeability in Equation 2.
### 3.4. Solid Deformation Modeling: Darcy-Brinkman-Biot

The Darcy-Brinkman-Biot model proposed by Carrillo and Bourg (2019) describes hydro-mechanical coupling in two-scale deformable porous media. In addition to the micro-continuum fluid momentum and continuity equations, their model considers a mass balance and a momentum equation for the solid phase as well. For plastic and elastic solids, it writes,

$$ \frac{\partial \phi_s \rho_s}{\partial t} + \nabla \cdot (\phi_s \rho_s \bar{v}_s) = 0, \qquad (13) $$

and,

$$ -\nabla \cdot \sigma = \phi_s \nabla \cdot \tau + \phi_s \mu k^{-1} (\bar{v} - \bar{v}_s), \qquad (14) $$

where $\rho_s$ is the solid density, $\bar{v}_s$ is the solid velocity, $\tau$ is the Terzaghi stress tensor defining the effective stress due to both fluid pressure and confining pressure, and $\sigma$ is the plastic (or elastic) stress tensor. Note that the last term is a mutual momentum transfer due to the viscous resistance (i.e., drag) of the fluid on the solid surface.

According to poromechanics theory, the Terzaghi stress tensor is defined as,

$$ \tau = \sigma_{\text{conf}} - \bar{p} l, \qquad (15) $$

where $\sigma_{\text{conf}}$ is the confining pressure and $\bar{p}$ is the fluid pressure. For modeling clay swelling, Carrillo and Bourg (2019) also add a swelling pressure.

For compressible plastic solid, the stress tensor is,

$$ \sigma = \phi_s \mu_s^{\text{eff}} \left( \nabla \bar{v}_s + \nabla^t \bar{v}_s - \frac{2}{3} (\nabla \cdot \bar{v}_s) l \right), \qquad (16) $$

where the effective solid viscosity $\mu_s^{\text{eff}}$ can be described using various non-Newtonian plastic viscosity models.

## 4. Applications in Geosciences

In this section, cutting-edge micro-continuum simulations of coupled processes in porous media are reviewed. They result from combinations of the base elements introduced in Section 3. First, we highlight the ability of micro-continuum models to account for subvoxel porosity in image-based simulations. Then, we review recent work on multi-scale reactive transport modeling and solid deformation.

### 4.1. Digital Rock Physics and Micro-Continuum Models

Digital Rock Physics—a new scientific discipline that emerged from the development of High-Performing Computing combined with high-resolution imaging—aims at the evaluation of rock properties based on flow simulations (Andrä et al., 2013a, 2013b; Blunt et al., 2013). For example, by volume averaging the velocity solution of a Stokes flow in a porous rock sample imaged with micro-CT, one can get an estimate of the absolute permeability. Standard approaches, however, are limited because some pores are not always resolved in the image (Sadeghnejad et al., 2021). Micro-continuum models push the limit of current Digital Rock Physics capabilities by considering sub-voxel porosity and coupled processes.

#### 4.1.1. Digital Rock Physics With Subvoxel Porosity

Sub-resolution porosity (SRP) is ubiquitous in Digital Rock Physics. SRP has two origins: (a) it results from the trade-off between image resolution and field-of-view, and in this case, the imaging resolution is larger than the smallest pores, (b) the imaging instrument resolution is larger than the smallest pores that are therefore not visible in the image. These situations are common in heterogeneous rock samples including carbonates. In both cases, there is a large uncertainty in classifying voxels into pores or solids. An alternative approach consists in performing a ternary segmentation including a solid phase, resolved pores, and unresolved pores (Bulteys
et al., 2016). Differential imaging between a scan in a dry state and a scan in which the sample is saturated by a contrast fluid improves the identification of unresolved pores. Micro-continuum models have—by construction— the ability to simulate flows in images with SRP. In this context, flow and transport are described using Navier- Stokes-based physics in the resolved region, and are modeled using porous media concepts in unresolved pores. The choice of appropriate sub-grid models is therefore key to accurately capturing processes in the unresolved regions.

Single-phase micro-continuum models are fairly flexible to describe digital rock systems including two characteristic length scales. They are mature enough to be included in DRP commercial packages such as GeoDict. They can account for the impact of SRP on the estimation of the rock sample absolute permeability (Apourvari & Arns, 2014; Kang et al., 2019; Knackstedt et al., 2006; Singh, 2019; Soulaine et al., 2016; Xu, Guan, et al., 2022). The assessment of the sub-grid permeability, *k*, that quantifies the flow resistance of the microstructure that is not visible in the image is key. It can be probed using higher-resolution imaging or measurements. For example, Menke et al. (2022) zoomed in on the sub-voxel using nanotomography and computed the sub-voxel properties. In the absence of knowledge about the sub-voxel hydraulic properties, it has to be evaluated using the information within reach: local porosity (based on the image gray-scale), the instrument resolution, and the pore-size distribution. Kozeny-Carman permeability-porosity relationship, $k = \frac{d^2\phi^3}{180(1-\phi)^2}$, offers an appealing sub-grid model including the cell porosity, $\phi$, and the imaging instrument resolution, *d* (Abu AlSaud et al., 2020; Scheibe et al., 2015; Soulaine et al., 2016). Using such a law, Soulaine et al. (2016) observe that SRP in a Berea sandstone imaged with a 3.16 µm³/voxel resolution has a disproportionately large impact on permeability—up to a factor 2 for only 2% of voxels containing SRP if they are neglected—implying that SRP forms key percolation pathways for fluid flow. Similar results were obtained later by Cid et al. (2021) who measured deviations up to 100% in volcanic-hosted geothermal systems. Nevertheless, assigning the image resolution to *d* overestimates local permeability values. Other approaches describe the sub-voxel microstructure in each unresolved voxel as a cylinder with pore throat diameter, $d_t$. In such cases, the local permeability is described by Hagen-Poiseuille law, $k = \frac{d^2\phi}{64}$. The microtube diameter value is constrained by correlating the pore-size distribution obtained experimentally using mercury injection capillary pressure with the mapping of the local porosity field in the unresolved region (Choi et al., 2020; Yang, Kang, & Yun, 2021). Even though further work still needs to be done to reconstruct accurately the unresolved pores using the pore-size distribution and the grayscale images, micro-continuum models have demonstrated their abilities to compute absolute permeability including the instrument resolution.

The long-standing objective of DRP is to compute unsaturated porous media properties including capillary pressure curve and relative permeabilities. The two-phase micro-continuum models introduced in Section 3.1 enable DRP with sub-voxel porosity to characterize two-phase flow properties in heterogeneous rock samples. Carrillo et al. (2022), and later Lanetc et al. (2022), investigated the effects of micro-CT sub-resolution porosity on a rock's relative permeabilities, residual saturations, and fluid breakthrough times (see Figure 4). Sub-grid models for two-phase flow involve much more parameters than for single-phase flow. Indeed, in addition to the local porosity and permeability, unsaturated micro-continuum models require data for sub-voxel relative permeabilities, sub-voxel capillary pressure, as well as wettability properties at the interface between resolved and unresolved pores. If they are not properly constrained, sub-grid models including relative permeability and capillary pressure curve described by van Genuchten (1980) or Brooks and Corey (1964) models involved too many parameters. Carrillo et al. (2022) performed a sensitivity analysis on the sub-grid model parameters and showed how sample-scale relative permeabilities react to changes in the porosity, permeability, and wettability of the SRP. They found out that, even at low saturations, SRP can act as a persistent connector between otherwise-isolated fluid clusters. Multi-scale DRPs confirm that the influence of SRP cannot be disregarded without incurring significant errors in numerical predictions or experimental analyses of multiphase flow in heterogeneous porous media. The numerical implementation of the two-phase micro-continuum model used by Carrillo et al. (2022) and Lanetc et al. (2022) was not able to deal with capillary-dominated regimes. Therefore, their studies were limited to relative permeabilities and did not focus on the capillary pressure curve. Recent improvements that filtered the spurious velocities due to the computation of the fluid-fluid interface curvature should leverage this limitation (Liu et al., 2022). The a priori characterization of two-phase sub-grid models remains an open
<figure><img src="image_6.png" /><figcaption>**Figure 4.** The impact of sub-voxel on rock properties in Digital Rock Physics, adapted from Carrillo et al. (2022). (a) Full 2-D view of the IXCT scan of Estailades carbonate rock sample by Bultreys (2016), which is 7 mm in diameter with a resolution of 3.1 µm per voxel. Black is open pore space, dark gray corresponds to domains that contain SRP, and the lightest color is solid calcite. (b) 500 × 500 voxel cropped sample. (c) Corresponding segmented image. (d) Spatial distribution of the SRP (red), pore space (blue), and solid rock (transparent) within the extracted 3-D rock representation (200 × 200 × 200 cells). (e) The corresponding SRP accounts for 21% of the voxels. (f) The associated open pore space, which accounts for 40% of the voxels. (g) Sensitivity of drainage and imbibition relative permeability curves to SRP absolute permeability, from $k_0 = 10^{-17}$ to $10^{-12}$ m$^2$. (h) Sensitivity of drainage and imbibition relative permeability curves to the external wettability of rock and SRP domains, from water-wetting ($\theta_r = 30°, \theta_p = 30°$), to mixed-wetting ($\theta_r = 150°, \theta_p = 30°$ and $\theta_r = 30°, \theta_p = 150°$), to oil-wetting ($\theta_r = 150°, \theta_p = 150°$).</figcaption></figure>

challenge. Wang et al. (2022) used time-lapse X-ray microtomography to reconstruct the local capillary pressure from the dynamic water saturation profile in unresolved voxels. Their approach, developed in the context of Dual-Pore Network Modeling could help to constrain the sub-voxel capillary pressure curve in two-phase micro-continuum simulations. Two-phase micro-continuum DRP has an incredible potential to improve the understanding of the physics of flow in heterogeneous rock and to guide the development of more accurate and predictive upscaled petrophysical models.

### 4.1.2. Digital Rock Physics for Characterizing Heat Transfer

Recently, Maes and Menke (2022) modeled conjugate heat transfer during injection of cold water into a micro-CT image of Bentheimer sandstone which is analogous to a geothermal reservoir using a micro-continuum-based approach. Standard approaches for modeling heat transfer between a fluid and a solid phase usually deal with two regions meshed independently from each other. Energy equations depend on advection and diffusion in the fluid region and on diffusion only in the solid regions. Both equations are connected through coupled boundary conditions at the fluid-solid interface (continuity of fluxes and temperatures). Setting up such simulations remains challenging for complex pore morphologies. Instead, Maes and Menke (2022) used the Cartesian grid of the microtomography images and considered the solid grains as low-permeability low-porosity porous media (see Section 2.3.2) in which the temperature evolution depends only on thermal diffusion. The interfacial conditions are handled directly by the micro-continuum energy equation. They performed simulations with various fluid properties at various Reynolds numbers in order to investigate the evolution of the Nusselt numbers which describe the heat exchange between solid and fluid. Finally, they obtained different trends for the heat exchange coefficient whether the system was in the conduction-dominated regime or in the convection-dominated regime.
## 4.2. Solid Matrix Deformation

The Darcy-Brinkman-Biot model proposed by Carrillo and Bourg (2019) pushes the capabilities of micro-continuum models to describe processes involving the deformation of soft porous materials such as clays, hydrogels, and biofilms. The model is somehow a combination of DBS and the Biot poromechanics theory. The model was successfully validated against experimental data and numerical models, for example, compaction of a clay plug in a flow of water, and standard oedometric measurements of clay swelling driven by a change in salinity. It has been used to model the deformable aspect of biofilm under fluid shear stress (Kurz et al., 2022). Lately, this model has been extended to multiphase flow (Carrillo & Bourg, 2021a) opening further possibilities to describe natural porous systems. For example, the two-phase Darcy-Brinkman-Biot approach was used to investigate viscous and capillary fracturing during drainage in viscoplastic porous media that occur if viscous and/or capillary stresses are sufficient to overcome the solid's structural forces (Carrillo & Bourg, 2021b).

## 4.3. Reactive Transport Modeling

Over the last decades, reactive transport modeling has become an essential tool for the study of subsurface processes involving flow, transport and geochemical reactions (Deng et al., 2022; Steefel et al., 2005). This discipline arises from the junction of two scientific communities, namely geochemistry, and transport in porous media. They consist of computational models that describe the coupled physical, chemical, and biological processes that interact with each other over a broad range of spatial and temporal scales. The transport of chemical species and the evolution of the solid mineral volume fractions due to surface reactions is described by the equations introduced in Section 3.3. Because they have the ability to displace fluid-rock interfaces on regular grids, micro-continuum models have demonstrated great capacities to model geochemical processes at different scales of interest. In particular, the micro-continuum framework allows for the use of geochemical packages developed for continuum-scale without major modifications (Pavuluri et al., 2022; Soulaine et al., 2021): (a) they use DBS momentum instead of Darcy's law, and (b) appropriate sub-grid models for the reactive surface area per unit volume enable hybrid-scale simulations. In the following section, we review recent applications of micro-continuum-based reactive transport modeling for simulating hydro-geochemical processes.

### 4.3.1. Pore-Scale Modeling of Mineral Dissolution and Precipitation

Micro-continuum approaches with penalization—that is, low-porosity low-permeability media, see Ssection 2.3.2—offer an appealing framework to solve pore-scale reactive transport with moving fluid-solid interfaces including mineral precipitation and dissolution. At the pore-scale, the surface reactions occur only at the resolved fluid-solid interface. It is described by kinetic laws proportional to the fluid-solid surface area, $A_e$. For example, in Equation 11, we have $\dot{m}_{f,i} = A_e k_i (C_i - C_{eq})$. At the pore-scale, the local surface area, $A_e$, is an output of the simulation computed using the Volume of Solid technique (Maes et al., 2022; Soulaine et al., 2017). The fluid-solid interface is, then, moved according to the solid mass balance equation, Equation 10, with $\dot{m}_{s,i} = \beta \dot{m}_{f,i}$ where $\beta$ is a stoechiometric coefficient.

This approach has been used successfully to simulate the dissolution of a calcite crystal in a flowing acidic solution. The results are in close agreement with high-fidelity microfluidic experiments (Soulaine et al., 2017) and also with other numerical approaches (Molins et al., 2020). A couple of extensions of this model have been developed. For example, Li et al. (2023) proposed a similar approach but based on a phase-field model instead of the porosity field. Soulaine et al. (2018) model calcite dissolution including the production of CO₂ gas bubbles. The geochemical reactions used so far were very simplistic. Soulaine et al. (2021) coupled a micro-continuum solver (OpenFOAM-based) with the geochemical package PHREEQC and obtained similar results, paving, therefore, the path to simulations with multiple minerals and more comprehensive reaction networks.

Today, the pore-scale modeling of mineral precipitation with a micro-continuum approach is less mature than the simulation of dissolution processes. Yang, Stack, et al. (2021) proposed such a model but without nucleation theory. Deng et al. (2021) have nucleation theory but do not consider advection. In both cases, there is still the need to compare the results with high-resolution experimental data set in well-controlled environments. Reactive microfluidics combined with Raman spectroscopy could be an option to acquire such reference data (Poonoosamy et al., 2020).
4.3.2. Fracture-Matrix Interactions

Fractures are key flow pathways that dominate fluids migration and solute transport in many geological for-
mations (Berkowitz, 2002). Experimental studies show that geochemical reactions may enlarge or seal fractures
depending on a number of factors including the flow regime, the composition of fluids, the mass transfer between
the fracture and the host porous matrix, and the alteration of mechanical properties (Fitts & Peters, 2013).
Moreover, the dissolution of the minerals that surround the fracture leads to the development of an altered layer
whose porosity varies with mineral composition, flow rate, and fluid composition (Noiriel et al., 2013). As a result
of the evolution of the fracture morphology, the hydraulic properties are subject to modifications. This feedback
between fluid flow and geochemical reactions can enhance the flow conductivity for deep geothermal energy, or
inversely, can threaten the integrity of storage systems. For example, in the context of CO2 sequestration in the
subsurface, the acidizing of the system can reorganize the pore-space and produce leakage pathways to the surface
(DePaolo & Cole, 2013).

Micro-continuum approaches are powerful tools to shed new light on fracture-matrix interactions. Indeed, hybrid-scale micro-continuum models solve Navier-Stokes-based physics within the fracture and Darcy-based flow and transport in the porous matrix. Such models are able to capture the development of the weathered zone at the vicinity of the fracture using appropriate surface area model within the matrix (e.g., A_e is computed using a power-law) (Noiriel & Soulaine, 2021). Zhang et al. (2022) coupled a micro-continuum framework with the geochemical modeling capability of CrunchFlow to investigate the development of the altered porous layer resulting from the preferential dissolution of the fast-reacting minerals and their subsequent chemical-physical alteration of the fractures. It will be interesting to use such approaches to study the thickness of the altered layer according to the predominance of diffusion, advection, and surface reactions described by the Péclet and Damköhler numbers. Hu et al. (2023) and Luo et al. (2012) used reactive micro-continuum models including the Boussinesq approximation to investigate the role of gravity in dissolving horizontal fractures. They were able to obtain scallop patterns that are characteristic of many subsurface systems. The role of the mass transfer in the altered layer in the emergence of scallop patterns, however, is still not clear, and more in-depth analysis is needed. Laboratory experiments using analog materials such as the caramel blocks used by Cohen et al. (2020) will help to verify the robustness and accuracy of the approach.

Not only fractures may be altered geochemically, but they also deform under stress. Pyrak-Nolte and Nolte (2016) demonstrated that a relationship exists between flow and stiffness. Coupled flow-geomechanical simulations of fracture-matrix systems are key to understand and characterize better this relationship. He et al. (2022) investigate the stress influence on fluid leakage caused by matrix-fracture interaction using loosely-coupled flow-geomechanical model in which the fluid flow is described by the DBS equation and a Brown-Scholz model to describe the rock deformation based on the fracture aperture. A micro-continuum approach based on the Darcy-Brinkman-Biot model (see Section 3.4) will enable the investigation of flow-geomechanics processes in fractured rocks using a fully-coupled model and provides more insights not only into the complex feedback between flow and mechanics in fracture-matrix systems but also with chemical alterations.

4.3.3. Wormholes Formation in Acidic Environments

The micro-continuum approach is particularly well-suited to simulate the formation and growth of wormholes without involving complex re-meshing strategies (Golfier et al., 2002; Ormond & Ortoleva, 2000). The development of wormholes occurs during acid stimulation. that is, an engineering process to increase the productivity of clogged wells or to increase the conductivity in deep geothermal reservoirs. The modeling of such structures is particularly challenging as they are intrinsically multi-scale and dynamic. Similar to the hybrid-scale micro-continuum model used to simulate fracture-matrix interactions, the evolution of the porous matrix is governed by the transport equations for chemical species and solid mineral volume fractions described in Section 3.3. If the matrix is fully dissolved within a cell, the flow micro-continuum equation tends to Navier-Stokes.

Micro-continuum models capture the different dissolution patterns involved in reactive transport, namely: facial dissolution, conical dissolution, one-dominant wormhole, ramified wormholes, and uniform dissolution (Golfier et al., 2002). An example of the emergence of wormholes obtained by a reactive micro-continuum model is shown in Figure 5. In this example, a first-order kinetic combined with A_e modeled as a power-law is used to simulate mineral dissolution. The dissolution regimes result from the complex interplay between advection, diffusion, and surface reaction. Although the transition between each regime is still a matter of controversy, micro-continuum
<figure><img src="image_7.png" /><figcaption>**Figure 5.** Simulation results of dissolution instabilities in a carbonate core. At early times, instabilities begin to develop when the reactant preferentially penetrates the regions with the biggest pores (the cells with the highest porosity and permeability) and erodes the walls to form the flow channels. When a finger is slightly longer than its surrounding neighbors, the pressure drop in this finger is smaller than in the surrounding fingers, which favors a preferential flow through the finger and reduces the growth of the other surrounding smaller fingers. The process continues until only one dominant wormhole takes over and breakthroughs. Time is normalized by the breakthrough time.</figcaption></figure>

models are in qualitatively good agreement with experiments (Daccord & Lenormand, 1987; Fredd & Fogler, 1998). Comparison, however, always relies on post-mortem analysis of the experiment. The recent advances in time-lapse microtomography will provide three-dimensional images of the evolution of dissolution patterns (Ott & Oedai, 2015). Reactive microfluidics is also a promising technique to verify the robustness of the model (Osselin et al., 2016; Rembert et al., 2023; Song et al., 2014).

Future research directions include the simulation of wormholes in multi-mineral rocks with comprehensive geochemistry, wormholing in unsaturated porous media, and coupling with geomechanics.

### 4.3.4. Modeling of Biofilm Growth

Micro-continuum models have an interesting potential for modeling biofilm growth and the subsequent pore-clogging in pore-scale micro-structures (Yan et al., 2017). In such cases, the biofilm is seen as a deformable solid phase that grows and moves with respect to the supply of nutrients (Tian & Wang, 2019). The development of micro-continuum models for biofilm growth is still in its infancy. Kurz et al. (2022) demonstrate using microfluidics and micro-continuum simulations that microbial growth and its competition with shear stresses lead to intermittent pathways in bioclogged porous media. In their simulations, the biofilm was described as a viscoplastic porous continuum through a Darcy-Brinkman-Biot model (Carrillo & Bourg, 2019). The biofilm was growing by the supply of nutrients and was allowed to be transported within the microfluidic device and to leave through the outlet. This model is the first of its kind. Further development will include a more comprehensive multispecies reactive transport, a more realistic biological model as well as multiphase flow to describe the growth of biofilm in the vadose zone.

### 4.3.5. Microbially-Induced Calcite Precipitation

Among the scientific literature on micro-continuum models for hydro-geochemical processes, the work of Minto and co-workers is of particular interest (Minto et al., 2018, 2019). They report leading-edge simulations of
<figure><img src="image_8.png" /><figcaption>**Figure 6.** Microbially-induced calcite precipitation from Minto et al. (2019). (a) Greyscale X-ray CT data from a 6-mm-diameter core subsampled from a larger 40-cm-diameter block of MICP-treated beach sand. The data have been segmented based on X-ray attenuation to distinguish sand grains from CaCO₃. (b) Modeled velocity in a section (flow is from left to right) and (c) resulting velocity-dependent bacterial attachment coefficient. (d) CaCO₃ precipitation in model with velocity-dependent bacterial attachment, and, for comparison, (e) a model without.</figcaption></figure>

microbially-induced calcite precipitation. Microbially induced calcium carbonate precipitation (MICP) is a promising technique that could be used for soil stabilization, permeability control in porous and fractured media, for sealing leaky hydrocarbon wells, and for immobilizing contaminants (Phillips et al., 2016).

They investigate optimum MICP treatment strategies by developing a micro-continuum model that includes bacteria transport and attachment, urea hydrolysis, tractable CaCO₃ precipitation, and modification to the porous media in terms of porosity and permeability. The transport of bacteria is modeled by an advection-dispersion-reaction equation that describes the evolution of the bacteria concentration including their attachment/detachment to the mineral surfaces. As other micro-continuum models, their approach describes processes both at the pore-scale using image-based simulations (see Figure 6) and at the field-scale. Their model—conceptually validated against experimental data—n narrows down the range of possible injection strategies to determine the most promising for implementation at a field scale. Their results indicate that phased injection strategies may lead to the most uniform precipitation in a porous medium.

Importantly, the authors emphasize the importance of an accurate description of bacterial attachment processes— at least the attachment rate should be a function of fluid velocity. Today, the modeling of attachment/detachment of all kinds of colloidal particles (nanometals, fine particles, bacteria, viruses, asphaltenes...) in the subsurface remains a challenge that needs to be investigated further.

**5. Conclusion**

Micro-continuum models are a versatile and powerful approach to simulate multi-scale coupled processes in porous media. The governing equations are rooted in the elementary physical principles and combined with appropriate sub-grid models for describing processes in unresolved porosity. A key interest of the method is its ability to model pore-scale and field-scale processes within the same framework. State-of-the-art micro-continuum models handle two-phase flow, hydro-bio-geochemical processes, and poromechanics.

Despite its early age micro-continuum approach for pore-scale processes has already demonstrated its strength in image-based simulations and coupled physics. The growing interest in the micro-continuum approach is emphasized by the increasing number of scientific publications using such a concept.

Micro-continuum models are still in their infancy, and active research in the area intends to improve their robustness and extend their capabilities. Open challenges and current lines of development include numerical efficiency, advanced sub-grid models, and additional coupled physics. In this regard, an important challenge consists in the derivation of sound resolved-unresolved interfacial conditions for scalar transport and multiphase
flow. Efforts should be dedicated to the development of coupled multiphysics micro-continuum models including the transport of clouds of inorganic colloidal and nano-particles, bacteria, and viruses, as well as the consideration of electrical potential gradients and electrical currents. Accurate modeling of colloidal particles attachment-detachment at the solid walls within the micro-continuum framework remains an open question. Importantly, hybrid-scale micro-continuum models need to be verified and benchmarked using experimental reference data.

## Data Availability Statement

The result presented in Figure 5 is obtained with the open-source package porousMedia4Foam (Soulaine, 2023). The case setup is located in the folder: tutorials-pM4F/dbsFoam/simpleFirstOrderKineticMole/dissolutionWormholing of the Zenodo repository.

## Acknowledgments

The author has been funded by the European Union (ERC, COCONUT, 101043288).

## References

Abu AlSaud, M., Gmira, A., Al-Enezi, S., & Yousef, A. (2020). Pore-scale simulation of fluid flow in carbonates using micro-CT scan images. In *International petroleum technology conference*. IPTC-19832-MS. Society of Petroleum Engineers. https://doi.org/10.2523/iptc-19832-ms

Ahmadi, N., Heck, K., Rolle, M., Helmig, R., & Mosthaf, K. (2021). On multicomponent gas diffusion and coupling concepts for porous media and free flow: A benchmark study. *Computational Geosciences*, 25(5), 1493-1507. https://doi.org/10.1007/s10596-021-10057-y

Alazmi, B., & Vafai, K. (2001). Analysis of fluid flow and heat transfer interfacial conditions between a porous medium and a fluid layer. *International Journal of Heat and Mass Transfer*, 44(9), 1735-1749. https://doi.org/10.1016/s0017-9310(00)00217-9

Andrä, H., Combaret, N., Dvorkin, J., Glatt, E., Han, J., Kabel, M., et al. (2013a). Digital rock physics benchmarks Part I: Imaging and segmentation. *Computers & Geosciences*, 50, 25-32. https://doi.org/10.1016/j.cageo.2012.09.005

Andrä, H., Combaret, N., Dvorkin, J., Glatt, E., Han, J., Kabel, M., et al. (2013b). Digital rock physics benchmarks Part II: Computing effective properties. *Computers & Geosciences*, 50, 33-43. https://doi.org/10.1016/j.cageo.2012.09.008

Angot, P., Bruneau, C.-H., & Fabrie, P. (1999). A penalization method to take into account obstacles in incompressible viscous flows. *Numerische Mathematik*, 81(4), 497-520. https://doi.org/10.1007/s002110050401

Apourvari, S. N., & Arns, C. H. (2014). An assessment of the influence of micro-porosity for effective permeability using local flux analysis on tomographic images. In *International petroleum technology conference, Doha, Qatar*.

Arns, C., Bauget, F., Limaye, A., Sakellariou, A., Senden, T., Sheppard, A., et al. (2005). Pore-scale characterization of carbonates using X-ray microtomography. *SPE Journal*, 10(4), 475-484. https://doi.org/10.2118/90368-pa

Auriault, J.-L. (2009). On the domain of validity of Brinkman's equation. *Transport in Porous Media*, 79(2), 215-223. https://doi.org/10.1007/s11242-008-9308-7

Balhoff, M. T., Thompson, K. E., & Hjortsø, M. (2007). Coupling pore-scale networks to continuum-scale models of porous media. *Computers & Geosciences*, 33(3), 393-410. https://doi.org/10.1016/j.cageo.2006.05.012

Beavers, G. S., & Joseph, D. D. (1967). Boundary conditions at a naturally permeable wall. *Journal of Fluid Mechanics*, 30(1), 197-207. https://doi.org/10.1017/S0022112067001375

Berkowitz, B. (2002). Characterizing flow and transport in fractured geological media: A review. *Advances in Water Resources*, 25(8-12), 861-884. https://doi.org/10.1016/S0309-1708(02)00042-8

Blunt, M. J., Bijeljic, B., Dong, H., Gharbi, O., Iglauer, S., Mostaghimi, P., et al. (2013). Pore-scale imaging and modelling. *Advances in Water Resources*, 51, 197-216. https://doi.org/10.1016/j.adwatres.2012.03.003

Brackbill, J. U., Kothe, D. B., & Zemach, C. (1992). A continuum method for modeling surface tension. *Journal of Computational Physics*, 100(2), 335-354. https://doi.org/10.1016/0021-9991(92)902040-y

Brinkman, H. C. (1947). A calculation of the viscous force exerted by a flowing fluid on a dense swarm of particles. *Applied Scientific Research, AI*(1), 27-34. https://doi.org/10.1007/bf02120313

Brooks, R., & Corey, A. (1964). Hydraulic properties of porous media. In *Hydrology papers* (p. 37). Colorado State University.

Bultreys, T. (2016). Estaillades carbonate 2. Retrieved from http://www.digitalrocksportal.org/projects/58

Bultreys, T., Boever, W. D., & Cnudde, V. (2016). Imaging and image-based fluid transport modeling at the pore scale in geological materials: A practical introduction to the current state-of-the-art. *Earth-Science Reviews*, 155, 93-128. https://doi.org/10.1016/j.earscirev.2016.02.001

Bultreys, T., Van Hoorebeke, L., & Cnudde, V. (2015). Multi-scale, micro-computed tomography-based pore network models to simulate drainage in heterogeneous rocks. *Advances in Water Resources*, 78, 36-49. https://doi.org/10.1016/j.adwatres.2015.02.003

Carman, P. C. (1937). Fluid flow through granular beds. *Transactions of the Institution of Chemical Engineers*, 15, 150-166.

Carrillo, F. J., & Bourg, I. (2021a). Modeling multiphase flow within and around deformable porous materials: A Darcy-Brinkman-Biot approach. *Water Resources Research*, 52(2), e2028WR028734. https://doi.org/10.1029/2028WR028734

Carrillo, F. J., & Bourg, I. C. (2019). A Darcy-Brinkman-Biot approach to modeling the hydrology and mechanics of porous media containing macropores and deformable microporous regions. *Water Resources Research*, 55(10), 8096-8121. https://doi.org/10.1029/2019wr024712

Carrillo, F. J., & Bourg, I. C. (2021b). Capillary and viscous fracturing during drainage in porous media. *Physical Review E*, 103(6), 063106. https://doi.org/10.1103/physreve.103.063106

Carrillo, F. J., Bourg, I. C., & Soulaine, C. (2020). Multiphase flow modeling in multiscale porous media: An open-source micro-continuum approach. *Journal of Computational Physics*, X(8), 100073. https://doi.org/10.1016/j.jcpx.2020.100073

Carrillo, F. J., Soulaine, C., & Bourg, I. C. (2022). The impact of sub-resolution porosity on numerical simulations of multiphase flow. *Advances in Water Resources*, 161, 104094. https://doi.org/10.1016/j.adwatres.2021.104094

Cassie, A., & Baxter, S. (1944). Wettability of porous surfaces. *Transactions of the Faraday Society*, 40, 546-551. https://doi.org/10.1039/tf9444000546

Choi, C.-S., Lee, Y.-K., & Song, J.-J. (2020). Equivalent pore channel model for fluid flow in rock based on microscale X-ray CT imaging. *Materials*, 13(11), 2619. https://doi.org/10.3390/ma13112619

Cid, H. E., Carrasco-Núñez, G., Manea, V. C., Vega, S., & Castaño, V. (2021). The role of microporosity on the permeability of volcanic-hosted geothermal reservoirs: A case study from Los Humeros, Mexico. In *Geothermics* (Vol. 90). https://doi.org/10.1016/j.geothermics.2020.102020.102020.
Cohen, C., Berhanu, M., Derr, J., & du Pont, S. C. (2020). Buoyancy-driven dissolution of inclined blocks: Erosion rate and pattern formation. *Physical Review Fluids, 5*, 053802. https://doi.org/10.1103/physrevfluids.5.053802

Daccord, G., & Lenormand, R. (1987). Fractal patterns from chemical dissolution. *Nature, 325*(6099), 41–43. https://doi.org/10.1038/325041a0

Deng, H., Gharasoo, M., Zhang, L., Dai, Z., Hajizadeh, A., Peters, C. A., et al. (2022). A perspective on applied geochemistry in porous media: Reactive transport modeling of geochemical dynamics and the interplay with flow phenomena and physical alteration. *Applied Geochemistry, 146*, 105445. https://doi.org/10.1016/j.apgeochem.2022.105445

Deng, H., Tournassat, C., Molins, S., Claret, F., & Steefel, C. I. (2021). A pore-scale investigation of mineral precipitation driven diffusivity change at the column-scale. *Water Resources Research, 57*(5), e2020WR028483. https://doi.org/10.1029/2020wr028483

DePaolo, D. J., & Cole, D. R. (2013). Geochemistry of geologic carbon sequestration: An overview. *Reviews in Mineralogy and Geochemistry, 77*(1), 1–14. https://doi.org/10.2138/rmg.2013.77.1

Fitts, J. P., & Peters, C. A. (2013). Caprock fracture dissolution and CO₂ leakage. *Reviews in Mineralogy and Geochemistry, 77*(1), 459–479. https://doi.org/10.2138/rmg.2013.77.13

Forchheimer, P. (1901). Wasserbewegung durch Boden. *Z. Ver. Deutsch. Ing., 45*, 1782–1788.

Fredd, C. N., & Fogler, H. S. (1998). Influence of transport and reaction on wormhole formation in porous media. *AIChE Journal, 44*(9), 1933–1949. https://doi.org/10.1002/aic.690440902

Golfier, F., Zarcone, C., Bazin, B., Lenormand, R., Lasseux, D., & Quintard, M. (2002). On the ability of a Darcy-scale model to capture wormhole formation during the dissolution of a porous medium. *Journal of Fluid Mechanics, 457*, 213–254. https://doi.org/10.1017/s0022112002007735

He, X., AlSinan, M., Zhang, Z., Kwak, H., Hoteit, H., & Oct (2022). Micro-continuum approach for modeling coupled flow and geomechanical processes in fractured rocks. In *SPE annual technical conference and exhibition. SPE-210453-MS*. Society of Petroleum Engineers. https://doi.org/10.2118/210453-ms

Hirt, C. W., & Nichols, B. D. (1981). Volume-of-fluid (VOF) method for the dynamic of free boundaries. *Journal of Computational Physics, 39*(1), 201–225. https://doi.org/10.1016/0021-9991(81)90145-5

Horgue, P., Prat, M., & Quintard, M. (2014). A penalization technique applied to the volume-of-fluid method: Wettability condition on immersed boundaries. *Computers & Fluids, 100*, 255–266. https://doi.org/10.1016/j.compfluid.2014.05.027

Hu, R., Li, K., Zhou, C.-X., Wang, T., Yang, Z., & Chen, Y.-F. (2023). On the role of gravity in dissolving horizontal fractures. *Journal of Geophysical Research: Solid Earth, 128*(3), e2022JB025214. https://doi.org/10.1029/2022jb025214

Kang, D. H., Yang, E., & Yun, T. S. (2019). Stokes-Brinkman flow simulation based on 3-D-CT images of porous rock using grayscale pore voxel permeability. *Water Resources Research, 55*(5), 4448–4464. https://doi.org/10.1029/2018wr024179

Khadra, K., Angot, P., Parneix, S., & Caltagirone, J.-P. (2000). Fictitious domain approach for numerical modelling of Navier-Stokes equations. *International Journal for Numerical Methods in Fluids, 34*(8), 651–684. https://doi.org/10.1002/1097-0363(20001230)34:8<651::AID-FLD61>3.0.CO;2-D

Klinkenberg, L. (1941). The permeability of porous media to liquids and gases. In *Drilling and production practice*. American Petroleum Institute.

Knackstedt, M. A., Arns, C. H., Ghous, A., Sakellariou, A., Senden, T. J., Sheppard, A. P., et al. (2006). 3D imaging and characterization of the pore space of carbonate core; implications to single and two phase flow properties. In *SPWLA 47th annual logging symposium*. Society of Petrophysicists and Well-Log Analysts.

Kozeny, J. (1927). *Über kapillare leitung der wasser in boden* (Vol. 136, pp. 271–306). Royal Academy of Science. Proc. Class I.

Kurz, D. L., Secchi, E., Carrillo, F. J., Bourg, I. C., Stocker, R., & Jimenez-Martinez, J. (2022). Competition between growth and shear stress drives intermittency in preferential flow paths in porous medium biofilms. *Proceedings of the National Academy of Sciences of the United States of America, 119*(30). https://doi.org/10.1073/pnas.2122202119

Lanetc, Z., Zhuravljov, A., Jing, Y., Armstrong, R. T., & Mostaghimi, P. (2021). Coupling of transient matrix diffusion and pore network models for gas flow in coal. *Journal of Natural Gas Science and Engineering, 88*, 103741. https://doi.org/10.1016/j.jngse.2020.103741

Lanetc, Z., Zhuravljov, A., Shapoval, A., Armstrong, R. T., & Mostaghimi, P. (2022). Inclusion of microporosity in numerical simulation of relative permeability curves. In *International petroleum technology conference. IPTC-21975-MS*. Society of Petroleum Engineers. https://doi.org/10.2523/iptc-21975-ms

Leal, A. (2015). Reaktoro: An open-source unified framework for modeling chemically reactive systems. Retrieved from https://reaktoro.org

Li, H., Wang, F., Wang, Y., Yuan, Y., Feng, G., Tian, H., & Xu, T. (2023). Phase-field modeling of coupled reactive transport and pore structure evolution due to mineral dissolution in porous media. *Journal of Hydrology, 619*, 129363. https://doi.org/10.1016/j.jhydrol.2023.129363

Lisitsa, V., Khachkova, T., Krutko, V., & Avdonin, A. (2023). Simulation of two-phase flow in models with micro-porous material. In (Eds), O. Gervasi, B. Murgante, A. M. A. C. Rocha, C. Garau, F. Scorza, Y. Karaca, et al., *Computational science and its applications – ICCSA 2023 Workshops* (pp. 3–18). Springer Nature Switzerland.

Liu, Z., Yang, J., Xu, Q., & Shi, L. (2022). Improved micro-continuum approach for capillary-dominated multiphase flow with reduced spurious velocity. *Physics of Fluids, 34*(12), 122108. https://doi.org/10.1063/5.0127603

Luo, H., Quintard, M., Debenest, G., & Laouafa, F. (2012). Properties of a diffuse interface model based on a porous medium theory for solid-liquid dissolution problems. *Computational Geosciences, 16*(4), 913–932. https://doi.org/10.1007/s10596-012-9295-1

Maes, J., & Menke, H. P. (2022). GeoChemFoam: Direct modelling of flow and heat transfer in micro-CT images of porous media. *Heat and Mass Transfer, 58*(11), 1937–1947. https://doi.org/10.1007/s00231-022-03221-2

Maes, J., Soulaine, C., Menke, H. P., & July (2022). Improved volume-of-solid formulations for micro-continuum simulation of mineral dissolution at the pore-scale. *Frontiers in Earth Science, 10*. https://doi.org/10.3389/feart.2022.917931

Meakin, P., & Tartakovsky, A. M. (2009). Modeling and simulation of pore-scale multiphase fluid flow and reactive transport in fractured and porous media. *Reviews of Geophysics, 47*(RG3002), 1–47. https://doi.org/10.1029/2008rg000263

Menke, H. P., Gao, Y., Linden, S., Andrew, M. G., & July (2022). Using nano-XRM and high-contrast imaging to inform micro-porosity permeability during Stokes-Brinkman single and two-phase flow simulations on micro-CT images. *Frontiers in Water, 4*. https://doi.org/10.3389/frwa.2022.935035

Minto, J. M., Lunn, R. J., & Mountassir, G. E. (2019). Development of a reactive transport model for field-scale simulation of microbially induced carbonate precipitation. *Water Resources Research, 55*(8), 7229–7245. https://doi.org/10.1029/2019wr025153

Minto, J. M., Tan, Q., Lunn, R. J., Mountassir, G. E., Guo, H., & Cheng, X. (2018). Microbial mortar’-restoration of degraded marble structures with microbially induced carbonate precipitation. *Construction and Building Materials, 180*, 44–54. https://doi.org/10.1016/j.conbuildmat.2018.05.200

Molins, S., Soulaine, C., Prasianakis, N., Abbasi, A., Poncet, P., Ladd, A., et al. (2020). Simulation of mineral dissolution at the pore scale with evolving fluid-solid interfaces: Review of approaches and benchmark problem set. *Computational Geosciences, 25*(4), 1–34. https://doi.org/10.1007/s10596-019-09903-x
Molins, S., Trebotich, D., Arora, B., Steefel, C. I., & Deng, H. (2019). Multi-scale model of reactive transport in fractured media: Diffusion limitations on rates. *Transport in Porous Media, 128*(2), 701–721. https://doi.org/10.1007/s11242-019-01266-2

Muskat, M. (1949). *Physical principles of oil production*. McGraw-Hill.

Neale, G., & Nader, W. (1974). Practical significance of Brinkman's extension of Darcy's law: Coupled parallel flows within a channel and a bounding porous medium. *Canadian Journal of Chemical Engineering*, 52(4), 475–478. https://doi.org/10.1002/cjce.5450520407

Noiriel, C., Gouze, P., & Madé, B. (2013). 3D analysis of geometry and flow changes in a limestone fracture during dissolution. *Journal of Hydrology*, 486, 211–223. https://doi.org/10.1016/j.jhydrol.2013.01.035

Noiriel, C., Luquot, L., Madé, B., Raimbault, L., Gouze, P., & van der Lee, J. (2009). Changes in reactive surface area during limestone dissolution: An experimental and modelling study. *Chemical Geology*, 265(1–2), 160–170. [CO2] geological storage: Integrating geochemical, hydrodynamical, mechanical and biological processes from the pore to the reservoir scale. https://doi.org/10.1016/j.chemgeo.2009.01.032

Noiriel, C., & Soulaine, C. (2021). Pore-scale imaging and modelling of reactive flow in evolving porous media: Tracking the dynamics of the fluid-rock interface. In *Transport in porous media*. https://doi.org/10.1007/s11242-021-01613-2

Ormond, A., & Ortoleva, P. (2000). Numerical modeling of reaction-induced cavities in a porous rock. *Journal of Geophysical Research*, 105(B7), 16737–16747. https://doi.org/10.1029/2000jb900116

Osselin, F., Kondratiuk, P., Budek, A., Cybulski, O., Garstecki, P., & Szymczak, P. (2016). Microfluidic observation of the onset of reactive-infiltration instability in an analog fracture. *Reactive-infiltration instability. Geophysical Research Letters*, 43(13), 6907–6915. https://doi.org/10.1002/2016gl069261

Ott, H., & Oedai, S. (2015). Wormhole formation and compact dissolution in single-and two-phase CO-brine injections. *Geophysical Research Letters*, 42(7), 2270–2276. https://doi.org/10.1002/2015GL063582

Parkhurst, D. L., & Wissmeier, L. (2015). PhreeqcRM: A reaction module for transport simulators based on the geochemical model PHREEQC. *Advances in Water Resources*, 83, 176–189. https://doi.org/10.1016/j.adwatres.2015.06.001

Pavuluri, S., Tournassat, C., Claret, F., & Soulaine, C. (2022). Reactive transport modelling with a coupled OpenFOAM®-PHREEQC platform. *Transport in Porous Media*, 145, 475–504. https://doi.org/10.7185/gold2021.4620

Perez, S., Moonen, P., & Poncet, P. (2022). On the deviation of computed permeability induced by unresolved morphological features of the pore space. *Transport in Porous Media*, 141(1), 151–184. https://doi.org/10.1007/s11242-021-01713-z

Phillips, A. J., Cunningham, A. B., Gerlach, R., Hiebert, R., Hwang, C., Lomans, B. P., et al. (2016). Fracture sealing with microbially-induced calcium carbonate precipitation: A field study. *Environmental science & technology*, 50(7), 4111–4117. https://doi.org/10.1021/acs.est.5b05559

Poonoosamy, J., Soulaine, C., Burmeister, A., Deissmann, G., Bosbach, D., & Roman, S. (2020). Microfluidic flow-through reactor and 3D Raman imaging for in situ assessment of mineral reactivity in porous and fractured porous media. In *Lab-on-a-chip*. https://doi.org/10.1039/d0lc00360c

Pyrak-Nolte, L. J., & Nolte, D. D. (2016). Approaching a universal scaling relationship between fracture stiffness and fluid flow. *Nature Communications*, 7(1), 10663. https://doi.org/10.1038/ncomms10663

Rembert, F., Stolz, A., Soulaine, C., & Roman, S. (2023). A microfluidic chip for geoelectrical monitoring of critical zone processes. *Lab on a Chip*, 23(15), 3433–3442. https://doi.org/10.1039/d3lc00377a

Rusche, H. (2003). *Computational fluid dynamics of dispersed two-phase flows at high phase fractions*. PhD thesis. Imperial College London (University of London).

Sadeghnejad, S., Enzmann, F., & Kersten, M. (2021). Digital rock physics, chemistry, and biology: Challenges and prospects of pore-scale modelling approach. *Applied Geochemistry*, 131, 105028. https://doi.org/10.1016/j.apgeochem.2021.105028

Scheibe, T. D., Perkins, W. A., Richmond, M. C., McKinley, M. I., Romero-Gomez, P. D. J., Oostrom, M., et al. (2015). Pore-scale and multiscale numerical simulation of flow and transport in a laboratory-scale column. *Water Resources Research*, 51(2), 1023–1035. https://doi.org/10.1002/2014WR015959

Sharaborin, E. L., Rogozin, O. A., & Kasimov, A. R. (2021). The coupled volume of fluid and brinkman penalization methods for simulation of incompressible multiphase flows. *Fluids*, 6(9), 334. https://doi.org/10.3390/fluids6090334

Singh, K. (2019). How hydraulic properties of organic matter control effective liquid permeability of mudrocks. *Transport in Porous Media*, 129(3), 761–777. https://doi.org/10.1007/s11242-019-01305-y

Song, W., de Haas, T. W., Fadaei, H., & Sinton, D. (2014). Chip-off-the-old-rock: The study of reservoir-relevant geological processes with real-rock micromodels. *Lab on a Chip*, 14(22), 4382–4390. https://doi.org/10.1039/C4LC00608A

Soulaine (2023). porousMedia4Foam: December 2, 2023 release (version v11legacy) [Software]. Zenodo. https://doi.org/10.5281/zenodo.10233603

Soulaine, C., Creux, P., & Tchelepi, H. A. (2019). Micro-continuum framework for pore-scale multiphase fluid transport in shale formations. In *Transport in porous media*.

Soulaine, C., Gjetvaj, F., Garing, C., Roman, S., Russian, A., Gouze, P., et al. (2016). The impact of sub-resolution porosity of X-ray tomography images on the permeability. *Transport in Porous Media*, 113(1), 227–243. https://doi.org/10.1007/s11242-016-0690-2

Soulaine, C., Pavuluri, S., Claret, F., & Tournassat, C. (2021). porousMedia4Foam: Multi-scale open-source platform for hydro-geochemical simulations with OpenFOAM. *Environmental Modelling & Software*, 145, 105199. https://doi.org/10.1016/j.envsoft.2021.105199

Soulaine, C., & Quintard, M. (2014). On the use of a Darcy–Forchheimer like model for a macro-scale description of turbulence in porous media and its application to structured packings. *International Journal of Heat and Mass Transfer*, 74(0), 88–100. https://doi.org/10.1016/j.ijheatmasstransfer.2014.02.069

Soulaine, C., Roman, S., Kovscek, A., & Tchelepi, H. A. (2017). Mineral dissolution and wormholing from a pore-scale perspective. *Journal of Fluid Mechanics*, 827, 457–483. https://doi.org/10.1017/jfm.2017.499

Soulaine, C., Roman, S., Kovscek, A., & Tchelepi, H. A. (2018). Pore-scale modelling of multiphase reactive flow. Application to mineral dissolution with production of CO. *Journal of Fluid Mechanics*, 855, 616–645. https://doi.org/10.1017/jfm.2018.655

Soulaine, C., & Tchelepi, H. A. (2016). Micro-continuum approach for pore-scale simulation of subsurface processes. *Transport in Porous Media*, 113(3), 431–456. https://doi.org/10.1007/s11242-016-0701-3

Steefel, C. I. (2009). *CrunchFlow: Software for modeling multicomponent reactive flow and transport. user's manual*. Earth Sciences Division Lawrence Berkeley National Laboratory Berkeley.

Steefel, C. I., Beckingham, L. E., & Landrot, G. (2015). Micro-continuum approaches for modeling pore-scale geochemical processes. *Reviews in Mineralogy and Geochemistry*, 80(1), 217–246. https://doi.org/10.2138/rmg.2015.80.07

Steefel, C. I., DePaolo, D. J., & Lichtner, P. C. (2005). Reactive transport modeling: An essential tool and a new research approach for the Earth sciences. *Earth and Planetary Science Letters*, 240(3–4), 539–558. https://doi.org/10.1016/j.epsl.2005.09.017
Tam, C. K. (1969). The drag on a cloud of spherical particles in low Reynolds number flow. *Journal of Fluid Mechanics*, **38**(03), 537–546. https://doi.org/10.1017/s0022112069000322

Tian, Z., & Wang, J. (2019). Lattice Boltzmann simulation of biofilm clogging and chemical oxygen demand removal in porous media. *AIChE Journal*, **65**(9). https://doi.org/10.1002/aic.16661

van Genuchten, M. T. (1980). A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. *Soil Science Society of America Journal*, **44**(5), 892–898. https://doi.org/10.2136/sssaj1980.03615995004400050002x

Veyskarami, M., Michalkowski, C., Bringedal, C., & Helmig, R. (2023). Droplet formation, growth and detachment at the interface of a coupled free-flow-porous medium system: A new model development and comparison. In *Transport in porous media*. https://doi.org/10.1007/s11242-023-01944-2

Voller, V. R. (2009). Numerical methods for phase-change problems. In *Handbook of numerical heat transfer* (pp. 593–622). John Wiley & Sons, Inc. https://doi.org/10.1002/9780470172599.ch19

Voller, V. R., Cross, M., & Markatos, N. C. (1987). An enthalpy method for convection/diffusion phase change. *International Journal for Numerical Methods in Engineering*, **24**(1), 271–284. https://doi.org/10.1002/nme.1620240119

Wang, S., Ruspini, L. C., Øren, P.-E., Offenwert, S. V., & Bultreys, T. (2022). Anchoring multi-scale models to micron-scale imaging of multiphase flow in rocks. *Water Resources Research*, **58**(1), e2021WR030870. https://doi.org/10.1029/2021wr030870

Weishaupt, K., Joekar-Niasar, V., & Helmig, R. (2019). An efficient coupling of free flow and porous media flow using the pore-network modeling approach. *Journal of Computational Physics*, **X**(1), 100011. https://doi.org/10.1016/j.jcpx.2019.100011

Weishaupt, K., Koch, T., & Helmig, R. (2022). A fully implicit coupled pore-network/free-flow model for the pore-scale simulation of drying processes. *Drying Technology*, **40**(4), 697–718. https://doi.org/10.1080/07373937.2021.1955706

Wenzel, R. N. (1936). Resistance of solid surfaces to wetting by water. *Industrial and Engineering Chemistry*, **28**(8), 988–994. https://doi.org/10.1021/ie50320a024

Whitaker, S. (1986a). Flow in porous media I: A theoretical derivation of Darcy's law. *Transport in Porous Media*, **1**(1), 3–25. https://doi.org/10.1007/BF01036523

Whitaker, S. (1986b). Flow in porous media II: The governing equations for immiscible, two-phase flow. In *Transport in porous media* (Vol. 1(2), 105–125).

Xu, Q., Dai, X., Yang, J., Liu, Z., & Shi, L. (2022). Image-based modelling of coke combustion in a multiscale porous medium using a micro-continuum framework. *Journal of Fluid Mechanics*, **932**, A51. https://doi.org/10.1017/jfm.2021.1039

Xu, Q., Guan, L., Zhang, W., Shi, L., Shao, H., Wang, G., & Long, W. (2022). Multiscale digital rock imaging and modeling for measuring the heterogeneous carbonate and conglomerate permeability at the laboratory plug scale. *Energy and Fuels*, **36**(18), 11025–11039. https://doi.org/10.1021/acs.energyfuels.2c02179

Yan, Z., Liu, C., Liu, Y., & Bailey, V. L. (2017). Multiscale investigation on biofilm distribution and its impact on macroscopic biogeochemical reaction rates. *Water Resources Research*, **53**(11), 8698–8714. https://doi.org/10.1002/2017wr020570

Yang, E., Kang, D. H., & Yun, T. S. (2021). Ternary segmentation and estimation of permeability for porous rocks based on 3D X-ray computed tomographic images by hidden Markov random field and Brinkman-force lattice Boltzmann model. *Journal of Hydrology*, **599**, 126377. https://doi.org/10.1016/j.jhydrol.2021.126377

Yang, F., Stack, A. G., & Starchenko, V. (2021). Micro-continuum approach for mineral precipitation. *Scientific Reports*, **11**(1), 1–14. https://doi.org/10.1038/s41598-021-82807-y

Yang, G., Coltman, E., Weishaupt, K., Terzis, A., Helmig, R., & Weigand, B. (2019). On the Beavers-Joseph interface condition for non-parallel coupled channel flow over a porous structure at high Reynolds numbers. *Transport in Porous Media*, **128**(2), 431–457. https://doi.org/10.1007/s11242-019-01255-5

Zhang, Q., Deng, H., Dong, Y., Molins, S., Li, X., & Steefel, C. (2022). Investigation of coupled processes in fractures and the bordering matrix via a micro-continuum reactive transport model. *Water Resources Research*, **58**(2), e2021WR030578. https://doi.org/10.1029/2021WR030578