Assessment of an efficient numerical solution of the Richards'equation for bare soil

A new numerical scheme has been proposed by Ross (2003) to solve the 1D Richards' equation (1931). This noniterative solution uses the description of soil properties proposed by Brooks and Corey (1964). It allows the derivation of an analytical expression for the Kirchhoff potential used in the calculation of water fluxes. The degree of saturation is used as the dependent variable when the soil is unsaturated and the Kirchhoff potential is used in case of saturation. A space and time discretisation scheme leads to a tridiagonal set of linear equations that is solved non-iteratively. We propose in this paper an extensive test of this numerical method, evaluated only on a single case by Ross (2003). The tests are conducted in two steps. First, the solution is assessed against two analytical solutions. The first one (Basha, 1999) provides the water content profile when simplified soil characteristics such as the exponential law of Gardner (1958) are used. The Ross (2003) solution is compared to this solution on eight different soils that were fitted to this law. Analytical solution with the Brooks and Corey (1964) models are not available at the moment for the moisture profile but some exist for cumulative infiltration. Therefore the second analytical solution, used in this study, is the one developed by Parlange et al. (1985) and Haverkamp et al. (1990). The various parameters of this solution are written for the Brooks and Corey (1964) models and the cumulative infiltration calculated by the Ross (2003) solution is compared to this analytical solution. The second part of the test is a comparison with a reference model: the SiSPAT (Simple Soil Plant Atmosphere Transfer) model, which provides a reference iterative solution of the soil water flow equation. A wide range of simulations is performed, with various soil types (homogeneous or not), various climate forcing and several initial conditions. It allows the comparison of various variables such as deep drainage, soil moisture profile, surface ponding, and evaporation under non uniform initial moisture content and time varying climate forcing. We show that the model provides robust and accurate solutions as compared with the analytical solutions and with the SiSPAT model. This study also shows that a finer discretization than the one proposed by Ross (2003) is necessary close to the soil surface to accurately model the cumulative infiltration, especially for clayey soils.