Pore network modeling

        Variably-saturated flow and multi-component reactive adsorptive transport

Drainage and Flow Simulations


To simulate a drainage process, the nonwetting phase is considered to enter the network through an external reservoir which is connected to the inlet side of the network. The displaced phase escapes through the outlet face on the opposite side. We assume that the progress of the displacement is controlled by the capillary forces. This forms the basis for the invasion-percolation algorithm used to model drainage. At every stage of the process, nonwetting phase invades all accessible pore bodies and throats with the lowest entry capillary pressure.

where Pw and Pn are the pressure of the wetting and the nonwetting phases, respectively, r* is the mean radius of curvature and Ƴwn is the surface tension. For a capillary tube of radius r, we have r*=r/cos(ϑ) (Young-Laplace’s equation), in which ϑ is the contact angle between fluid interface and capillary wall.

The invading fluid enters and fills an available pore throat only when the injection pressure is equal to or larger than the entry capillary pressure of the pore. The wetting phase is hydraulically connected though filaments within edges of angular pore bodies and pore throats. The capillary pressure is increased incrementally so that fluid-fluid interfaces will move only a short distance before coming to rest in equilibrium at the opening of slightly smaller pore throats.

PoreFlow uses the generated pore networks to simulate drainage process and to calculate pressure and velocity distributions in such a network. A new approach is formulated for a more accurate modeling of flow under variably saturated conditions. The continuous pore space is descretized into “pore elements”, each with its own pressure, concentration, and adsorbed mass, etc. The discretization depends on saturation state of a pore. In the case of a saturated pore, one pressure and one concentration values are assigned for each pore element representing either a saturated pore body or a pore throat. In the presence of the nonwetting phase, a given pore body or pore throat could be invaded and filled partially by the nonwetting phase, forcing the wetting phase to flow only along the edges. Under such conditions, pore throats are not necessarily the narrowest constriction along the flow path, anymore. In fact, resistance to the flow within the edges of drained pore bodies may be comparable to, or even larger than, the resistance to the flow within the pore throats. Also, the lower conductivity of a drained pore body reduces the connectivity and therefore the flow/solute flux among saturated pore throats connected to a drained pore body. Under such conditions, a drained pore body is divided into corner units, each of them being a pore element.

(a) Example of two drained pore-bodies connected to each other using a drained pore throat with angular cross section. Each drained pore body is refined into 8 corner units (pore elements), and (b) a pore body which is invaded by the nonwetting phase through one of its throats, and, as a result, reduces the connectivity of the neighboring saturated pore throats (such as pore throats number 2, 3, 4).


A corner unit is composed of a corner together with half pore body edges connected to it. Fluid flow between these elements occurs through edges within drained pore body.

In a drained pore throat with angular cross section, there will be flow of wetting phase only through the edges of the pore throats. Since, in an angular pore, edges may have different angles, they have different flow rates, and presence of the nonwetting phase limits mixing between them. To take this effect into account, each edge of a drained pore throat is treated as a pore element with its own flow rate.

After performing the drainage step, knowing the curvature of interfaces and local saturations in each pore, PoreFlow calculates conductances of edges of the drained pores throats and drained pore bodies as well as conductances of saturated pores.

Imposing a pressures difference across the network, we assume that the discharge through a given pore throat can be prescribed by Hagen-Poiseuille equation:

where qi j;tot is the total volumetric flow rate through pore throat ij, gi j;k is conductance of kth edge of pore throat ij, and pi and pj are pressures at pore elements i and j, respectively.

For incompressible, steady-state flow, the sum of discharges into and out a pore body, or a pore-body corner unit in the case of a drained pore body, must be zero. For the pore bodies, the continuity equation may be written as:

where NCU,iedge shows the number of edges through which corner unit i, within a drained pore body, is connected to other corner units, n, within the same pore body. zi is the coordination number of pore i.

Solving for pore body pressures and flow velocities in all pore throats, considering the network as an REV, the average pore velocity v can be determined.

where Qtot is the total discharge through the network which is the sum of fluxes through all pore throats at the inlet or outlet boundary of the network, L is the network length, and Vf is total fluid volume. The relative permeability of the network to the wetting phase at a given saturation and capillary pressure is calculated using Darcy’s law.

where µ is viscosity, k is intrinsic permeability, and ΔP is the pressure difference between inflow and outflow reservoirs. Repetition of this process at consecutively larger imposed capillary pressures results in a graph of capillary pressure versus saturation and relative permeability versus saturation.

Relative permeability curves obtained using PoreFlow and study by Al-Kharusi and  Blunt (2008) together with the measured values during primary drainage experiment (after Raoof and Hassanizadeh, 2011).



Raoof, A. and Prof. M. Hassanizadeh (2011), A new formulation for pore-network modeling of two-phase flowWater Resour. Res., doi:10.1029/2010WR010180


Al-Kharusi, A. S., and M. J. Blunt (2008), Multiphase flow predictions from carbonate pore space images using extracted network models, Water Resour. Res., 44, W06S01, doi:10.1029/2006WR005695.

 Pore network modeling

        Variably-saturated flow and multi-component reactive adsorptive transport