View source for Finite water-content vadose zone flow method

The '''finite water-content vadose zone flux method'''<ref name="TO2008">Talbot, C.A., and F. L. Ogden (2008), A method for computing infiltration and redistribution in a discretized moisture content domain, ''Water Resour. Res.'', 44(8), doi: 10.1029/2008WR006815.</ref><ref name="Ogden_etal2015a">Ogden, F. L., W. Lai, R. C. Steinke, J. Zhu, C. A. Talbot, and J. L. Wilson (2015), A new general 1-D vadose zone solution method, ''Water Resour.Res.'', 51, doi:10.1002/2015WR017126.</ref> represents a one-dimensional alternative to the numerical solution of [[Richards' equation]]<ref>Richards, L. A. (1931), Capillary conduction of liquids through porous mediums, ''J. Appl. Phys.'', 1(5), 318–333.</ref> for simulating the movement of water in [[Vadose zone|unsaturated]] soils.  The finite water-content method solves the advection-like term of the [[Soil Moisture Velocity Equation]], which is an [[ordinary differential equation]] alternative to the Richards  [[partial differential equation]].  The Richards equation is difficult to approximate in general because it does not have a [[Closed-form expression|closed-form]] analytical solution except in a few cases.<ref name="Ross_Par_1994">Ross, P.J., and J.-Y. Parlange (1994).  Comparing exact and numerical solutions of Richards equation for one-dimensional infiltration and drainage.  Soil Sci. Vol 1557, No. 6, pp. 341-345.</ref>  The finite water-content method, is perhaps the first generic replacement for the numerical solution of the [[Richards' equation]].  The finite water-content solution has several advantages over the [[Richards equation]] solution.  First, as an [[ordinary differential equation]] it is explicit, guaranteed to converge <ref name="Yu_etal2012">Yu, H., C. C. Douglas, and F. L. Ogden, (2012), A new application of dynamic data driven system in the Talbot–Ogden model for groundwater infiltration, ''Procedia Computer Science'', 9, 1073–1080.</ref> and computationally inexpensive to solve.  Second, using a [[finite volume]] solution methodology it is guaranteed to conserve mass.  The finite water content method readily simulates sharp wetting fronts, something that the Richards solution struggles with.<ref>Tocci, M. D., C. T. Kelley, and C. T. Miller (1997), Accurate and economical solution of the pressure-head form of Richards' equation by the method of lines, ''Adv. Wat. Resour''., 20(1), 1–14.</ref>  The main limiting assumption required to use the finite water-content method is that the soil be homogeneous in layers.

[[File:Fwc domain.png|thumbnail|Finite Water-Content discretization.  The porous medium is divided into ''n'' uniform "bins" of <math>\Delta \theta</math> water content.]]

