Published in Soil Sci. Soc. Am. J. 68:1-6 (2004).
© 2004 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
DIVISION S-1SOIL PHYSICS
Improvements to Estimating Unsaturated Soil Properties from Horizontal Infiltration
John S. Tyner* and
G. O. Brown
The University of Tennessee, 2506 E.J. Chapman Drive, Knoxville, TN 37996-4531
* Corresponding author (jtyner{at}utk.edu).
 |
ABSTRACT
|
|---|
Current methods to determine unsaturated soil properties are expensive, difficult, and time-consuming. A procedure is presented that quickly estimates the unsaturated hydraulic conductivity, diffusivity, water retention, and sorptivity functions from a BruceKlute test in conjunction with a semi-analytic analysis. Data collection was performed using a computer controlled
-ray attenuation system to measure moisture content. A computer controlled syringe pump maintained the inlet boundary conditions. Improvements to existing analytical procedures include the independent optimization of van Genuchten soil hydraulic parameters n and
, and a method to optimize the residual moisture content from a single measurement of water potential versus water content. The method was applied to a sandy loam and its strengths and weaknesses are discussed.
 |
INTRODUCTION
|
|---|
KNOWLEDGE OF SOIL HYDRAULIC PROPERTIES is required to estimate fluid flow in unsaturated porous media. Commonly, numerical models are employed to estimate fluid flow in a variety of settings. Although much effort has been expended to increase the sophistication of numerical models, lack of sufficient model input is arguably a larger problem for many modelers. Because direct measurement of unsaturated hydraulic functions is difficult and expensive, properties are often estimated by applying the Mualem (Mualem, 1976) or Burdine (Burdine, 1953) model to a water retention curve described by Brooks and Corey (1966) or van Genuchten (1980). Many practitioners have adopted the additional step of estimating soil water retention from grain-size distribution, bulk density, and other easily measured soil properties (Arya and Paris, 1981).
The method presented within this manuscript describes an improved procedure to estimate soil hydraulic functions using data collected from a BruceKlute test (Bruce and Klute, 1956). We employed a computer-controlled
-ray attenuation system to measure volumetric moisture content and a computer-controlled syringe pump to maintain the Bruce-Klute inlet boundary condition. An improved data analysis procedure is also presented, which permits independent optimization of van Genuchten soil hydraulic parameters n and
. Additionally, the procedure optimizes an additional parameter, residual moisture content, from a single measurement of water potential versus water content.
 |
THEORY
|
|---|
Van Genuchten Parameters
This method uses van Genuchten's expression for water retention (van Genuchten, 1980), although similar versions utilizing Brooks and Corey (Brooks and Corey, 1966) or numerous other versions could be used. Van Genuchten describes the water retention curve by
 | [1] |
where
is the normalized water content,
is the matric potential, and
, m, n,
s, and
r are fitting parameters. The parameters,
s and
r, are the saturated water content and residual water content, respectively. By letting
 | [2] |
and applying the Mualem model (Mualem, 1976) to relate water retention and unsaturated hydraulic conductivity, van Genuchten obtains
 | [3] |
where Ks is the saturated hydraulic conductivity. Using the definition of diffusivity (Klute, 1952),
 | [4] |
van Genuchten also provides an expression for D(
)
 | [5] |
where h is head. Summarizing, van Genuchten's expressions for unsaturated hydraulic functions contain six parameters,
, m, n,
s,
r, and Ks. As is common practice, we apply Eq. [2] to reduce the number of required parameters to five.
BruceKlute Test
Bruce and Klute (1956) provide a convenient method to measure D(
). The method represents horizontal, semi-infinite, one-dimensional flow as
 | [6] |
where x is the horizontal distance from the inlet, and t is the elapsed time. The initial and boundary conditions for the test are
 | [7] |
where
n is the antecedent water content and
o is the inlet water content. These conditions are traditionally referred to as the constant concentration boundary condition. The Boltzmann transformation (Boltzmann, 1894) reduces Eq. [6] and [7] to the ordinary differential equation,
 | [8] |
subject to,
 | [9] |
