VZJ Download to Citation Manager
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
QUICK SEARCH:   [advanced]


     

Published online 26 May 2006
Published in Vadose Zone J 5:628-640 (2006)
DOI: 10.2136/vzj2005.0061
© 2006 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via HighWire
Right arrow Citing Articles via ISI Web of Science (3)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Schwärzel, K.
Right arrow Articles by van Genuchten, M. Th.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Schwärzel, K.
Right arrow Articles by van Genuchten, M. Th.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Schwärzel, K.
Right arrow Articles by van Genuchten, M. Th.
Related Collections
Right arrow Vadose Zone Processes and Chemical Transport
Right arrow Soil Hydrology

ORIGINAL RESEARCH

Estimation of the Unsaturated Hydraulic Conductivity of Peat Soils

Laboratory versus Field Data

K. Schwärzela,*, J. Simunekb, H. Stoffregenc, G. Wessolekc and M. Th. van Genuchtend

a Institute of Soil Science and Site Ecology, Univ. of Technology, Dresden, Germany
b Dep. of Environmental Sciences, Univ. of California, Riverside, CA 92521 USA
c Institute of Ecology, Dep. of Soil Science and Soil Protection, Technical Univ. of Berlin, Germany
d USDA-ARS, George E. Brown, Jr. Salinity Lab., Riverside, CA 92521, USA
* Corresponding author (Kai.Schwaerzel{at}Forst.TU-Dresden.de)

Received 10 May 2005.



   ABSTRACT
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
INTRODUCTION
 THEORY
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
As compared with mineral soils, there are few in situ measurements of the unsaturated hydraulic properties of peat soils available. We used parameter estimation (inverse) methods to estimate the water retention and hydraulic conductivity functions of drained peat soils from both laboratory and field data. The laboratory data were obtained on small cores using the traditional evaporation method, while the field data were obtained by means of evaporation experiments using groundwater lysimeters with and without vegetation. The field experiments without vegetation produced highly uncertain parameters and were limited to a relatively small pressure head range. Better results were obtained for the lysimeter with a grass cover that caused the peat soil to dry out more. However, a physically realistic minimum in the objective function for the plant-covered lysimeter could be found only when prior information about several parameters was included in the optimization. Good agreement was obtained between the laboratory and field measurements. The hydraulic functions were subsequently tested by comparing forward simulations with independently measured pressure heads and water contents of an additional lysimeter experiment under grass. The dynamics of the drying process was described well using the optimized soil hydraulic properties.

Abbreviations: VGM, van Genuchten–Mualem • WRC, water retention curve


   INTRODUCTION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
THEORY
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
UNDERSTANDING the processes that control the retention and flow of water in peat soils is critical to effective management of such soils from both agricultural and ecological perspectives. In contrast to mineral soils, much less is known about the unsaturated soil hydraulic properties of peat soils, especially the hydraulic conductivity. This is in part due to the unique nature of the physical and hydraulic properties of peat soils, such as volume changes during dewatering (Schwärzel et al., 2002; Price, 2003). Only a relatively few measurements of the unsaturated hydraulic conductivity of peat soils have been reported in the literature (e.g., Rijtema, 1965; Renger et al., 1976; Loxham and Burghardt, 1986; Schouwenaars and Vink, 1992; Baird, 1997; Silins and Rothwell, 1998; Schlotzhauer and Price, 1999; Schindler et al., 2003; Naasz et al., 2005). Most of these studies involved direct laboratory measurements of the unsaturated hydraulic conductivity using steady-state methods on small core samples.

Many laboratory and field methods are currently available for direct determination of the soil water retention and unsaturated hydraulic conductivity functions (Dane and Topp, 2002). In addition to these direct methods, inverse solution techniques are now also increasingly used to determine the hydraulic properties from transient field or laboratory measurements (e.g., Dane and Hruska, 1983; Parker et al., 1985; Eching and Hopmans, 1993; Simunek et al., 1998b; Abbaspour et al., 2000; Jacques et al., 2002; Sonnleitner et al., 2003). The application of inverse methods assumes a certain parametric model for the hydraulic properties with as yet unknown parameters. A flow experiment is then typically conducted to estimate these parameters by minimization of deviations between measured and predicted flow variables (Kool and Parker, 1987). An excellent review of the principles and advantages of inverse methods was given by Hopmans et al. (2002). In this study we use inverse procedures to determine the unsaturated hydraulic properties of peat soils.

In general, laboratory measurements are quick and precise, but often lead to soil hydraulic properties that are not representative of field conditions. Differences between laboratory and field measurements have long been noted for mineral soils (e.g., Sonnleitner et al., 2003), but such differences may be especially important for peat soils because of their unique physical properties. Moreover, standard field methods such as the instantaneous profile method (Vachaud and Dane, 2002) or the plane-of-zero-flux method (Arya, 2002) are generally not applicable to peat soils, which are usually situated in areas with relatively shallow water tables. While the use of lysimeters may overcome some of these limitations, their installation and operation is often very expensive. In contrast to conventional lysimeters, we developed a relatively low-cost lysimeter that can be installed more easily, and with minimal technical effort (Schwärzel and Bohl, 2003). Field lysimeters of this type, without vegetation, were used in the Schwärzel and Bohl (2003) study to determine the unsaturated soil hydraulic functions of peat soils using traditional direct methods. Results showed that direct estimation of the hydraulic conductivities from the lysimeter data was very tedious, uncertain (due to very small pressure head gradients), and limited to a relatively small pressure head interval. In contrast, measurements of the hydraulic properties in the laboratory were quicker and more accurate over a much wider range of water contents. Still, questions arose as to what extent the relatively small size of the laboratory samples and soil disturbance during sample collection affected the hydraulic properties of the peat soils.

