The Kinetic Theory of water flow in porous media - Lausanne seminars

 Work in progress

Understanding the Mathematics of the Richards as a limit paper

In Part 1 we built the toolkit on finite matrices: a Laplacian-like operator with \(\ker = \mathrm{span}\{\mathbf{1}\}\), the Fredholm alternative as the source of macroscopic equations, and the pseudo-inverse as the source of transport coefficients. Now we let the matrix indices become continuous and watch the same algebra, verbatim, derive Richards' equation. This post is a plain-language companion to the second of our two Physical Review E manuscripts, where Richards' equation is obtained as the hydrodynamic (Chapman–Enskog) limit of a kinetic theory of unsaturated soil water.

From nodes to pore classes

In the kinetic theory the state of the soil water at a macroscopic point \(x\) and time \(t\) is not a single number \(\theta(x,t)\) but a whole distribution: the pore-occupancy function \(g(r, x, t)\), telling us how water is apportioned among pore classes of radius \(r\). The water content is recovered as a moment,

$$ \theta(x,t) = \int g(r,x,t)\, \mu(dr), $$

with \(\mu\) the pore-size measure of the medium. The “nodes” of Part 1 have become the continuum of pore radii \(r\); a “vector” \(\mathbf{f}\) has become a function \(f(r)\); the dot product has become an integral, \(\langle f, h \rangle = \int f(r)\, h(r)\, \mu(dr)\). Nothing else changes.

Water is exchanged between pore classes — capillary rearrangement, film flow, local equilibration — and this exchange is encoded by a linear(ized) operator \(\mathcal{I}\) acting on functions of \(r\). Schematically, and up to the details spelled out in the papers,