where
= x/
. Integration of Eq. [8] with respect to
yields
 | [10] |
Philip (1969) showed that the inlet flux, qo, is given by
 | [11] |
where qo is the water flux at the inlet (x = 0) and S(
n,
o) is the sorptivity. S(
n,
o) is defined as
 | [12] |
and is a constant for given values of
n and
o.
Equation [10] provides a direct method to calculate D(
); however, its application is problematic for three reasons.
- The slope of
(
) near the inlet and wetting front are difficult to measure accurately due to the very small and large slopes, respectively.
- The calculated diffusivities are tabular making their application within numerical models more complicated.
- The calculated D(
) values are not necessarily consistent with the measured
(
) data from the BruceKlute test. That is, there is no measure of the ability for D(
) to accurately predict the measured
(
).
Several papers have suggested methods to address difficulties I and II by initially fitting a function to the measured
(
) data and then numerically determining the slope and area of the curve (McBride and Horton, 1985; Meyer and Warrick, 1990). Others have also addressed difficulty III by employing solutions that ensure the estimated D(
) values will predict the measured
(
) curve (Clothier et al. 1983; Warrick, 1994). Given that D(
) at high water contents is large and sensitive to the difficult interpretation of d
/d
, small errors in the interpretation of d
/d
can lead to large errors of D(
) near saturation.
McWhorter's and Philip and Knight's Solutions
While developing equations to describe two-phase constant concentration boundaries, McWhorter (1971) introduced fractional flow, F(
), where
 | [13] |
q(
) is the flux density, and q0 is the flux density at x = 0. Inspection of Eq. [13] reveals that F(
) encompasses values from one at
=
0 to zero at
=
n. Using the same definition for F(
), Philip (1973) derived a similar expression for fractional flow given as
 | [14] |
Philip and Knight (1974) provide expressions for S(
n,
o) and
(
) as
 | [15] |
 | [16] |
and an exact, quasi-analytical solution of Eq. [14] using the iterative procedure
 | [17] |
where i represents the iteration of F(
), and ß is a variable of integration. Equation [17] is a special case of McWhorter's (1971) solution for horizontal flow of two viscous fluids. The rightmost expression is used during computations of F(
). As F(
)i converges to F(
)i+1, Eq. [17] converges to the exact solution
 | [18] |
rapidly for all D(
) represented by a monotonically increasing function, including D(
) defined using van Genuchten soil parameters. Brown and McWhorter (1990) found that the first guess,
 | [19] |
provides a stable convergence to a final estimate even for non-monotonic diffusivity functions. An expression to calculate S(
n,
0) is obtained by the substitution of Eq. [12] into Eq. [14] followed by differentiation with respect to
(Brown, 1987).
 | [20] |
Equation [6] can also be solved using a purely numeric solution (Simunek et al., 2000).
Clothier Et Al.'s and Warrick's Procedures
Clothier et al. (1983) showed the difficulty of evaluating Eq. [10] and then provided a method to ensure the estimated D(
) is consistent with measured
(
). They first fit a free-hand curve through measured
(
) data, and using Eq. [10], calculated D(
). Next, F(
) and S(
n,
o) were determined using Eq. [14] and [15], respectively. Finally,
(
) was estimated using Eq. [16], and then compared with the measured
(
). The estimate for S(
n,
o) using Eq. [15] was found to be 11% larger than S(
n,
o) calculated from Eq. [12] using the measured data. They emphasized that this disparity "may be fortuitously small." To avoid the discrepancy between measured and predicted
(
), they proposed estimating D(
) using analytical expressions derived in Philip (1960) that allow expressions for D(
) to be fit directly to the measured
(
), thereby ensuring proper scaling of measured and predicted S(
n,
o). A drawback of this method is that it does not utilize familiar soil parameters such as Brooks and Corey (1966) or van Genuchten (1980), thereby limiting the utilization of such parameters beyond expressions for D(
).
McBride and Horton (1985) stated their method to determine D(
) from a BruceKlute test compares favorably to that of Clothier et al. (1983). This statement was founded on the fact that their
(
) function fits measured
(
) data better than that of Clothier et al. (1983), although they neglected to show that D(
) calculated by their method accurately predicts the measured
(
) as is done analytically in the Clothier et al. solution. There exists a multitude of
(
) expressions to fit measured
(
) data. Since the ultimate goal of measuring hydraulic properties is to predict flow, arguments that one
(
) expression fits measured
(
) data slightly better than another are largely inconsequential unless they can also be shown to better predict
(
) from the estimates of D(
).
Using several common D(
) functions including Eq. [5], Warrick (1994) optimized soil parameters using the Philip (1969) finite difference solution of the constant concentration condition. This solution is similar to that of Philip and Knight (1974) in that it allows D(
) to predict
(
). By assuming
s and
r were known and equal to
0 and
n, respectively, the predicted
(
) was fit to the measured
(
) by simultaneously optimizing the van Genuchten m and a lumped parameter A(m),
 | [21] |
