Numerical modeling of variably saturated flow in porous media

Richards' equation (RE) with associated constitutive relations is a simple mathematical representation of saturated/unsaturated (variably saturated) flow in porous media. Highly nonlinear constitutive relations make analytical solutions impossible, except for special cases. The nonlinearity also leads to computationally expensive numerical solutions. Many numerical methods have been proposed, yet there remains room for improvement in terms of efficiency and robustness. Effective numerical solutions to RE are of substantial practical importance in providing solutions to difficult classes of groundwater flow problems.; For problems involving porous media with non-uniform grain size distributions, constitutive relations described using the van Genuchten and Mualem models lead to numerical convergence difficulties. These difficulties arise because of discontinuities in the derivative of specific moisture capacity and relative permeability as a function of capillary pressure. Convergence difficulties are illustrated using standard numerical approaches to simulate such problems. Constitutive relations, interblock permeability, nonlinear algebraic system approximation methods, and two time integration approaches are investigated. An integral permeability approach approximated by Hermite polynomials is recommended and shown to be robust and economical.; Sharp fronts with rapid changes in fluid saturations over short distance and time scales often exist in variably saturated flow. Standard solution approaches are often inefficient and unreliable for such highly nonlinear problems, yet alternative approaches such as transformation methods may provide for more efficient solutions. Four existing transforms and one new transform are summarized and compared for solving RE within a common framework over a wide range of media properties. Results show that transformation methods can significantly improve solution efficiency and robustness compared to standard solution approaches, and that the new transform compares favorably with existing transforms.; An effective variable time-step approach has recently been introduced using a differential algebraic equation (DAE) implementation of the method of lines (MOL). The DAE/MOL approach gives solutions to the pressure-head form of RE that are accurate and more efficient than constant step-size or empirically-based variable step-size methods. Numerical results are presented which show that a combination of transformation and DAE/MOL solution methods provides the most accurate, efficient and robust solutions to RE.