The objectives of this study hence were to use parameter inverse estimation methods to determine the hydraulic properties of peat soils and to investigate possible differences between laboratory and in situ field measurements. Transient evaporation experiments were for this purpose performed on small laboratory soil cores and using field lysimeters with and without a grass cover. An important aspect of this study was to investigate the suitability of transient evaporation experiments on lysimeters with a grass cover. The estimated hydraulic functions were subsequently tested by comparing predictions with measured pressure heads and water contents of an additional lysimeter experiment under grass.


   THEORY
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 THEORY
MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Governing Equations
Water flow simulations were conducted with the Hydrus-1D model (Simunek et al., 1998a), which numerically solves the Richards equation:

Formula 1[1]
where {theta} is the volumetric water content (L3 L–3), t is time (T), z is the vertical coordinate (positive upward) (L), h is the pressure head (L), K is the hydraulic conductivity (L T–1), and S represents the root water uptake rate (T–1). Solving Eq. [1] requires a set of initial and boundary conditions. The initial condition in our study was given in terms of measured pressure head, hi:

Formula 2[2]

An atmospheric boundary condition was specified at the soil surface:

Formula 3[3]
and

Formula 4[4]
where E is the maximum potential rate of infiltration or evaporation under the current atmospheric conditions (L T–1), and hA and hS are the minimum and maximum pressure heads, respectively, allowed at the soil surface (L). Following Simunek et al. (1998a) values for hA and hS were assumed to be –1.0 x 106 and 0 cm, respectively. We additionally applied a zero flux condition to the bottom boundary of the lysimeter:

Formula 5[5]
where L is the depth of the lysimeter.

Soil Hydraulic Properties
The unsaturated soil hydraulic properties were described with the van Genuchten–Mualem (VGM) model (van Genuchten, 1980; Mualem, 1976):

Formula 6[6]

Formula 7[7]

Formula 8[8]

Formula 9[9]
where {Theta} is effective saturation, {theta}R is the residual volumetric water content, {theta}S is the saturated volumetric water content, KS (L T–1) is the saturated hydraulic conductivity, and {alpha} (L–1), n, m, and {lambda} are empirical parameters.

Weiss et al. (1998) previously demonstrated the flexibility of the VGM model in describing the water retention data of peat soils. They obtained optimal results using only three VGM parameters ({theta}S, {alpha}, and n), with {theta}R being fixed at zero. In our preliminary studies we could not find improvements in predictions of K(h) by varying the {lambda} parameter in Eq. [8]. The use of alternative parametric expressions (Burdine, 1953; Brooks and Corey, 1964) also did not lead to better optimization results. For these reasons we reduced the number of unknown parameters in Eq. [6] by fixing the values of {theta}R and {lambda} to zero and 0.5, respectively, thus enhancing the likelihood of uniqueness and stability of the inverse solution. We note that in our inverse analysis we considered all parameters (including KS and {theta}S) to be fitting parameters that define the scale and the shape of the soil hydraulic functions without any clear physical meaning. Their mutual correlation and possible nonidentifiability hence were not a major concern to us, as long as they reproduced the measured data used for calibration and subsequent validation.

Root Water Uptake
In this study the sink term, S, in Eq. [1] is defined as the volume of water removed from a unit volume of soil per unit time due to plant water uptake. Feddes et al. (1978) defined S as

Formula 10[10]
where ß (0 ≤ ß ≤ 1) is a prescribed dimensionless function representing root water uptake in response to water stress, and Sp is the potential root water uptake rate distribution (T–1). When the potential root water uptake rate is uniformly distributed over the root zone, Sp becomes

Formula 11[11]
where Lr is the depth (L) of the root zone and Tp is the potential transpiration rate (L T–1). Equation [11] may be generalized by introducing a nonuniform distribution of the potential root water uptake rate over a root zone of arbitrary shape (Fig. 1 ):

Formula 12[12]
where b(z) is the normalized water uptake distribution (L–1). The spatial variation of the potential extraction term, Sp, over the root zone is obtained by normalizing any arbitrarily measured or prescribed root distribution function, b'(z), as follows

Formula 13[13]
where Lr is the region occupied by the root zone (L).


Figure 1
View larger version (40K):
[in this window]
[in a new window]
  		Fig. 1. Schematic of the potential root water uptake distribution function.

 
The potential root-water uptake function Sp (Eq. [11]) corresponds to the transpiration rate Tp of a growing canopy (Fig. 1) such that

Formula 14[14]
in which Tp depends on weather conditions, plant height, and canopy resistance against water loss. The actual root water uptake rate S is obtained by substituting Eq. [12] into Eq. [10], while the actual transpiration rate Ta (L T–1) follows by integrating S over the root zone:

Formula 15[15]

Inverse Parameter Estimation Procedure
Parameter estimation in this study was accomplished by minimizing the following objective function, {Phi} (Simunek et al., 1998a):

Formula 16[16]
in which the first term represents deviations between measured, qj*, and calculated, qj, space–time variables, such as pressure heads or water contents. In this term mq is the number of different sets of measurements, nqi the number of measurements in a particular measurement set, qj*(z, ti) represents a specific measurement at time ti for the jth measurement set, qj(z,ti,b) are the corresponding model predictions for the vector of estimated parameters b({theta}R, {theta}S, {alpha}, n, KS, {lambda}), and {nu}i and wij are weights associated with a particular measurement set or data point, respectively. The second term of Eq. [16] represents differences between independently measured, pj*, and predicted, pj, soil hydraulic properties, such as retention or conductivity data at particular pressure heads. The terms variables, mp, npj, pj*({theta}i), pj({theta}i,b), {upsilon}'j and w'ij, in this second term have similar meanings as in the first term, but now for the soil hydraulic properties. The last term of Eq. [16] represents a penalty function for deviations between prior knowledge of the soil hydraulic parameters, bj*, and their final estimates, bj, with nj being the number of parameters with prior knowledge and Formula 16j representing preassigned weights.

In this study we assume that the weighting coefficients wij and w'ij in Eq. [16] are equal to 1. The weighting coefficients {nu}i and {upsilon}'j that equalize the contribution of different data types to the objective function {Phi} are given by

Formula 17[17]
which normalizes the objective function in terms of the variances {sigma}j2 of measurement type j.

Minimization of the objective function {Phi} was accomplished using the Levenberg–Marquardt nonlinear optimization method (Marquardt, 1963; Simunek et al., 1998a). Since the Levenberg–Marquardt method is a local search method, the final optimized parameters may depend on the initial parameter values (e.g., Abbaspour et al., 2004). The use of a well-defined flow problem combined with prior information about the optimized parameters generally enhances the likelihood of uniqueness of the inverse solution (Russo et al., 1991). Although uniqueness of the inverse problem cannot always be guaranteed, solving the problem several times with different initial parameter estimates should ensure that the search algorithm locates the global minimum of the objective function.

Compared with global optimization methods (e.g., the genetic algorithm method used by Vrugt et al., 2001) in conjunction with Hydrus) that have been used recently to determine empirical parameters of lumped hydrological (Vrugt et al., 2002) or other black box models (e.g., Barhen et al., 1997), the Levenberg–Marquardt method provides useful statistical information about the optimized parameters, such as their mutual correlation and confidence intervals. This information should be indicative whether or not the experimental data have enough information content for unique determination of the optimized parameters (Simunek and Hopmans, 2002).

Forward Simulations
Based on an additional lysimeter experiment under grass, forward simulations were conducted to test the accuracy of estimated hydraulic functions. Differences between various simulations were quantified by comparing predicted results Pi with the observed data Oi using the RMSE given by (Loague and Green, 1991):

Formula 18[18]
where N is the number of values and O is the arithmetic average of the N measurements.


   MATERIALS AND METHODS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 THEORY
 MATERIALS AND METHODS
RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
The Rhinluch Field Site
Field measurements were performed in Rhinluch, a fen area of about 87 000 ha, 60 km northwest of Berlin, Germany. Peat soils in the area, drained more than 250 yr ago, are characterized as a Eutri Folic Histosol, with an average thickness of about 120 cm. Below the peats are glacifluvial sands (mostly fine sand) and limnic sediments such as detritus or calcareous mud. The upper peat layers are strongly decomposed and altered pedogenically because of intensive land use and drainage. Humified and very humified peats are mostly present near the soil surface. Deeper layers consist primarily of sedge (Carex spp.) or reed-peats (Phragmites spp.), or of a mixture of these two peats. Vegetation at the study site was dominated by reed canary grass (Phalaris arundinacea L.).

Laboratory Methods
The unsaturated soil hydraulic properties were determined in the laboratory during transient conditions using the evaporation method (Wendroth et al., 1993) and during steady state conditions using the hanging water column method and a pressure apparatus (Dane and Topp, 2002). Undisturbed cores (3 replicates from each layer) with a height of 10 cm and an inside diameter of 5.64 cm were used for the evaporation method. Five tensiometers with cups of 4 cm length and 0.3 cm outside diameter were installed at depths of 1, 3, 5, 7, and 9 cm and connected to pressure transducers with a precision of ±1 cm. Once placed on a ceramic plate, a hanging water column (with deionized water) was used to create a negative pressure equal to the initial soil water pressure head of –20 cm at the bottom of the cores. After hydraulic equilibrium was reached, the lid-covered samples were moved from the ceramic onto impermeable plates for the evaporation experiments. After an additional day of equilibration, immediately before the evaporation experiments, tensiometer readings were compared and corrected, assuming that hydraulic equilibrium had been reached. The evaporation process was then started, with pressure heads recorded every 30 min. Water losses with time were measured by weighing the cores several times each day. We also estimated water losses from water contents determined with five mini TDR probes (Easy Test, Lublin, Poland, 1993) consisting of two wires (5.3 cm long and spaced 0.6 cm apart) installed at depths of 1, 3, 5, 7, and 9 cm. The error of the laboratory TDR probes for measuring water content was about 0.01 m3 m–3 (Plagge et al., 1994). At the end of the experiments, tensiometers and TDR probes were removed and residual water contents determined by oven drying the samples at 105°C and weighing.