Since Eq. [5] approaches infinity as
0 approaches
s, it is assumed some type of numerical procedure was used to approximate D(
) at saturation, possibly the estimation method presented in Warrick et al. (1985). Essentially Warrick's method starts with five unknown parameters (Ks,
, m,
s, and
r), assumes two are known (
s, and
r) and optimizes two others [m and A(m)]. Selection of a D(
) function whose parameters are related to the water retention curve and unsaturated hydraulic conductivity function is quite appealing. Given independent knowledge of
s,
r, and either Ks or
, a soil's hydraulic functions can be fully described. However, measurement of Ks,
,
s, and
r requires laboratory procedures that are difficult, expensive, and time-consuming using more traditional techniques such as permeameters, hanging water columns, and pressure plates.
 |
MATERIALS AND METHODS
|
|---|
Improved Method
Like Warrick (1994), we choose to fit van Genuchten's description of D(
), Eq. [5], to the measured
(
). We prefer to solve the constant concentration condition using Philip and Knight (1974) due to both its relative ease and more analytical basis compared with Philip (1969). Additionally, our method offers two distinct improvements over previous methods. First, n is determined independently of Ks/
instead of using a lumped parameter such as Eq. [21]. Second, we provide a method to estimate
r from a single measurement of
versus
.
Applying the Boltzmann transformation to the collected data results in the calculation of
(
). The
axis is then normalized to
 | [22] |