$$ (\mathcal{I} f)(r) = \int \kappa_s(r, r')\, \big[ f(r') - f(r) \big]\, \mu(dr'), $$

with a symmetric pair conductance built as the harmonic mean of the single-class conductances,

$$ \kappa_s(r, r') := \frac{2\, \kappa(r)\, \kappa(r')}{\kappa(r) + \kappa(r')}. $$

Compare this with \((L\mathbf{f})_i = \sum_j A_{ij}(f_j - f_i)\) (sign flipped): \(\mathcal{I}\) is a weighted graph Laplacian on a continuum of nodes, with \(-\mathcal{I}\) playing the role of \(L\). The harmonic mean is not decoration — it is the series-resistor composition rule of Part 1, guaranteeing that exchange between two classes is throttled by the less conductive of the two, and it makes \(\kappa_s\) manifestly symmetric, hence \(\mathcal{I}\) self-adjoint.

The three properties, revisited

Every structural fact from Part 1 now reappears with physical meaning attached.

(i) \(\ker \mathcal{I} = \mathrm{span}\{\mathbf{1}\}\) — one collision invariant. The exchange operator annihilates constants because pairwise exchange conserves total water: \(\int (\mathcal{I}f)\, \mu(dr) = 0\) identically, by the antisymmetry of the integrand. In the Boltzmann theory of gases the collision operator has a five-dimensional kernel (mass, three momenta, energy), and the hydrodynamic limit correspondingly produces five balance equations — the compressible Euler/Navier–Stokes system. In soil water the exchange between pore classes conserves only mass: the kernel is one-dimensional, and the hydrodynamic limit produces exactly one balance equation. That equation is Richards'. The dimension of a kernel dictates the size of your macroscopic PDE system — I find this one of the cleanest structural insights the kinetic viewpoint offers.

(ii) \(\mathcal{I}\) is self-adjoint and negative semi-definite — an H-theorem. The continuum version of the sum-of-squares identity of Part 1 reads

$$ \langle f, \mathcal{I} f \rangle = -\tfrac{1}{2} \iint \kappa_s(r,r')\, \big[f(r') - f(r)\big]^2\, \mu(dr)\, \mu(dr') \ \le\ 0, $$

with equality iff \(f\) is constant (on a “connected” pore network, in the sense that \(\kappa_s\) does not decompose the pore space into non-communicating blocks). Exchange strictly dissipates any non-uniformity: this is the H-theorem of the model, and the reason equilibrium exists and is unique. The equilibrium itself is a packing state — pores fill in order of capillary strength, a Fermi-sea-like picture — but for the linear algebra all we need is that fluctuations around it relax under a self-adjoint, negative semi-definite \(\mathcal{I}\).

(iii) Spectral gap — the small parameter exists. Because the zero eigenvalue is simple and isolated, there is a gap \(\lambda_1 > 0\), hence a fastest-conserved and slowest-decaying separation of time scales. Its ratio to the macroscopic time defines the Damköhler number, \(\mathrm{Da} := \tau_{eq}/\tau_{mac}\), and \(\mathrm{Da} \ll 1\) is the regime in which a hydrodynamic description can be honest. When \(\mathrm{Da}\) is not small — coarse structured soils, preferential flow, rapid forcing — the fast modes never fully slave to \(\theta\) and Richards' equation degrades. The linear algebra even tells you how it degrades: through the modes just above the gap.

The Chapman–Enskog march, order by order

Write the kinetic equation schematically as

$$ \partial_t g + (\text{transport in } x) = \frac{1}{\mathrm{Da}}\, \mathcal{I} g, $$

and expand \(g = g^{(0)} + \mathrm{Da}\, g^{(1)} + \dots\). The algebra of Part 1 now executes itself.

Order \(\mathrm{Da}^{-1}\): \(\mathcal{I} g^{(0)} = 0\), so \(g^{(0)}\) lies in the kernel — it is the local equilibrium distribution, parametrized by the single conserved moment \(\theta(x,t)\). The population of pore classes is enslaved to the water content.

Order \(\mathrm{Da}^{0}\): an equation of the form \(\mathcal{I} g^{(1)} = \mathcal{S}[g^{(0)}]\), with \(\mathcal{S}\) collecting the transport terms. This is exactly the singular problem \(L\mathbf{u} = \mathbf{b}\) of Part 1. The Fredholm alternative demands \(\langle \mathbf{1}, \mathcal{S} \rangle = 0\): projecting the transport terms onto the conserved direction. That projection is the continuity equation,

$$ \partial_t \theta + \nabla \cdot \mathbf{q} = 0. $$

The macroscopic balance law is not assumed; it is the solvability condition of a singular linear problem.

The flux and the conductivity: granted solvability, the correction is \(g^{(1)} = \mathcal{I}^{+} \mathcal{S}\) — the pseudo-inverse at work — and inserting it into the flux moment yields a Buckingham–Darcy law, \(\mathbf{q} = -K(\theta)\, \nabla (\psi(\theta) + z)\)-type, in which the hydraulic conductivity emerges as a bracket of the form

$$ K \ \sim\ -\,\big\langle \Phi,\ \mathcal{I}^{+} \Phi \big\rangle, $$

with \(\Phi(r, r') = -\Phi(r', r)\) the antisymmetric driving potential of the exchange. Readers of Part 1 will recognize the structure immediately: it is the Green–Kubo / effective-resistance formula, the continuum sibling of \(R_{ij}\) built from \(L^+\). The conductivity of a soil is, quite literally, the inverse Kirchhoff-type resistance of its pore-class network, evaluated on the mode forced by gravity and capillarity. This is where the empirical shapes of \(K(\theta)\) — Mualem, van Genuchten and relatives — acquire the status of approximations to a spectral object.

The variational subtlety, honestly told. In an earlier version of the manuscript we characterized this bracket by a single-field quadratic functional, Cercignani-style. That is legitimate when the bracket is a genuine quadratic form \(\langle \chi, \mathcal{I}\chi \rangle\); it silently fails when the object of interest is a bilinear pairing between two different functions. The fix is a two-field, primal–dual functional \(\mathcal{I}[\tilde\chi, \tilde\chi^*]\), stationary in each argument separately, whose stationary value returns the bilinear bracket without smuggling in a false symmetry. On a 3×3 matrix this distinction is invisible; in function space it decides whether your variational bound on \(K(\theta)\) is a theorem or wishful thinking. I mention it because it is a good example of how the finite-dimensional intuition of Part 1, taken too casually, can bite.

Coda: why bother

One can use Richards' equation for a lifetime without this machinery. The point of the derivation is not to re-obtain a 1931 result; it is that every object in the equation now has an address. \(\theta\) is the kernel coordinate; the continuity equation is a Fredholm solvability condition; \(K(\theta)\) is a pseudo-inverse bracket over the pore-class network; and the validity of the whole enterprise is a statement about a spectral gap, quantified by \(\mathrm{Da}\). When the equation fails — and we all know soils where it does — the linear algebra tells you which assumption broke, and what the next term in the expansion looks like.

References and further reading

• S. Chapman and T. G. Cowling, The Mathematical Theory of Non-Uniform Gases, 3rd ed., Cambridge University Press, 1970. (The original Chapman–Enskog method.)
• C. Cercignani, The Boltzmann Equation and Its Applications, Springer, 1988. (Linearized collision operator, Fredholm alternative, variational principles for transport coefficients — the template we adapt from gases to soils.)
• H. Grad, “Asymptotic theory of the Boltzmann equation,” Physics of Fluids 6:147–181, 1963. (The role of the collision invariants and the hydrodynamic projection.)
• F. Golse, “The Boltzmann equation and its hydrodynamic limits,” in Handbook of Differential Equations: Evolutionary Equations, Vol. 2, Elsevier, 2005. (A modern, rigorous survey of hydrodynamic limits.)
• L. A. Richards, “Capillary conduction of liquids through porous mediums,” Physics 1:318–333, 1931.
• E. Buckingham, “Studies on the movement of soil moisture,” USDA Bureau of Soils Bulletin 38, 1907.
• Y. Mualem, “A new model for predicting the hydraulic conductivity of unsaturated porous media,” Water Resources Research 12:513–522, 1976; M. Th. van Genuchten, Soil Science Society of America Journal 44:892–898, 1980. (The empirical \(K(\theta)\) shapes reinterpreted here as spectral approximations.)
• [PRE-1] R. Rigon et al., A kinetic theory of unsaturated soil water — manuscript, Physical Review E (submitted). (insert final title/DOI)
• [PRE-2] R. Rigon et al., Richards' equation as a Chapman–Enskog hydrodynamic limit — manuscript, Physical Review E (submitted). (insert final title/DOI)