The water content was determined by using a site-specific calibration between the dielectric constant, e, and the water content (standard error of calibration = 0.022 m3 m–3, Schwärzel and Bohl, 2003):

Formula 19[19]

Unsaturated hydraulic conductivities were estimated from the evaporation experiments in two ways. First, when TDR measurements were available, water fluxes were calculated from the water content changes and pressure head gradients. For these calculations, water contents for each TDR probe were fitted with a quadratic polynomial function of time, and measured pressure heads with a quadratic polynomial function of depth at each measurement time. Hydraulic conductivities as a function of the pressure head were next calculated using Darcy's Law. Second, when TDR measurements were not available, hydraulic conductivities were calculated using the approach of Wind (1968). In this approach, water contents were first calculated from measured pressure heads using a retention curve (the VGM model in our case) that was optimized such that water content changes in the soil column became equal to the measured total water losses. These water contents were subsequently used to calculate hydraulic conductivities as described above.

In addition to the above direct methods for estimating the soil hydraulic functions, an inverse parameter optimization method (Simunek et al., 1998b) was applied to all laboratory measurements. The tensiometer readings as a function of time and the total water volume at the end of the experiment were used for the calibration. Water volume measurements were used to calculate evaporation rates needed for the upper boundary condition. A zero-flux condition was used as the bottom boundary condition. As the initial condition we assumed hydraulic equilibrium with the initially imposed pressure head (–20 cm) at the bottom of the soil sample. No prior information was used in the inverse analysis of the laboratory evaporation data.

Undisturbed soil core samples (5 replicates from each layer) with a height of 4 cm and a diameter of 5.5 cm were used to determine retention curves assuming steady-state conditions. Samples were first placed on a ceramic plate and saturated. The dewatering process (desorption) was performed using ceramic plates connected to a hanging water column down to a pressure of –100 cm. A pressure apparatus was used below this pressure. Volume changes during the dewatering process were quantified with a caliper rule at the end of every pressure step. The peat soils used in our study started to noticeably shrink at pressure heads below –30 m (5% loss of volume). More details about this problem were given by Schwärzel et al. (2002).

Field Methods
Three nonweighing groundwater table lysimeters (1 by 1 by 1 m) were installed in two replicates and equipped with fully automated tensiometers using pressure transducers with a precision of ±1 cm and manually operated TDR probes (Easy Test, Lublin, Poland, 1993), made from two wires (10 cm long and spaced 1.6 cm apart, and having measurement precision of 0.01 m3 m–3), at depths of 10, 20, 30, 40, 50, 60, 70, and 80 cm. All sensors were horizontally installed. TDR readings were converted into water content values using the site-specific calibration equation (Eq. [19]). Groundwater levels in the lysimeters and the surrounding soil were monitored using slotted PVC pipes connected via a flexible tube to a differential pressure transducer a precision of 1 cm. Following the approach of Feddes (1971), the water level in the lysimeter was either adjusted to the surrounding groundwater level (the reference level) or fixed at a certain depth. Water level differences recorded by the control unit either triggered a bilge pump in the lysimeter pipe, or released a valve in a water supply container. The volume of water automatically added or removed from the lysimeter was manually measured by recording water levels in the container. Hysteresis in the imposed groundwater level was observed to be at most 2 cm.

The average actual evapotranspiration rate, ETa, (L T–1) for a certain time interval between t1 and t2 was calculated as follows:

Formula 20[20]
where P (L T–1) is the precipitation rate, Q (L T–1) is the net inflow rate (inflow minus outflow) into the lysimeter needed to keep the groundwater table at the same level, and L (L) is the depth of the lysimeter. Additional information about the groundwater lysimeters and the field site was given in Schwärzel and Bohl (2003).

Two of the three lysimeters were operated with a grass cover. One of these was used for parameter estimation and the other one to test the estimated parameter in forward simulations. The third lysimeter was operated in 1994 and 1995 without plants and with roofing for rain protection. This third lysimeter was fully saturated at the beginning of the vegetation period, at which time the water level control mechanism was disabled by disconnecting the control unit. As a result, the soil monolith slowly dried out in about 70 d, similarly to the laboratory evaporation experiment. Bare lysimeter data were analyzed using the standard evaporation method (Hillel et al., 1972; Arya, 2002), as well as inverse procedures. For the latter approach, the evaporation rate, q0, (L T–1) from the soil surface was calculated from TDR-measured soil water content changes in the soil profile as follows:

Formula 21[21]
where W (L) is the total depth of water in the soil profile given by

Formula 22[22]