where
wf is the largest value of
such that
>
n.
An estimate of D(
) is calculated using Eq. [5] with an initial guess for n, unity for the ratio Ks/
,
s is set equal to the porosity, and
r is initially estimated based on soil type or prior knowledge. Next, the initial D(
) and F(
)1 are entered into Eq. [17] followed by iteration until convergence of F(
) is achieved. Philip and Knight (1974) and personal experience demonstrate that four iterations are adequate for most soils. A theoretical water profile,
(
), is predicted using Eq. [15] and Eq. [20], and is then normalized to
(
) using Eq. [21]. By normalizing
(
) to
(
), the influence of the lumped parameter Ks/
is removed, which allows for optimization of n independent of Ks/
. Next, we optimize the value of n by minimizing the weighted absolute value of differences between measured and predicted
(
).
The functional independence of
(
) and Ks/
can be demonstrated by noting from Eq. [2] and Eq. [5] that the predicted D(
) is linearly related to Ks/
and is a complex function of n. Substituting Eq. [5] into Eq. [18] reveals that F(
) is independent of Ks/
, since it resides once in the numerator and denominator. However, n cannot be reduced from Eq. [18] revealing that F(
) and n are dependent. Equation [15] shows that S(
n,
o) is a function of Ks/
and n since D(
) is only present in the numerator. Thus, when Eq. [20] is substituted into the numerator and denominator of the right side of Eq. [21], S(
n,
o) is cancelled and
(
) is shown to be a function of n, but not Ks/
.
Using the previously determined n, we compare the predicted
(
) with the measured
(
) to optimize Ks/
, again using the weighted absolute value of differences. We can also check to verify that the measured S(
n,
o), given by Eq.[12], is equal to S(
n,
o) delivered to the column by the syringe pump. Agreement of the measured and delivered S(
n,
o) provides evidence that measurement of
(
) was conducted accurately.
The only parameter yet to be determined is
r. The values for
, m, and n are all functions of the initial estimate of
r. Since
r is a fitting parameter and is only defined in terms of fitting Eq. [1] to measured
(
) data, it cannot be measured directly. The
r is determined by first collecting a single measurement of
versus
at a moderately dry water content. Next, we optimize
r such that the prediction of
(
) by Eq. [1] is in agreement with the single measurement of
versus
. Using this second estimate for
r,
and n are re-optimized using Philip and Knight (1974) followed by another optimization of
r. This iterative process is repeated until subsequent estimates of
r converge.
Restating the procedure:
- Start with six unknown hydraulic parameters (
, m, n,
s,
r, and Ks) and apply Eq. [2] to determine m in terms of n.
- Measure or estimate
s based on the porosity.
- Conduct BruceKlute test and normalize the data by Eq. [21].
- Temporarily assume a value for
r based on soil type or other measure.
- Use Philip and Knight's (1974) solution to obtain a best-fit estimate of n.
- Using the non-normalized
(
) data and the previously determined n, conduct an additional solution of Philip and Knight's (1974) solution to obtain a best-fit estimate of Ks/
.
- Independently measure Ks and calculate
.
- Measure a single point on the
(
) curve. Using the estimates for
and n, estimate
r by applying Eq. [1].
- Using the new value for
r, iterate through Steps 5 to 8, with the exception of the measurement of Ks and a point from
(
), until convergence of
r is achieved.
Laboratory Procedures
The Slaughterville sandy loam (Coarse-loamy, mixed, superactive, thermic Udic Haplustoll), collected near Perkins, OK, was selected for testing. Soil was uniformly packed with a uniform water content,
n, into a 0.2 m long, 0.035-m diam. clear acrylic column. A clear column allows for visual verification of the water front location during testing, but is not necessary. Larger diameter columns are allowable for fine-textured soils; however, coarse-textured soils may exhibit nonvertical wetting fronts due to gravitational effects. Packing was performed by compacting 0.01-m lifts with a 0.01-m diam. tamping rod. After each lift was placed, the upper 0.01 m of packed soil was loosened, followed by the placement of the next lift. Dry bulk density was calculated from column volume and the dry mass of packed soil. Porosity was estimated from dry bulk density and an assumed particle density of 2650 kg m3. Saturated water content was assumed to equal the porosity, and
r was initially estimated based on soil type (Carsel and Parish, 1988). Saturated hydraulic conductivity, Ks, was measured using a constant head permeameter on a similarly prepared hand packed soil. Field samples may be used if they are homogenous and if water content is allowed to equilibrate before testing.
Water was injected into the soil column with a custom-built programmable syringe pump at a rate following Eq. [11] (Brown and Allred, 1992). Use of a syringe pump, instead of a traditional Mariotte flask, allows the inlet boundary condition to be maintained at any desired water content. In particular, if it is desired to maintain a relatively small
0 while testing a fine-grained soil, the large negative pressure required makes a Mariotte flask impractical. Use of Mariottes also tend to cause d
/d
to become concave upward near the inlet. Bruce and Klute (1956) discuss this behavior and that it leads to calculation of a maximum D(
) at less than saturation. Brown and Allred (1992) conducted a direct comparison of the two types of inlets and their results show that a concave upward d
/d
is only present with the Mariotte flask induced boundary condition. Equation [5] predicts that D(
)
as
1, which is consistent with a concave downward d
/d
near the inlet as produced by the syringe pump induced boundary condition.
During imbibition, the volumetric water content was measured at various distances from the inlet over time using
-ray attenuation procedures (Gardner, 1986). Photons from an Am241 source were detected with a Bicron 2 x 2 NaI(Ti) detector (Bicron, Newberry, OH) in conjunction with an Ortec Ace-2K multichannel analyzer (Ortec, Oak ridge, TN). Although the required
(
) curve can be measured from a single location, measurement from multiple locations over time allows confirmation that infiltration is following the Boltzmann transformation by verifying that all of the
(
) data lay atop a single curve. Alternatively, one can dissect the column immediately after imbibition and measure the water content gravimetrically as described in Bruce and Klute (1956). Measurement of
at a single time or location does not enable confirmation of the Boltzmann transformation and fewer useful data can be collected. Additionally, measuring
gravimetrically hinders the use of field-collected samples since the columns must be quickly dissected following imbibition. Quick dissection of unconsolidated material would be difficult given the types of column sheathes often used. Consolidated material might also be difficult to dissect quickly because of the hardness of the sample itself. To meet the semi-infinite boundary condition of Eq. [9], data collection must cease before the water front advancing to the end of the column.
 |
