Porous flow in a viscous medium

Governing equations

At the pressure and temperature conditions of the creeping mantle, « Fluids » (hydrous melts or aqueous fluids) are expected to migrate by porous flow though an interconnected network of channels formed along the grain edges [1].

In 1984, D. McKenzie [2] derived the governing equations that describe the transport of a liquid through a porous viscously deformable matrix based on the continuum mechanics theory. The liquid characterized by a very low viscosity relatively to that of the matrix is called the fluid phase and the matrix the solid phase.

These equations are :
  – Conservation of mass for both phases :
    \qquad \dfrac{\partial (\rho_f \phi)}{\partial t} + \nabla \cdot (\rho_f \phi \boldsymbol{v}_f) = \Gamma
    \qquad \dfrac{\partial (\rho_s \phi)}{\partial t} + \nabla \cdot (\rho_s (1-\phi) \boldsymbol{v}_s) = -\Gamma
  – Conservation of momentum :
    \phi(\boldsymbol{v}_f - \boldsymbol{v}_s) = -\dfrac{K}{\mu}(\boldsymbol{\nabla} P + \rho_f \boldsymbol{g})
    \nabla \cdot (2 \eta \dot{\boldsymbol{\epsilon}}) - \boldsymbol{\nabla} P + \nabla (\xi \nabla \cdot \boldsymbol{v_s}) = (\rho_f \phi + \rho_s (1- \phi) \boldsymbol{g})
where subscripts f and s denote values for, respectively, the fluid and the solid. \rho is density, \boldsymbol{v} is velocity, \phi is porosity or fluid fraction, and \Gamma is rate of mass transfer between the phases. K is permeability, \mu is fluid viscosity, p total pressure, \dot{\boldsymbol{\epsilon}} deformation rate, and \boldsymbol{g} gravity. \eta and \xi are, respectively, the solid shear and bulk viscosities.

Numerical strategy

The McKenzie’s equations are re-written in a form more suitable for a numerical resolution [see 3].

Two systems are defined assuming that densities of both phases are constant :
  – an « incompressible » Stokes-like solid system :
    \nabla \cdot (2 \eta \dot{\boldsymbol{\epsilon}}) - \boldsymbol{\nabla}p = \phi \Delta \rho \boldsymbol{g}
    \nabla \cdot \boldsymbol{v_s} = \dfrac{{\cal P}}{\xi}
  – a « compressible » fluid system.
    \dfrac{\partial \phi}{\partial t} +  \boldsymbol{v}_s  \cdot \boldsymbol{\nabla} \phi = (1-\phi)  \dfrac{{\cal P}}{\xi} + \dfrac{\Gamma}{\rho_s}
    \dfrac{{\cal P}}{\xi} - \nabla \cdot \left[ \dfrac{K}{\mu} (\boldsymbol{\nabla} {\cal P} + \boldsymbol{\nabla} p + \Delta \rho \boldsymbol{g} ) \right] =  \Gamma \dfrac{\Delta \rho}{\rho_f\rho_s}
where \Delta \rho is density contrast. {\cal P}  and p are, respectively, compaction and dynamic pressures.


[1] Wark, D. A., et al. « Reassessment of pore shapes in microstructurally equilibrated rocks, with implications for permeability of the upper mantle. » Journal of Geophysical Research: Solid Earth 108.B1 (2003).
[2] McKenzie, D. « The generation and compaction of partially molten rock. » Journal of Petrology 25.3 (1984): 713-765.
[3] Katz, R. F., et al. « Numerical simulation of geodynamic processes with the Portable Extensible Toolkit for Scientific Computation. » Physics of the Earth and Planetary Interiors 163.1 (2007): 52-68.


Web Analytics