Experimental Limitations
The accuracy of the directly calculated unsaturated hydraulic conductivities and the inversely estimated hydraulic parameters depends on the precision of the measurement devices (i.e., on the instrumental error). The calculated unsaturated hydraulic conductivities are especially sensitive to hydraulic pressure gradients close to –1 cm cm–1. The influence of the measurement precision of the tensiometers on the calculated hydraulic gradients, {partial}h/{partial}z (L L–1), and consequently on the unsaturated hydraulic conductivities, were quantified as follows (Stoffregen, 1998):

Formula 23[23]
where {delta}h is the measurement precision of the pressure transducer (±1 cm), and h1 and h2 are measured pressure heads at depths of z1 and z2, respectively. Upper and lower limits of the calculated hydraulic conductivity at a particular pressure gradient may be determined using Eq. [21] and [23]. Figure 2 shows the ratios of the calculated unsaturated hydraulic conductivity (K) without considering the measurement error, to either the upper limit (Kmax) or the lower limit (Kmin) of the unsaturated hydraulic conductivities with the measurement error. The actual value of the hydraulic conductivity should be located between the upper and lower curves. Figure 2 demonstrates that tensiometer measurement errors can have a large effect on the calculated hydraulic conductivities, especially when small hydraulic gradients are present.


Figure 2
View larger version (12K):
[in this window]
[in a new window]
  		Fig. 2. Effect of tensiometer measurement precision on calculated unsaturated hydraulic conductivities.

 
Another source of error may be due to imprecise estimates of the evaporation rate as calculated from the TDR readings. The influence of this error on the accuracy of the unsaturated hydraulic conductivity was quantified as follows (Stoffregen, 1998):

Formula 24[24]
where {delta}q0 is the measurement error in the calculated evaporation rate. Equation [24] shows that a relative error in the calculated flux leads to the same relative error in calculated unsaturated hydraulic conductivity. The smaller the water content change, the higher the relative error will be in the evaporative flux, and hence the higher the errors in the hydraulic conductivity.


   RESULTS AND DISCUSSION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 THEORY
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
CONCLUSIONS
 REFERENCES
 
Laboratory Measurements
Directly measured and inversely estimated laboratory water retention and hydraulic conductivity functions of the peat soil are compared in Fig. 3 . Final optimized parameters obtained from the inverse approach are listed in Table 1. Since unsaturated hydraulic conductivities were calculated directly at various measurement locations within the soil core, the analysis leads to several data sets, as shown in Fig. 3. By comparison, the inverse method leads to only one effective conductivity function for the entire column. Notice the relatively good agreement in Fig. 3 between the directly and inversely estimated hydraulic conductivities for pressure heads less than –100 cm. Deviations between the direct and inverse estimates at the higher pressure heads are most likely due to the very small pressure head gradients during the early stages of the experiments (as discussed above). These small pressure gradients caused relatively high uncertainty in the directly estimated hydraulic conductivities near saturation (Simunek et al., 1998b) (Fig. 2).


Figure 3
View larger version (27K):
[in this window]
[in a new window]
  		Fig. 3. Hydraulic conductivity (left) and water retention (right) functions for the peat soil as determined in the laboratory (3 replicates for each horizon) using inverse parameter estimation and direct methods. Open circles represent directly estimated hydraulic conductivities, squares are directly measured retention data, {theta}(h); solid lines represent inversely estimated K(h) and {theta}(h) functions; and dashed lines are {theta}(h) functions fitted to directly measured retention data and K(h) functions predicted from the fitted retention curves and the directly measured KS.

 

View this table:
[in this window]
[in a new window]
  		Table 1. Inverse parameter estimation results for the laboratory experiments. "No." indicates replicate number (soil core), while {Phi}h and {Phi}{theta} are the minimized values of the pressure head and water content terms of the objective function, respectively.{dagger}

 
Despite the high correlation between the measured and simulated flow variables (Table 1), much uncertainty was found for the parameters KS and {alpha} because of their considerable mutual correlation (0.959–0.991). The uncertainties in KS and {alpha} were caused mainly by small pressure head gradients at the beginning of the transient evaporation experiment. Our findings are consistent with those by Romano and Santini (1999), who performed numerical evaporation experiments to explore the sensitivity of pressure heads to selected soil hydraulic parameters. Romano and Santini (1999) concluded that measurements of the pressure head later during the experiments produce steeper hydraulic gradients and hence better identifiable soil hydraulic parameters.

Variability in hydraulic conductivities was found to be relatively small within each of the three soil layers, while being in a relatively narrow range at intermediate and high tensions. The K(h) curves of the various horizons deviated only in the near-saturated region. Major differences in conductivities among the different horizons occurred for pressure heads larger than about –100 cm. Unsaturated hydraulic conductivities near saturation decreased with increasing depth. Conductivities of the third horizon were up to an order of magnitude smaller than those of the first horizon. This may have been due to the use of heavy machinery during tillage, leading to compaction of layers below the root zone, as well as due to natural peat compression often observed for peat soil profiles (e.g., Price, 2003). Figure 3 shows that the reduced unsaturated hydraulic conductivities correspond with lower saturated water contents and that the compacted subsurface layers exhibit more variability in both {theta}(h) and K(h). Figure 3 also shows that soil compaction affects the shape of the soil water retention curve (WRC), with the measured curve acquiring a bimodal-type distribution especially at the 15- to 25-cm soil layer. The retention curve of this layer could not be described with the standard VGM model.