RESULTS AND DISCUSSION
|
|---|
The hand-packed soil had a dry bulk density of 1480 kg m3 and a porosity of 0.44 m3 m3, which was assumed to be static throughout the test. Saturated water content was set equal to porosity, and
r was initially estimated at 0.065. Initial and inlet water contents,
n and
0, were set to measured values of 0.035 and 0.37, respectively. Saturated hydraulic conductivity, Ks, was found to be 6.4 x 105 m s1. A sorptivity of 8.9 x 106 m2 s1/2 was selected and entered into Eq. [11] to calculate the syringe pump flux as a function of time. Data was collected for 55 min and water content,
, was measured at 0.10 and 0.11 m from the inlet; the resulting
(
) plot is shown in Fig. 1
. After normalizing
(
) to
(
), n was optimized to a value of 1.87 (Fig. 2)
. Using this estimate for n, the
(
) curve was used to optimize Ks/
to a value of 2.30 x 105 m2 s1, which leads to an
of 2.78 m1 (Fig. 3)
. A single measurement of
= 153 m versus
= 0.058 was collected using pressure plates. Equation [1] was solved for
r, which led to
r = 0.056. The results converged after four additional iterations.
Table 1 presents the method's estimates of n,
, and
r for each iteration. Figures 2 and 3 present the fits of Iteration 1a. The fits of the predicted
(
) to measured
(
) in subsequent iterations (2a5a) are similar, but the values of n and Ks/
change markedly as
r converges. Figure 4
presents the measured
(
) and the final iteration (5a) of the predicted
(
). Note the similarity of the predicted
(
) for Fig. 3 (Iteration 1a) and Fig. 4 (Iteration 5a). Similar results were achieved using an initial
r = 0.00 (Iterations 1b8b). Using the estimated hydraulic parameters from the final iteration, K(
) and D(
) are can be calculated using Eq. [3] and Eq. [5], respectively.
The sorptivity, S(
n,
o), can be estimated for any chosen value of
n and
o using the Philip and Knight (1974) solution in conjunction with the previously determined hydraulic parameters. Figure 5 maps lines of equal sorptivity as a function of inlet and initial water contents. S(
n,
o) varies from zero to infinity with a three order of magnitude variation within the region of common interest. The infinite S(
n,
o) predicted is a result of the infinite D(
) predicted by Eq. [5] as
o approaches
s. Of course, an infinite S(
n,
o) is physically infeasible. The results imply that this type of test is best performed at inlet water contents less than saturation.
 |
CONCLUSIONS
|
|---|
This research describes an improved method to effectively optimize hydraulic parameters using data from a one-dimension infiltration test. Data collection was performed using a computer controlled
-ray attenuation system to measure moisture content. A computer controlled syringe pump maintained the inlet boundary conditions. Using a BruceKlute test to measure a specific wetting profile, van GenuchtenMualem's description of D(
) is optimized using Philip and Knight's (1974) constant concentration solution to predict the previously measured wetting profile. The van Genuchten soil parameters defining D(
) in conjunction with an independently measured Ks allow for the estimation of K(
),
(
), and S(
n,
o). The procedure optimizes n and
independently, and a procedure is also provided to increase the accuracy of the predicted hydraulic parameters using a single measurement of
versus
, which enables the optimization of
r.
Received for publication May 14, 2002.
 |