The finite water-content vadose zone flux method is derived from the same starting point as the derivation of [[Richards' equation]].  However, the derivation employs a [[hodograph transformation]]<ref>Philip, J. R. 1957. The theory of infiltration: 1. The infiltration equation and its solution, ''Soil Sci'', 83(5), 345–357.</ref> to produce an advection solution that does not include soil water diffusivity, wherein <math>z</math> becomes the dependent variable and <math>\theta</math> becomes an independent variable:<ref name="Ogden_etal2015a"/>

:<math>\left(\frac{dz}{dt}\right)_\theta = \frac{\partial K(\theta)}{\partial \theta} 
\left[ 1- \left (\frac{\partial \psi(\theta)}{\partial z}\right) \right]\ 
</math>

where:
:<math>K</math> is the unsaturated [[hydraulic conductivity]] [L T<sup>−1</sup>],
:<math>\psi</math> is the capillary [[pressure head]] [L] (negative for unsaturated soil),
:<math>z</math> is the vertical coordinate [L] (positive downward), 
:<math>\theta</math> is the [[water content]], (−) and
:<math>t</math> is [[time]] [T].

This equation was converted into a set of three [[ordinary differential equations]] (ODEs) <ref name="Ogden_etal2015a"/> using the Method of Lines<ref>{{Cite journal|doi=10.4249/scholarpedia.2859|title=Method of lines|year=2007|last1=Griffiths|first1=Graham|last2=Schiesser|first2=William|last3=Hamdi|first3=Samir|journal=Scholarpedia|volume=2|issue=7|page=2859|bibcode=2007SchpJ...2.2859H|doi-access=free }}</ref> to convert the [[partial derivatives]] on the right-hand side of the equation into appropriate [[finite difference]] forms.  These three ODEs represent the dynamics of infiltrating water, falling slugs, and capillary groundwater, respectively.

==Derivation==
A superior derivation was published<ref>Ogden, F.L., M.B. Allen, W.Lai, J. Zhu, C.C. Douglas, M. Seo, and C.A. Talbot, 2017.  The Soil Moisture Velocity Equation, ''J. Adv. Modeling Earth Syst.''https://doi.org/10.1002/2017MS000931</ref> in 2017, showing that this equation is a diffusion-free version of the [[Soil Moisture Velocity Equation]].

One way to solve this equation is to solve it for <math>q(\theta,t)</math> and <math>z(\theta,t)</math> by integration:<ref>Wilson, J. L. (1974), Dispersive mixing in a partially saturated porous medium, PhD dissertation, 355pp., Dept. of Civil Engrg., Mass. Inst. Tech., Cambridge, MA.</ref>

: <math>  \int \frac{\partial q}{\partial \theta} \, d\theta = \int \frac{\partial z}{\partial t} \, d\theta</math>

Instead, a finite water-content discretization is used and the integrals are replaced with summations:

: <math> \sum_{j=1}^N \left[\frac{\partial q}{\partial \theta}\right]_j \Delta \theta = \sum_{j=1}^N \left[\frac{\partial z}{\partial t}\right]_j \Delta \theta</math>

where <math>N</math> is the total number of finite water content bins.

Using this approach, the conservation equation for each bin is:

: <math> \left[\frac{\partial q}{\partial \theta}\right]_j = \left[\frac{\partial z}{\partial t}\right]_j.</math>

The method of lines is used to replace the partial differential forms on the right-hand side into appropriate finite-difference forms.  This process results in a set of three ordinary differential equations that describe the dynamics of infiltration fronts, falling slugs, and groundwater capillary fronts using a finite water-content discretization.

==Method essentials==

The finite water-content vadose zone flux calculation method replaces the [[Richards' equation]] [[Partial differential equation|PDE]] with a set of three [[ordinary differential equations]] (ODEs).  These three ODEs are developed in the following sections.   Furthermore, because the finite water-content method does not explicitly include soil water diffusivity, it necessitates a separate capillary relaxation step.  Capillary relaxation <ref name="Moeb_etal2012">Moebius, F., D. Canone, and D. Or (2012), Characteristics of acoustic emissions induced by fluid front displacement in porous media, ''Water Resour. Res.'', 48(11), W11507, doi:10.1029/2012WR012525.</ref> represents a free-energy minimization process at the pore scale that produces no advection beyond the REV scale.

===Infiltration fronts===
[[File:Inf fronts.png|thumbnail|Infiltration fronts in finite water-content domain.]]
With reference to Figure 1, water infiltrating the land surface can flow through the pore space between <math>\theta_d</math> and <math>\theta_i</math>.  In the context of the [[method of lines]], the partial derivative terms are replaced with:

: <math>\frac{\partial K(\theta)}{\partial \theta}=\frac{K(\theta_d)-K(\theta_i)}{\theta_d-\theta_i}.</math>

Given that any ponded depth of water on the land surface is <math>h_p</math>, the Green and Ampt (1911)<ref>Green, W. H., and G. A. Ampt (1911), Studies on soil physics, 1, The flow of air and water through soils, ''J. Agric. Sci.'', 4(1), 1–24.</ref> assumption is employed,

: <math>\frac{\partial \psi(\theta)}{\partial z}=\frac{|\psi(\theta_d)|+h_p}{z_j},</math>

represents the capillary head gradient that is driving the flow.  Therefore the finite water-content equation in the case of infiltration fronts is:

: <math>\left(\frac{dz}{dt}\right)_j= \frac{K(\theta_d)-K(\theta_i)}{\theta_d-\theta_i} \left(\frac{|\psi(\theta_d)|+h_p}{z_j}+1\right).</math>

===Falling slugs===
[[File:Falling slugs.png|thumbnail|Falling slugs in the finite water-content domain.  The water in each bin is considered a separate slug.]]
After rainfall stops and all surface water infiltrates, water in bins that contains infiltration fronts detaches from the land surface.  Assuming that the capillarity at leading and trailing edges of this 'falling slug' of water is balanced, then the water falls through the media at the incremental conductivity associated with the  <math>j</math>-th <math> \Delta\theta</math> bin:

: <math> \left(\frac{dz}{dt}\right)_j = \frac{ K(\theta_j) - K(\theta_{j-1}) }{\theta_j -\theta_{j-1}}</math>

===Capillary groundwater fronts===
[[File:Groundwater fronts.png|thumbnail|Groundwater capillary fronts in finite water-content domain.]]
In this case, the flux of water to the <math> j^\text{th}</math> bin occurs between bin ''j'' and ''i''.  Therefore in the context of the [[method of lines]]:

: <math>\frac{\partial K(\theta)}{\partial \theta}= \frac{K(\theta_j)-K(\theta_i)}{\theta_j - \theta_i}, </math>

and,

: <math> \frac{\partial\psi(\theta)}{\partial z} = \frac{|\psi(\theta_j)|}{H_j} </math>

which yields:

: <math>\left(\frac{dz}{dt}\right)_j= \frac{K(\theta_j)-K(\theta_i)}{\theta_j - \theta_i} \left(\frac{|\psi(\theta_j)|}{H_j}-1\right). </math>

The performance of this equation was verified for cases where the groundwater table velocity was less than 0.92 <math>K_s</math>,<ref name="Ogden_etal2015b">Ogden, F. L., W. Lai, R. C. Steinke, and J. Zhu (2015b), Validation of finite water-content vadose zone dynamics method using column experiments with a moving water table and applied surface flux, ''Water Resour. Res.'', 10.1002/2014WR016454.</ref> using a column experiment fashioned after that by Childs and Poulovassilis (1962).<ref>Childs, E. C., and A. Poulovassilis (1962), The moisture profile above a moving water table, ''J. Soil Sci''., 13(2), 271–285.</ref>  Results of that validation showed that the finite water-content vadose zone flux calculation method performed comparably to the numerical solution of Richards' equation.

===Capillary relaxation===
Because the hydraulic conductivity rapidly increases as the water content moves towards saturation, with reference to Fig.1, right-most bins in both capillary groundwater fronts and infiltration fronts can "out-run" their neighbors to the left.  In the finite water content discretization, these shocks<ref>Smith, R. E. (1983), Approximate soil water movement by kinematic characteristics,  ''Soil Sci. Soc. Am. J.'', 47(1), 3–8.</ref> are dissipated by the process of capillary relaxation, which represents a pore-scale free-energy minimization process that produces no advection beyond the REV scale<ref name="Moeb_etal2012"/>  Numerically, this process is a numerical sort that places the fronts in monotonically-decreasing magnitude from left-right.

==Constitutive relations==
The finite water content vadose zone flux method works with any monotonic [[water retention curve]]/unsaturated hydraulic conductivity relations such as Brooks and Corey<ref>Brooks, R.H., and A.T. Corey, 1964.  Hydraulic properties of porous media.  Hydrol. Pap. 3, Colo. State Univ., Fort Collins, Colorado, USA.</ref> Clapp and Hornberger <ref>Clapp R.B., and G.M. Hornberger, 1978.  Empirical equations for some soil hydraulic properties.  ''Water Resour. Res.''  14(4):601–604</ref> and van Genuchten-Mualem.<ref>van Genuchten, M. Th. (1980). "A Closed-Form Equation for Predicting the Hydraulic Conductivity of Unsaturated Soils" (PDF). ''Soil Sci. Soc. Am. J.'', 44 (5): 892-898. doi:10.2136/sssaj1980.03615995004400050002x</ref>  The method might work with hysteretic water retention relations- these have not yet been tested.

==Limitations==
The finite water content method lacks the effect of soil water diffusion.  This omission does not affect the accuracy of flux calculations using the method because the mean of the diffusive flux is small.  Practically, this means that the shape of the wetting front plays no role in driving the infiltration.  The method is thus far limited to 1-D in practical applications.  The infiltration equation <ref name="Ogden_etal2015a"/> was extended to 2- and quasi-3 dimensions.<ref name="Yu_etal2012"/>  More work remains in extending the entire method into more than one dimension.

==Awards==
The paper describing this method <ref name="Ogden_etal2015a"/> was selected by the Early Career Hydrogeologists Network of the International Association of Hydrogeologists to receive the "Coolest paper Published in 2015" award in recognition of the potential impact of the publication on the future of hydrogeology.

==See also==
* [[Richards' equation]]
* [[Infiltration (hydrology)]]
* [[Soil Moisture Velocity Equation]]

==References==
{{reflist|2}}

[[Category:Soil physics]]
[[Category:Hydrology]]
[[Category:Partial differential equations]]
[[Category:Ordinary differential equations]]