Results in the right column of Fig. 3 also show that within the measurement range of the tensiometers, retention data points obtained from the steady-state measurements agree reasonably well with retention curves obtained using the inverse parameter estimation technique. However, differences between the directly measured retention data and inversely estimated retention curves increased substantially at lower water contents. Extrapolation of the optimized WRCs beyond the measurement range covered by the experiment hence produced, as expected, a high level of uncertainty. Some overestimation of the directly measured retention data at lower pressure heads may have been caused by shrinkage during dewatering of the peat soil. Schwärzel et al. (2002) compared WRCs of slightly decomposed peat soils with and without considering shrinkage. They found substantially higher volumetric water contents at low (more negative) pressure heads when shrinkage was taken into account.

The independently measured retention data points fitted with the van Genuchten (1980) analytical model and the KS value measured in the laboratory with the constant head method were used to predict the unsaturated hydraulic conductivity function (Fig. 3). Even though the shapes of the K(h) curves were quite similar, the predicted K(h) curves (dashed lines in Fig. 3) did not correspond well with the other (direct and inverse) results, mostly because of a KS value that was too high. Good agreement between measured and estimated hydraulic conductivity curves was found only for the top layer.

The VGM parameters listed in Table 1 differ from those of Weiss et al. (1998) and Naasz et al. (2005). This is due to the fact that our peat soils were much more influenced by drainage processes, and hence pedogenetic changes, than the peat material of the other studies.

Field Measurements
Only very small changes in pressure heads and water contents were recorded during the later stages of the evaporation experiment for the bare lysimeter (see also Schwärzel and Bohl, 2003). The topsoil layer at that time was air dry with a very low hydraulic conductivity, which prevented further evaporation. Actual soil evaporation rates were calculated from weekly measured water content values. While pressure heads were measured hourly, only daily averages were used in the inverse parameter optimization. Figure 4 shows that the measured pressure heads and water contents of the bare lysimeters decreased substantially only in the upper parts of the profile in 1994 (top plots). Similar results were also found for 1995 (not shown).


Figure 4
View larger version (35K):
[in this window]
[in a new window]
  		Fig. 4. Measured (open circles) and fitted (dashed lines) pressure heads (left) and water contents (right) for the lysimeter experiments without (top) and with (bottom) plants.

 
The evaporation experiments of the bare lysimeter were first analyzed using the direct approach. Because water contents varied very little, we could only calculate unsaturated hydraulic conductivities of the first two major soil horizons between 10 and 20 cm and between 20 and 30 cm. For 1994 and 1995 we obtained 61 and 54 K(h) data pairs, respectively (Fig. 5). Because of very small hydraulic gradients and flow rates, the calculated values are quite scattered, causing much uncertainty in the directly determined K(h) data.


Figure 5
View larger version (21K):
[in this window]
[in a new window]
  		Fig. 5. Hydraulic conductivity (left) and water retention (right) functions of the peat soil as determined from the lysimeter experiments using inverse parameter estimation and direct methods.

 
To obtain more information about the hydraulic properties in the dry range, we used another lysimeter covered with vegetation. This lysimeter was operated during the unusually dry summer of 1994. Although at the beginning of the experiment the groundwater level in the lysimeter was 60 cm below the soil surface, after 2 wk the lysimeter was already without groundwater (below the measurement depth) because of the high evaporative demand. Figure 4 (bottom plots) shows that the decrease in pressure heads and water contents, especially in the two surface layers, was much more significant than for the bare lysimeter. The experiment with vegetation produced three weekly measured evapotranspiration rates that were scaled to 22 daily values using the calculated potential evapotranspiration, 24 water content values (again measured weekly), and 150 pressure head values (measured hourly and used as daily averages in numerical inversions) for use with the inverse method.

We next used inverse methods to determine the soil hydraulic parameters from both the bare and vegetated lysimeters. The soil profile in the numerical calculations was for this purpose divided into two different soil layers (0–30 cm and below 30 cm). A zero flux condition was imposed at the bottom of the soil profile in all cases. Measured evaporation fluxes were used as the upper boundary condition for the bare lysimeter, while evaporation from the lysimeters with plants was considered to be negligible because of the presence of complete plant cover. For the lysimeters with plants, root water uptake was calculated from weekly measured evapotranspiration rates, which were scaled to daily values using calculated daily grass reference evaporation rates (Allen et al., 1994). The initial condition was given in terms of measured pressure heads for all runs.

The VGM soil hydraulic parameters of the two soil horizons were fitted simultaneously to all available data. The residual water content {theta}R and the pore-connectivity exponent {lambda} in Eq. [8] were fixed at zero and 0.5, respectively. Hence, only eight parameters were optimized simultaneously for the two horizons. Three sets of initial estimates of the parameters were used to identify the global minimum.