REFERENCES
|
|---|
- Arya, L.M., and J.F. Paris. 1981. A physicoempirical model to predict the soil moisture characteristics from particle-size distribution and bulk density data. Soil Sci. Soc. Am. J. 45:10231030.[Abstract/Free Full Text]
- Boltzmann, L. 1894. Zur integration der diffusionsgleichung bei variabeln diffusions coefficienten. Ann. Phys. (Berlin) 53:959964.
- Brooks, R.H., and A.T. Corey. 1966. Properties of porous media affecting fluid flow. J. Irrig. Drain. Soc. Civ. Eng. Div. IR2:6187.
- Brown, G.O. 1987. Solute transport by a volatile solvent. Ph.D. diss., Dep. Civil Eng., Colorado State University, Fort Collins.
- Brown, G.O., and B. Allred. 1992. The performance of syringe pumps in unsaturated horizontal column experiments. Soil Sci. 154:243249.
- Brown, G.O. and D.B. McWhorter. 1990. Solute transport by a volatile solvent. J. Contam. Hydrol. 5:387402.
- Bruce, R.R., and A. Klute. 1956. The measurement of soil diffusivity. Soil Sci. Soc. Am. Proc. 20:458462.
- Burdine, N.T. 1953. Relative Permeability calculations from pore-size distribution data. J. AIME 198:7177.
- Carsel, R.F., and R.S. Parish. 1988. Developing joint probability distributions of soil and water retention characteristic. Water Resour. Res. 24:755769.
- Clothier, B.E., D.R. Scotter, and A.E. Green. 1983. Diffusivity and one-dimensional absorption experiments. Soil Sci. Soc. Am. J. 47:641644.[Abstract/Free Full Text]
- Gardner, W.H. 1986. Water Content. p. 493544. In A. Klute (ed.) Methods of soil analysis. Part 1. 2nd ed. Agron. Monogr. 9. SSSA, Madison, WI.
- Klute, A. 1952. Some theoretical aspects of the flow of water in unsaturated soils. Soil Sci. Soc. Proc. 16:144148.
- McBride, J.F., and R. Horton. 1985. An empirical function to describe measured water distributions from horizontal infiltration experiments. Water. Resour. Res. 21:15391544.
- McWhorter, D.B. 1971. Infiltration affected by flow of air. Hydrology Paper No. 49, Colorado State University, Fort Collins.
- Meyer, J.J., and A.W. Warrick. 1990. Analytical expression for soil water diffusivity derived from horizontal infiltration experiments. Soil Sci. Soc. Am. J. 54:15471552.[Abstract/Free Full Text]
- Mualem, Y. 1976. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12:513522.
- Philip, J.R. 1960. General method of exact solution of the concentration-dependent diffusion equation. Aust. J. Phys. 13:112.
- Philip, J.R. 1969. Theory of infiltration. Adv. Hydroscience 5:215296.
- Philip, J.R. 1973. On solving the unsaturated flow equation: 1. The flux concentration relation. Soil Sci. 116:328335.
- Philip, J.R., and J.H. Knight. 1974. On solving the unsaturated flow equation: 3. New quasi-analytical technique. Soil Sci. 117:113.
- Simunek, J., J.W. Hopmans, D.R. Nielson, and M.Th. Van Genuchten. 2000. Horizontal infiltration revisited using parameter estimation. Soil Sci. 165:708717.
- van Genuchten, M.Th. 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892898.[Abstract/Free Full Text]
- Warrick, A.W. 1994. Soil water diffusivity estimates from one-dimensional absorption experiments. Soil Sci. Soc. Am. J. 58:7277.[Abstract/Free Full Text]
- Warrick, A.W., D.O. Lomen, and S.R. Yates. 1985. A generalized solution to infiltration. Soil Sci. Soc. Am. J. 49:3438.[Abstract/Free Full Text]
This article has been cited by other articles:

|
 |

|
 |
 
D. Ma, Q. Wang, and M. Shao
Analytical Method for Estimating Soil Hydraulic Parameters from Horizontal Absorption
Soil Sci. Soc. Am. J.,
May 1, 2009;
73(3):
727 - 736.
[Abstract]
[Full Text]
[PDF]
|
 |
|