Figure 4 compares measured and fitted tensiometer and TDR readings for the lysimeters with and without plants. We obtained an excellent fit of the tensiometer readings in both cases and for all depths. Deviations between measured and fitted TDR readings were somewhat higher, especially at depths of 10 and 20 cm, for the plant-covered lysimeter. Notice that the measured water contents at the 10-cm depth were slightly overestimated and at the 20-cm depth underestimated. We believe that these deviations were most likely due to the use of a simplified potential root water uptake distribution function with a constant rooting depth of 30 cm. A better fit to the measured water contents may have been possible using a root uptake function that accounts for root growth into the subsoil and increased root intensity in the topsoil. We did not pursue this possibility because of a lack of independent root growth data.

Final results of the parameter optimization for both lysimeters and both years are presented in Table 2. The different sets of initial parameter estimates always generated the same set of optimized parameters for the topsoil of the bare lysimeter. Also, only small differences were found between the parameter values of the bare lysimeter for years 1994 and 1995. Despite the good fit of data and the small differences between parameter values for years 1994 and 1995, we found a high level of uncertainty for the parameters of the bare lysimeter. We always obtained strong correlations between the parameters n and {alpha} (–0.990 to –0.992), n and KS (–0.989 to –0.992), and {alpha} and KS (0.990–0.992). These strong correlations are a consequence of the much narrower range of measured pressure heads and water contents and suggest that the parameters {alpha}, n, and KS should not be identified simultaneously in such cases.


View this table:
[in this window]
[in a new window]
  		Table 2. Parameter estimation results obtained with three different initial estimates of the optimized parameters (Runs 1–3) for the field lysimeter experiments in 1994 and 1995 with and without plants. {Phi}h and {Phi}{theta} are the minimized values of the pressure head and water content terms of the objective function, respectively.

 
The inverse simulations for the subsoil of the lysimeter without plants did not always converge to the same optimized parameters, which suggests that local minima were present in the objective function {Phi} (Table 2). This difficulty occurred for both 1994 and 1995. These relatively poor results reflect the fact that soil moisture variations in the subsoil were very small and hence did not provide enough resolution to permit a unique parameter optimization process. The parameter n, in particular, produced a wide range of values when different sets of initial parameter values were used. As opposed to the topsoil, estimated parameters for the subsoil were only weakly correlated to each other. The large confidence intervals of the estimated parameters indicate that the optimized parameters are highly uncertain (data not shown). Despite the large confidence intervals for the soil hydraulic parameters of the bare lysimeter, differences between the unsaturated hydraulic conductivity functions obtained from the different experiments were relatively small. This is shown in Fig. 5 , which compares the hydraulic conductivity and water retention functions of the topsoil for all successful optimizations.

Tables 1 and 2 indicate relatively good agreement between the laboratory and field estimates of the parameters {alpha} and n. However, the saturated hydraulic conductivity, KS, estimated from the field data differed significantly from the laboratory estimates. This is attributed to the very narrow range of the field data as compared to the laboratory data, particularly a lack of information to define the field hydraulic properties at or near saturation.

For the plant-covered lysimeter, a minimum in the objective function was found only when prior information about several subsoil parameters, estimated from the independent laboratory water retention measurements and from field borehole tests (Bohl et al., 1994), was incorporated into the optimization process (Table 2). One reason for this is the fact that the soil was already slightly unsaturated at the beginning of the experiment, and hence little or no information was available about the near-saturated hydraulic properties. Consequently, the global minimum of the objective function was difficult or impossible to locate when different initial parameter values were used. Our findings are consistent with those of Russo et al. (1991), who reported that uniqueness and stability of the inverse solution often depend on the quality of prior information about the model parameters. Sonnleitner et al. (2003) similarly found ill-defined hydraulic properties close to saturation when full saturation was not reached at the beginning of their lysimeter experiments.

In contrast to the hydraulic conductivity functions, the water retention functions obtained from the different lysimeter experiments were similar only within the measurement range (Fig. 5). With increasing pressure heads, differences between the different retention functions increased significantly. This shows that the optimized soil hydraulic functions should be used only for the conditions and measurement range of the experiment. Extrapolation of results, especially toward lower pressure heads (drier conditions), could quickly lead to incorrect predictions.

Contrary to the transient evaporation experiments conducted in the laboratory, application of the inverse method to the transient field evaporation data was subject to several experimental difficulties. Nonuniqueness and nonidentifiability problems were especially evident for the evaporation experiment performed on the bare lysimeter. This was primarily due to the limited range in the experimental pressure head and water content data. Soil moisture and pressure head variations were simply too small to produce a well-defined global minimum in the objective function (Fig. 4).

An important source of error in field studies may result when the evaporation flux is estimated from TDR readings. Flow rates calculated from the water content data of the bare lysimeter were very small and well within the range of the measurement precision of the TDR probes (Fig. 4). Romano and Santini (1999) performed numerical evaporation experiments to explore the influence of the imposed evaporation rate on ill-posedness of the inverse problem. They suggested that the inverse method for fine-textured soils subject to higher evaporation rates might experience problems of ill-posedness. For coarse-textured soils higher evaporation rates did lead to better identifiability of the unknown parameters. Peat soils behave more similarly to coarse- rather than to fine-textured soils. Higher evaporation rates hence should produce a wider range in water contents and pressure heads. The wider measurement range can be achieved also when the lysimeter is covered with grass. However, to find the global minimum of the objective function, prior information for the subsoil still needed to be incorporated in the objective function. A well-defined global minimum without prior information is much more likely when the water content data show more resolution, especially in the topsoil. Figure 4 illustrates that much more water is lost to evaporation from the upper part of the soil profile, with soil water depletion being very low deeper in the profile. Additional TDR probes in the upper part of the profile hence would have provided more useful information for the inverse estimation procedure. This would have enabled better definition of the evaporation rate needed as an upper boundary condition for the inverse simulations.

Forward Predictions
The accuracy of the estimated hydraulic parameters was tested by simulating pressure heads and water contents at different depths for a second plant-covered lysimeter during an extended drying period. The groundwater level in the lysimeter for this experiment was adjusted to the water table of the surrounding area. Within 28 d of the experiment the groundwater level had declined from 30 to 70 cm below the soil surface. Water contents during the experiment at the 10-cm depth decreased from 0.70 to 0.29, at the 30-cm depth from 0.80 to 0.75, and at the 60-cm depth from 0.87 to 0.83. Pressure heads similarly decreased from –20 to below –700 cm at the 10-cm depth (i.e., below the measurement limit of the tensiometers), and from 32 to –12 cm at the 60-cm depth. Forward simulations were conducted with both the directly determined and inversely estimated soil hydraulic parameter sets. Table 3 gives details about the different simulation scenarios and invoked parameter sets.


View this table:
[in this window]
[in a new window]
  		Table 3. Various scenarios and associated parameter sets used in forward calculations of the additional field lysimeter experiment with vegetation.{dagger}

 
Tables 4 and 5 summarize the RMSE between predicted and observed variables for different scenarios. Figure 6 presents time series of simulated and measured pressure head and water content values. Simulations with parameter values obtained from the various inverse estimates reproduced the measured pressure heads reasonably well. Relatively large differences between observed and predicted values were found only for the 20-cm depth. Inspection of Fig. 6 shows that the WRC parameter set (Table 3) produced substantially lower (negative) pressure heads than the 30- and 40-cm depth data, and much higher water contents than the 10-cm depth data when the peat soil dried out. The RMSE for the WRC set was in most cases the highest as compared with the other scenarios (Table 4), except for the depth of 20 cm. The best agreement between measured and predicted pressure heads at all depths was found with the parameter set determined under field conditions (see "Field" in Tables 4 and 5).


View this table:
[in this window]
[in a new window]
  		Table 4. The RMSE values in percentage of observed versus predicted pressure heads at different depths for the various forward simulations listed in Table 3.

 

View this table:
[in this window]
[in a new window]
  		Table 5. The RMSE values in percentage of observed versus predicted water contents at different depths for the various forward simulations listed in Table 3.

 

Figure 6
View larger version (24K):
[in this window]
[in a new window]
  		Fig. 6. Observed (circles) and predicted (lines) pressure heads (left) and water contents (right) for the additional evaporation experiment using a lysimeter with plants. The different simulation scenarios are explained in Table 3.

 
The hydraulic properties for the Lab 1, Lab 2, and Lab 3 scenarios were obtained from transient laboratory measurements with three replicates. Only minor differences between predicted pressure heads were observed when different laboratory parameter values were used. An exception is the 20-cm depth, where the simulated pressure heads for the three scenarios diverged more from the data than at other depths. This is mostly likely due to the variability in the hydraulic properties of the underlying layer (Fig. 3, 22–32 cm).

Simulations with different parameter sets also matched the dynamics of water contents reasonably well, especially for the 10-cm depth. Only the WRC scenario for the 10-cm depth showed relatively large deviations during the drier conditions, while the other parameter sets closely matched the drying process. Large differences between observed and predicted water content values were found only at the 20-cm depth. The poor agreement between the predicted and observed water contents at the 20-cm depth could have been due to the use of inaccurate root distribution and/or root growth functions. In general, simulations with different parameter sets led to similar values of the water content. Table 5 shows that the RMSE values were about the same for all cases and all depths.


   CONCLUSIONS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 THEORY
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
REFERENCES
 
Inverse analysis of the unsaturated hydraulic parameters of the peat soils using both field and laboratory experiments produced reliable results. Compared with the field, measurements in the laboratory allowed for quicker and more cost-effective estimation of the soil hydraulic properties across a much wider range of soil water contents. Considerable parameter uncertainty for the bare field lysimeters was due to the narrow range in soil moisture conditions of that lysimeter. This uncertainty was reduced substantially by using a vegetative cover to produce a much broader range of pressure heads and water contents. Prediction of the unsaturated hydraulic conductivity function from retention data and the saturated hydraulic conductivity, KS, produced unsatisfactory results, mostly because of estimates of KS that were too high.

The accuracy of the estimated hydraulic functions was further tested by comparing forward predictions using data from an additional evapotranspiration experiment on a different lysimeter. The dynamics of the drying process were very well described using hydraulic properties determined from the transient field and laboratory experiments. In contrast, predictions obtained from measured retention data and KS overestimated the drying process. Our results showed only little difference between the laboratory and field data. Hence, for peat soils of the type used in our study we recommend measurements of the unsaturated hydraulic conductivity in the laboratory in combination with inverse parameter optimization.


   REFERENCES
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 THEORY
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES