J. Jpn. Soc. Soil Phys.
土壌の物理性No. 125, p.3 ∼ 15 (2013)

Estimating the unsaturated hydraulic conductivity of Andisols
using the evaporation method
RUDIYANTO1 , Nobuo TORIDE2 , Masaru SAKAI2 and Martinus Th. van GENUCHTEN3

Abstract: Parameters of the bimodal van Genuchten (VG)
hydraulic functions for two aggregated Andisols were inversely determined using the evaporation method. Initial
estimates of the water retention parameters were determined from separate retention measurements, which facilitated rapid convergence of the parameter optimization
process regardless of the number of optimized parameters. When the bimodal water retention parameters were
fixed according to the independently measured retention
data from near saturation to very low pressure heads down
to −105 cm, it was possible to estimate the unsaturated
hydraulic conductivity, K(h), by optimizing only two conductivity parameters (Ks , ℓ). Since the flat region of the bimodal retention curve at intermediate pressures is difficult
to measure precisely, however, we still recommend optimizing all bimodal VG parameters to yield the best overall
results. Including water retention data at very low pressure
heads in the dry range extended the applicable range of themodel predictions, at least down to pressure heads of approximately −104 cm.
Key Words : unsaturated hydraulic conductivity, water retention curve, aggregated soil, Andisol, evaporation
method

1. Introduction
Andisols generally are developed from volcanic ash
consisting of noncrystalline materials such as allophone,imogolite, Al-humus complexes, and ferrihydrite (Shoji etal., 1993; Nanzyo, 2002). Andisols cover 17 % of the land
surface in Japan, where they are widely used for agriculture. Water flow and solute transport processes in Andisols
are of considerable interest because of their unique physical and chemical properties. They typically have a very
low bulk density and a well-developed aggregated structure. Because of this, Andisols usually exhibit a composite (stepwise) water retention function reflecting distinct
1 Dept. of Civil and Environmental Engineering, Bogor Agricultural

Univ., Indonesia.
2 Graduate School of Bioresources, Mie Univ., Tsu, Mie 514-8507, Japan.
Corresponding authour:取出伸夫, 三重大学大学院生物資源学研究科
3 Dept. of Mechanical Engineering, Federal University of Rio de Janeiro,
UFRJ, Rio de Janeiro, RJ, Brazil.
2013 年 9 月 3 日受稿 2013 年 9 月 20 日受理

but interacting inter-aggregate and intra-aggregate pore regions (Miyamoto et al., 2003).
Several bimodal or multimodal functions have been proposed over the years to account for the additive effects ofinter- and intra-aggregate pore regions to the overall soil
hydraulic properties (Peters and Klavetter, 1988; Othmer
et al., 1991; Gerke and van Genuchten, 1993, Mallants
et al., 1997; Mohanty et al., 1997). One frequently used
multimodal formulation stems from Durner (1992, 1994)
who developed a composite retention function by summing
multiple van Genuchten (VG) models (van Genuchten,1980). The bimodal form of this function has been usedin several recent studies (e.g., Coppola, 2000; Peters andDurner, 2008; Schelle et al., 2010), including for thewater retention properties of Andisols (Miyamoto et al.,
2003; Hamamoto et al., 2009; Chamindu Deepagoda et al.,2012). Most of these studies were concerned with the water retention properties of soils. By comparison, very few
studies have applied the bimodal VG model to the unsaturated hydraulic conductivity function of aggregated soils.
When used also for the hydraulic conductivity, most applications were concerned with improved descriptions of
the unsaturated conductivity function near saturation in attempts to account for the effects of macropores or rockfractures (Peters and Klavetter, 1988; Zurmühl and Durner,
1998; Iden and Durner, 2007; Durner and Iden, 2011).
The bimodal VG model has considerable flexibility indescribing the hydraulic properties of aggregated media.
Unfortunately, optimization of the large number of parameters in the model against transient flow data (such as frommultistep outflow or evaporation methods) is a major challenge. Since several of the hydraulic parameters are often
correlated (Zurmühl and Durner, 1998), it is inherently difficult to find the global minimum of the objective functionin the optimization process. Zurmühl and Durner (1998)
showed that convergence of the parameter optimization
process depends very much on the initial estimates of theunknown parameters. In addition to the transient measurements, the objective function could also include independently measured soil water retention or unsaturated
hydraulic conductivity data points. Parameter uncertainty

4

土壌の物理性第 125 号(2013)

pressure heads. The validity of the estimated unsaturated
hydraulic conductivity beyond the pressure measurement
range, however, was not discussed in detail in their study.

Fig. 1 Water retention curves for the Kumamoto and Mie Andisols fitted with the bimodal VG model (solid line). Retention
data were obtained with a hanging water column (open circles),a pressure plate (crosses), and a WP4 dewpoint potentiameter(open squares).

(i.e., the reliable range of the optimized parameter values)
generally will be reduced by including independently measured soil hydraulic property data in the objective function,
but this may affect the goodness of fit between modeled
and experimental data (Šimůnek et al., 1998; Hopmans etal., 2002). More rapid minimization of the objective function may be achieved by using as initial estimates water
retention parameters fitted to observed retention data. For
example, Spohrer et al. (2006) succeeded in estimating 25
parameters of the bimodal VG model for four soil layersof a tropical Acrisol from transient field data using initialvalues of the fitted water retention parameters.
When the hydraulic properties are estimated inversely
using the evaporation method, the optimized hydraulicfunctions are generally assumed to be representative only
within the range of the tensiometer data. No guarantee
exists that the hydraulic functions can be extrapolated beyond the invoked pressure measurement range. Šimůnek
et al. (1998) and Hopmans et al. (2002) suggested the applicable range of the model predictions could be extended
by including independently measured hydraulic data in theobjective function. Sakai and Toride (2007b) estimatedthe unsaturated conductivity of the model of Fayer andSimmons (1995) for a dune sand, as well as the bimodal
VG model for an Andisol, using the evaporation method
in combination with water retention data from close tosaturation to very dry conditions down to pressure headsof approximately −105 cm. They showed that the observed pressure heads agreed well with the model predictions when an appropriate hydraulic function was used
to describe the water retention data over a wide range of

The objective of this study was to determine parameters of the bimodal VG model for Andisols over a wide
range of pressure heads as determined using the evaporation method. In addition to the pressure heads we also usedwater retention data in the objective function. The saturated hydraulic conductivity, Ks , and the pore-connectivity
factor, ℓ, in the VG hydraulic conductivity function, alongwith several water retention parameters, were optimized.
Different sets of optimizations were conducted by restricting the pressure head range of the water retention and tensiometer measurements in the objective function. We were
especially interested in the role of water retention data atvery at low (negative) pressure heads in the objective function. We furthermore used the measured pressure headprofiles of the samples at the end of the evaporation experiment to confirm the shape of the estimated hydraulicfunctions beyond the tensiometer measurement range.

2. Material and methods
2.1 Evaporation experiment
For the evaporation experiments we used two different
Japanese Andisols collected from the surface horizons of
two sites. One set of soil samples was obtained from the
NARO Kyushu Okinawa Agricultural Research Center inKumamoto, and one from the NARO Institute of Vegetable
and Tea Science in Mie, Japan. Collected soil samples
were sieved using a 2 mm mesh. The disturbed Andisols were packed uniformly in 16-cm long, 3.8-cm diameter acrylic columns to bulk densities, ρb , of 0.48 and 0.75
g cm−3 for the Kumamoto and Mie soils, respectively. The
saturated hydraulic conductivity, Ks , based on the fallinghead method was estimated to be approximately 200 cm
d−1 for the Kumamoto Andisol, and 1000 cm d−1 for the
Mie Andisol. Fig. 1 shows water retention curves as measured using a hanging water column for the pressure head(h) range −200 < h < −5 cm, using a pressure plate for−1.2 × 104 < h < −250 cm, and a WP4 dew point potentiameter (Decagon Devices, Pullman, WA) based on relative humidity measurements equilibrated with the soil water pressure for the range −105 < h < −3 × 103 cm. Although WP4 measurements are generally used for pressure
heads below −104 cm, we applied the WP4 potentiameterto the higher pressure heads up to approximately −3 × 103
cm with considerable care as suggested by Maček et al.
(2013). Volumetric water contents for the WP4 measurements were determined from the gravimetric water contents and the bulk density of the soil.
After slowly saturating the soil samples from the bottom, the water supply was closed and evaporation was al-

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method

low to start using a fan to blow air away from the soil
surface in a 25 ◦ C constant temperature room. Pressure
heads were monitored using five tensiometers connected
to pressure transducers at 1, 2, 3, 5 and 10 cm depths.
Cumulative amounts of evaporation were calculated frommeasurements of the soil column weights using an electrical balance connected to a data logger. The evaporation experiment for the Mie Andisol continued until thepressure head became less than −600 cm at 1 cm depth.
The evaporation period for the Kumomoto Andisol was allowed to continue longer, until the pressure head at 10 cm
depth became less than −500 cm. The valves connected tothe tensiometers above 10 cm depth were closed at h values of about −500 cm to prevent water leakage from thetensiometer cups (Durner and Or, 2005; Schindler et al.,2010). Volumetric water content profiles at the end of experiments were determined gravimetrically by sectioning
the soil columns. Pressure head profiles near the surface
at the end of experiments for the Kumamoto Andisol were
measured using the WP4 potentiameter.
2.2 Bimodal van Genuchten model
The bimodal water retention (Durner, 1992, 1994) andunsaturated hydraulic conductivity (Priesack and Durner,
2006) functions based on the van Genuchten-Mualem
model (van Genuchten, 1980; Mualem, 1976) were used
to describe the hydraulic properties of the two Andisols:(1)

contents, while S2 corresponds to the second subregion at
the lower water contents. The hydraulic functions given by
Eqs. (1) to (3) are further referred to here as the bimodal
VG model. Note that the hydraulic functions contain a total of nine parameters. Seven of these parameters (θs , θr ,α1 , n1 , w2 , α2 , n2 ) relate primarily to the water retention
curve, and two additional parameters (Ks , ℓ) to the unsaturated hydraulic conductivity function.
2.3 Parameter optimization
Parameters of the bimodal VG model for our
two Andisols were optimized using version 4.08 of
the HYDRUS-1D software package (http://www.pcprogress.com/en/Default.aspx?hydrus-1d) of Šimůnek etal. (2008). The calculations assumed applicability of the
Richards equation to one-dimensional vertical water flow
under isothermal condition as follows:()∂ θ (h)
∂
∂h
=
K (h)
+ K (h)
∂t
∂z
∂z

(4)

t is time [T].

2(2)h (z, 0) = hi (z))ℓ(∑ wi Si−K (h)

i=1

 2
[
(
)mi ] 2


1/mi


w
1
−
1
−
S
α
∑
i
i
i


2

∑ wi αi





i=1

The subscript i in these equations represents the number ofsubregions (i = 1, 2 in our application), Si is effective saturation [−] of the i-th subregion, h is the soil water pressure
head [L], θ is the volumetric water content [L3 L−3 ], θs andθr are the saturated and residual water contents [L3 L−3 ],respectively, ni [-], αi [L−1 ], and mi (= 1 − 1/ni ) are shapeparameters subject to αi > 0 and ni > 1, w1 and w2 are theweighing factors subject to 0 < wi < 1 and w1 + w2 = 1,
K is the hydraulic conductivity [LT−1 ], Ks is the saturated

(5)

)
∂h
+ 1 = qevap (L,t)
∂z

(6)

)
∂h
+1 = 0
∂z

(7)

(3)

i=1





sociated with the first subregion of θ (h) at the higher water

experiment are given by (e.g., Šimůnek et al., 1998):

Si = [1 + (αi |h|)ni ]−mi×coefficient [-]. The S1 variable (subject to α1 > α2 ) is as-

The initial and boundary conditions for the evaporation

where Si is given by:

K (h) =Ks

hydraulic conductivity [LT−1 ], and ℓ is a pore-connectivity

where z is the vertical coordinate [L] positive upward, and2θ (h) − θrS== ∑ wi Si
θs − θr
i=1

(

5

(
q (0,t) = −K (h)where hi (z) is the initial pressure head [L] distribution in
the column, qevap (t) is the time-variable evaporation rate
[LT−1 ] at the soil surface, and L is the z coordinate of thesoil surface [L]. The initial condition, hi (z), was linearlyinterpolated using the initial tensiometer readings, while
qevap (t) was described with a polynomial function fitted
to the observed evaporation rate as shown in Fig. 2. A
no-flow boundary condition was imposed at the bottom
boundary. We note that vapor flow was not included in theflow model, thus assuming that the effects of vapor flow

土壌の物理性

6

第 125 号(2013)

Fig. 2 Evaporation rate from the surface of the Kumamoto (a) and Mie (b) Andisol soils samples. Sown are the observed evaporation
rate data (open circles) and the fitted polynomial equation (solid line).

on the parameter optimization process using the evapora-sample after the experiment. Estimates of the remaining

tion method can be neglected (Sakai and Toride, 2007b).

water retention parameters (θr , α1 , n1 , w2 , α2 , n2 ) in Eqs.

The objective function, Φ, to be minimized during the(1) to (3) were determined first from the observed water re-parameter estimation process consisted of pressure heads

tention data in Fig. 1 using version 6.02 of the RETC code

measured at several depths as well as the independently(http://www.pc-progress.com/en/Default.aspx?retc) of van

measured water retention data shown in Fig. 1 (Hopmans

Genuchten et al. (1991). In order to reduce the number of

et al., 2002; Šimůnek et al., 2008):

optimized parameters, the value of θr was assumed to benhzero for our two Andisols since we found that optimiza-Φ = vh ∑ [hobs (ti ) − hfit (ti )]
i=1
nθ

2

tion of θr did not improve the fit of the data in Fig. 1. We(8)+vθ ∑ [θobs (hi ) − θfit (hi )]2

i=1

similarly found that the restrictions mi = 1 − 1/ni in thebimodal VG model did not affect the goodness of fit. Fig.
1 shows excellent visual matches of the bimodal VG soil

where nh and nθ are the number of observed pressure heads(in the column) and independently measured water retention data, respectively, vh = 1/(nh σh2 ) and vθ = 1/(nθ σθ2 )are weights for each data type (in which σ is the varianceof the observed data), while the subscripts obs and fit indicate observed and model fitted values, respectively, at timeti .

water retention functions to the data.
Following Spohrer et al. (2006) and Sakai and Toride(2007b), the fitted retention parameter values were used
next as initial estimates in the overall optimization of theevaporation experiment. In our study we optimized the hydraulic conductivity parameters Ks and ℓ since these twoparameters are generally difficult to measure directly. The

Šimůnek et al. (1998) recommended including the final
water content of the evaporation experiment in the objective function in order to anchor the retention curve alongthe θ axis. Instead of including the final water content in
the objective function, we fixed the value of θs by applying
a mass balance to the column. The value of θs was estimated from the final water content of the column and thecumulative amount of evaporation during the evaporation
experiment as determined from the loss of weight of theinitial value for Ks was fixed at the value obtained with thefalling head method, while ℓ was initially assumed to be
0.5 as suggested by Mualem (1976). The values of θs and
Ks in an evaporation experiment may be slightly different
from observed static water retention measurements and thefalling head method because of entrapped air during saturation or the length of soil column used (Hopmans andDane, 1986; Dane and Hopmans, 2002; Sakaguchi et al.,2005). The value of Ks can easily change with very small

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method
Table 1 Values of the optimized hydraulic parameters of the bimodal VG model, and their standard errors (SE)and coefficients of variation (CV), obtained with the evaporation method for cases 1, through 7 of the Kumamoto
and Mie Andisols.
Case

θr
θs
(cm3 cm−3 )α1(cm−1 )n1(-)

Ks
(cm d−1 )

ℓ
(-)

w2
(-)

α2
(cm−1 )

n2
(-)

0
-

0.78
-

-

-

0.057
0.045
0.001
0.022
0.044
0.001
0.023
0.055
0.001
0.018
-

1.88
1.92
0.03
0.015
1.96
0.03
0.017
1.87
0.04
0.020
-

0.5
0.751
0.107
0.142
0.054
0.022
0.418
0.436
0.214
0.492
0.614
0.231
0.377
2.301
0.227
0.099
0.007
0.016
2.167
0.001
0.012
15.250

0.544
0.544
0.003
0.005
0.549
0.003
0.006
0.541
0.008
0.014
-

0.000110.000130.00001
0.076
0.000250.00006
0.232
-

1.44
-

-

200
703
11
0.016
312
26
0.083
343
32
0.093
312
38
0.122
971
18
0.019
208
42
0.202
28
6
0.214

1.39
0.02
0.011
1.55
0.31
0.203
-

-

0
-

0.629
-

0.129
0.101
0.002
0.020
0.111
0.002
0.018
0.131
0.002
0.015
-

2.41
4.36
0.06
0.013
3.57
0.07
0.020
2.42
0.04
0.016
-

1000
969
66
0.068
968
50
0.052
555
34
0.061
507
176
0.347
850
223
0.262
617
39
0.063
391
22
0.056

0.5
0.040
0.016
0.389
0.044
0.006
0.145
0.011
0.004
0.381
0.330
0.608
1.844
0.620
0.378
0.610
0.407
0.059
0.145
0.018
0.014
0.765

0.509
0.532
0.003
0.005
0.545
0.003
0.005
0.513
0.002
0.005
-

0.000620.001250.00009
0.072
0.00072
0.00004
0.056
-

1.28
1.25
0.01
0.005
1.74
0.07
0.043
-

Kumamoto
initial
1
SE
CV
2
SE
CV
3
SE
CV
4
SE
CV
5
SE
CV
6
SE
CV
7
SE
CV
Mie
initial
1
SE
CV
2
SE
CV
3
SE
CV
4
SE
CV
5
SE
CV
6
SE
CV
7
SE
CV

Note: Entries indicated by- were not included in the optimization but fixed at the initial values

7

土壌の物理性

8

第 125 号(2013)

Fig. 3 Observed and fitted pressure heads as a function of time for the Kumamoto (a) and Mie (b) Andisols. Case 1 (long dash line),case 2 (short dash line) and case 3 (solid line).

variations in θs . Since Andisols have a large porosity (es-found that using pressure heads at two different depths inpecially the Kumamoto soil in our study), direct measure-the objective function produced smaller standard errors forment of saturation can be quite vulnerable also to small

Ks and ℓ than when data from only one depth were used.

measurement errors. Since we fixed θs based on water

Since Ks and ℓ were our primary concern, we will show be-mass balances of the sample, Ks is probably best treated

low results when using pressure heads at 1 and 3 cm depths

as fitting parameter in the evaporation method. Initial val-in the objective function given by Eq. (8). Still, obtaining

ues of all of the parameters used in the optimizations for

additional tensiometer measurements at other depths may

the two Andisols are listed in Table 1.

well be useful for backup information in case some of ten-

We optimized the Ks and ℓ values for three cases involv-

siometers failed to work properly (Hopmans et al., 2002).

ing different combinations with the water retention param-

Fig. 3 shows fitted and observed pressure heads as aeters: (1) only two hydraulic conductivity parameters (i.e.,function of time at the 1 and 3 cm depths for the Kumamoto

Ks , ℓ) while keeping the water retention parameters fixedand Mie Andisols. Estimated water retention and unsat-

according to the independently measured water retention

urated hydraulic conductivity functions corresponding tofunctions, (2) adding the water retention parameters of thethe three optimization cases of our study are shown in Fig.

first subregion (i.e., Ks , ℓ, α1 , n1 , w2 ) in the optimization,

4. The optimized parameter values and their standard er-

and (3) further including the retention parameters of therors, as well as the coefficients of variations (defined as thesecond subregion (i.e., Ks , ℓ, α1 , n1 , w2 , α2 , n2 ). The val-

optimized values divided by the standard errors) for the bi-

ues of θs and θr (= 0) remained fixed in all optimizations.

modal VG model are listed in Table 1.

3. Results and discussion

For the Kumamoto Andisol, all three optimization scenarios (i.e., case 1 with (Ks , ℓ) optimized, case 2 with

3.1 Parameter estimation(Ks , ℓ, α1 , n1 , w2 ) optimized, and case 3 with (Ks , ℓ,In earlier work, Šimůnek et al. (1998) showed that a

α1 , n1 , w2 , α2 , n2 ) optimized) gave almost identical re-single set of tensiometer readings near the sample surface

sults. All three cases produced excellent agreement withwas sufficient to yield accurate estimates of the soil hy-

the observed pressure heads at the 1 and 3 cm depths (Fig.

draulic parameters using the evaporation method. Our ex-3a). Predicted pressure heads for cases 2 and 3 were al-

periments confirmed this in that similar results were ob-most identical. The three cases also produced close agree-

tained irrespective of including pressure head data fromment with the observed water retention data (Fig. 4a). In

points deeper in the columns. Sakai and Toride (2007a)fact, as shown in Fig. 4b, the estimated hydraulic con-

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method

9Fig. 4 Estimated water retention (a, c) and unsaturated hydraulic conductivity (b, d) curves for the Kumamoto (top) and Mie (bottom)
Andisols: case 1 (long dash lines), case 2 (short dash lines) and case 3 (solid lines).

ductivities were almost identical for all three Kumamoto

sure heads being greater than the observations. Cases 2

optimizations. In addition, all of the optimized parametersand 3 improved the agreement between observed and fitted

had small standard errors as reflected by the narrow con-pressure heads. As the number of optimized parameter in-fidence intervals shown in Table 1. These results suggest

creases, the bimodal VG model clearly has more flexibility

that when reliable water retention data are available andin fitting the pressure head data. Improved fitting of theused in the objective function, it is possible to predict K(h)pressure head data could be achieved also by sacrificing

by optimizing only the two conductivity parameters (Ks , ℓ).

the close fit of the water retention data. However, com-

This implies that as long as the water retention parameters

pared to the differences between the fitted and observed

can be properly determined using observed water retention

pressure heads in Fig. 3b, the fitted water retention curves

data, the number of parameters in Mualem-type soil hy-

in Fig. 4c were very similar for all three cases. The bi-

draulic functions, even for the bimodal VG model as given

modal water retention curve of Fig. 4c has a relatively

by Eq. (3), should not pose a problem in the optimization.

flat region between about −50 cm and −103 cm. Since

Results for the Mie Andisol were different. Case 1 forthe soil water capacity (d θ /dh) is small in that flat region,

this soil failed to fit the sudden pressure drop at h = −50
cm after about 3.5 days in Fig. 3b, with the fitted pres-the abrupt pressure drop below h = −50 cm in Fig. 3bshould not change the water content very much. This ex-

10

土壌の物理性第 125 号(2013)

Fig. 5 Estimated water retention (a, c) and unsaturated hydraulic conductivity (b, d) curves for the Kumamoto (top) and Mie (bottom)
Andisols. Results were obtained with (solid lines) and without (dashed lines) including water retention data at the low pressure heads
in the objective function.

plains why the small change in the water retention curve

tent and pressure head relationships when the soil waterobtained by additionally optimizing the first subregion re-capacity is very small. Pressure plates are generally used

tention parameters (α1 , n1 , w2 ) still produced a significant

to measure water retention in this region. Nonequilibrium

improvement in the pressure head simulations. By com-conditions, poor contact between the soil sample and theparison, the pressure range of the flat region of the Mieplate (Cresswell et al., 2008; Bittelli and Flury, 2009), andAndisol is wider and the soil water capacity smaller thanhysteresis in the sample preparation often affect the mea-of the Kumamoto Andisol (Fig. 4a versus and 4c). We

surements. For this reason it may be better to optimize atemphasize that the discrepancies in the pressure head wereleast the first subregion parameters (α1 , n1 , w2 ) also. But asapparent when optimizing only the conductivity parame-long as water retention data over a wide range of pressureters (case 1), and then only for the Mie Andisol, but not forheads are available and can be used in the objective func-

the Kumamoto soil.

tion, as Durner (1994) also suggested, case 3 optimization

We also note that measurements of water retention data

should yield the best overall results.

in the flat part of the retention function are often subject to

We mention some concern about the value of the pore-errors. It is not easy to precisely determine the water con-connectivity coefficient, ℓ, which is often fixed at 0.5 fol-

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method

11Šimůnek et al. (1998) and Hopmans et al. (2002). To test
the need for retention data in the dry range, we comparedhydraulic functions optimized with (case 3) and without(case 4) including the independent water retention data at
the low pressure heads in the objective function. The optimized conditions for case 4 were the same as for case 3,
including the initial estimates, except that case 4 considered only water retention data for h > −103 cm. The estimated parameter values for case 4 are also listed in Table 1.
For both Andisols, good agreement between the observed
and fitted pressure heads in the columns was obtained, very
similar to case 3 as shown in Fig. 3.
Fig. 5 shows the observed and fitted water retention
curves and the estimated unsaturated hydraulic conductivFig. 6 Observed and predicted pressure head profiles after 6.08
days for the evaporation experiment of the Kumamoto Andisol.
Results were obtained with (solid lines, case 3) and without(dashed lines, case 4) including water retention data at the lowpressure heads in the objective function. Circles represent WP4data and the square a tensiometer datum.

lowing Mualem (1976). As pointed out in several studies(e.g., Wösten and van Genuchten, 1988; Kaveh and van
Genuchten, 1992; Schaap and Leij, 2000; Spohrer et al.,2006), ℓ can be quite variable depending upon soil type.
Using the case 3 optimization, ℓ was found to be 0.436
for the Kumamoto Andisol and 0.011 for the Mie Andisol.
More data are clearly needed to obtain a better definitionof possible ℓ values for Andisols.
Our parameter optimizations were found to converge
quickly, in all cases less than 10 iterations (5, 8 and 8 iterations for cases 1 to 3, respectively, of the Kumamoto
Andisol and 8, 9 and 9 iterations for cases 1 to 3 of the
Mie Andisol). This included case 3 which had the largest
number of optimized parameters. Our results suggests thatusing water retention parameters fitted to independently
observed retention data as initial estimates will facilitate
convergence.
3.2 Water retention data at low pressure heads
In addition to the pressure head measurements obtained
during the evaporation experiments for h > −600 cm, we
included in the objective function thus far also all of the independently measured water retention data between −105
and −5 cm (Fig. 1). This is to obtain reliable estimates ofK(h) in the low pressure head range as recommended by

ity functions for cases 3 and 4. The observed water retention data in the dry range were severely underestimatedusing case 4 relative to case 3. Since no information wasgiven for h < −103 cm for the case 4 optimizations, it is
not surprising that this case converged to different valuesof the second subregion parameter (α2 , n2 ). The discrepancies shown in Fig. 5 for case 4 are a reason why optimized hydraulic functions are generally assumed to bereliable only within the range of the tensiometer measurements (Šimůnek et al., 1998; Hopmans et al., 2002). Table
1 shows that most of the standard errors of the estimatedparameters for case 3 were smaller than those for case 4,
which indicates that including water retention data at thelower pressure heads in the objective function will reduceparameter uncertainty and lead to more reliable parametervalues.
Although observed pressure head data over a 5-day period were used for parameter optimization of the Kumamoto Andisol (Fig. 3a), the evaporation experiments
continued for 6.08 days until the pressure heads at 10 cm
depth reached approximately −500 cm. To confirm the accuracy of the hydraulic properties of cases 3 and 4 for the
Kumamoto Andisol, pressure heads were simulated for alonger time period using case 3 and 4 properties. Fig. 6
shows observed WP4 data near the surface with one tensiometer datum at the 10-cm depth, and predicted final
pressure head profiles after 6.08 days. The simulation using case 3 parameters agreed well with the observed pressure profile near the surface, whereas the prediction usingcase 4 parameters overestimated the pressure heads. These

12

土壌の物理性第 125 号(2013)

Fig. 7 Impacts of the pressure head range in the objective function on the estimated unsaturated hydraulic conductivity for the
Kumamoto (a) and Mie (b) Andisols: case 1 used all data, case 5 only h > −200 cm, case 6 only h > −100 cm, and case 7 onlyh > −50 cm.

results again indicate that including water retention data at

timizations: only data for which h > −200 cm (case 5),pressure heads below the tensiometer range (about −500

h > −100 cm (case 6), and h > −50 cm (case 7). The es-cm) will improve the estimates of the soil hydraulic pa-

timated parameter values for these three cases (5 to 7) are

rameters.

given in Table 1.

Since water retention at low pressure heads now can be

Fig. 7 shows the estimated K(h) plots for the Kumamoto

measured relatively easily and accurately using dew point

and Mie Andisols obtained with the different pressure headpotentiameters (e.g., WP4), we strongly recommended tomeasurement ranges. The conductivity function for caseinclude in the optimization also water retention data in

1 is the same as shown in Fig. 4. Except for the muchthe measurement range of the dew point potentiameterlower K(h) curve for case 7 of the Kumamoto Andisol,(−105 < h < −3 × 103 cm) to extend the applicable range

the plots neglecting the lower pressure data did not affectof the model predictions. On the other hand, estimates of

the K(h) results very much; they produced almost identi-

the hydraulic conductivity at the very low pressure heads

cal K(h) curves as for case 1. For the Mie Andisol, even(e.g., h < −104 ) may require further investigations because

case 7 did give a good estimate of K(h), likely becauseof the possible effects of film flow (Peters and Durner,h > −50 cm covered the entire region of the first slope of2008), diffusion of water vapor (Saito et al., 2006; Peters

the water retention curve (Fig. 1). We conclude that inde-and Durner, 2010) and the effect of shrinkage (Dorel et al.,pendently measured water retention data covering a wide
range of pressure heads provide very useful information to2000).
3.3 Pressure head measurement range

the parameter estimation process, leading to a much more

We next demonstrate the effect of using different ranges

robust optimization of the remaining parameters (Ks , ℓ inof the pressure head data measured in the columns during

this case). Note that the standard errors and coefficients ofthe evaporation experiments in the inverse analysis, whilevariation have a tendency to increase as the pressure range

fixing in the objective function the bimodal water retention

becomes narrower (Table 1). This means that one shouldparameters (θr , θs , α1 , n1 , w2 , α2 , n2 ) as determined fromstill use all available tensiometer data in the objective func-

the independently measured retention data. These scenar-tion.

ios are the same as case 1 in that only the conductivity
parameters (Ks , ℓ) are optimized, except that we used different ranges of the measured pressure heads in the op-

4. Conclusions
Parameters of the bimodal VG model for two aggregated Andisols were inversely determined using the evap-

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method

13

oration method. Independently measured water retention

Neglecting the lower pressure data did not affect the K(h)data from near saturation to very low pressure heads down

estimation process very much, leading to almost similar

to −105 cm were included in the objective function in addi-

K(h) functions. This confirms that collecting water reten-tion to soil pressure head data at two depths measured dur-tion data over a wide range of pressure heads will give very

ing the evaporation experiments. The saturated hydraulicuseful prior information to the parameter estimation pro-parameter, Ks , and the pore-connectivity factor, ℓ, alongcess.

with several sets of water retention parameters were opti-

References

mized. Since direct θs measurements are often quite variable for high-porosity Andisols, the saturated water content was determined from the final water content and themeasured cumulative amount of evaporation. When the

Bittelli, M. and Flury, M. (2009): Errors in water retention curves
determined with pressure plates and their effect on soil hydraulic functions. Soil Sci. Soc. Am. J., 73: 1453–1460, doi:
10.2136/sssaj2008.0082.

value of θs is fixed from these water mass balance consid-

Chamindu Deepagoda, T.K.K., Moldrup, P., Jensen, M.P., Jones,

erations, Ks should be estimated in the evaporation method.

S.B., de Jonge, L.W., Schjønning, P., Scow, K., Hopmans,

Since the initial estimates of the water retention parameters were determined from the independently measured
water retention data, parameter estimation succeeded toconverge quickly regardless of the number of optimized
parameters. When water retention data from near saturation to very low pressure heads are available and used in

J.W., Rolston, D.E., Kawamoto, K. and Komatsu, T. (2012):
Diffusion aspects of designing porous growth media for earthand space. Soil Sci. Soc. Am. J., 76: 1564–1578, doi: 10.2136/
sssaj2011.0438.
Cresswell, H.P., Green, T.W. and McKenzie. N.J. (2008): The adequacy of pressure plate apparatus for determining soil water
retention. Soil Sci. Soc. Am. J., 72: 41–49, doi: 10.2136/sssaj2006.0182.

the objective function, it is possible to predict K(h) by op-

Coppola, A. (2000): Unimodal and bimodal descriptions of hy-

timizing only two conductivity parameters (Ks , ℓ). Since

draulic properties for aggregated soils. Soil Sci. Soc. Am. J.,

the flat region of the bimodal water retention curve is difficult to measure precisely, we recommend optimizing allof the bimodal VG parameters except θs and θr to yield thebest overall fit.
In order to demonstrate the role of water retention data
at low pressure heads in the objective function, we com-

64: 1252–1262, doi: 10.2136/sssaj2000.6441252x.
Dane, J.H. and Hopmans, J.W. (2002): Water retention storage. In
Dane, J.H. and Topp, G.C. (ed.) Methods of soil analysis. Part
4. pp. 675–680. SSSA Book Ser. 5. SSSA, Madison, WI.
Dorel, M., Roger-Estrade, J., Manichon, H. and Delvaux, B.
(2000): Porosity and soil water properties of Caribbean volcanic ash soils. Soil Use Manage., 16: 133–140, doi: 10.1111/
j.1475-2743.2000.tb00188.x.

pared hydraulic functions optimized with and without wa-

Durner, W. (1992): Predicting the unsaturated hydraulic con-ter retention data at low pressure heads. Although almost

ductivity using multi-porosity water retention curves. In vansimilar matches to observed pressure heads were obtained

Genuchten, M.Th., Leij, F.J. and Lund, L.J. (ed.) Indirect

with or without the low pressure data in the dry range, thesecond subregion parameter (α2 , n2 ) converged to different

methods for estimating the hydraulic properties of unsaturated
soils: Proc. Intl. Workshop, Riverside, pp. 185–202. CA. 11–
13 Oct. 1989. U.S. Salinity Lab., Riverside, CA.

values if the low pressure data were omitted from the ob-

Durner, W. (1994): Hydraulic conductivity estimation for soils

jective function. Predictions based on the hydraulic func-

with heterogeneous pore structure. Water Resour. Res., 30:

tions optimized with the low pressure data agreed well with

211–223, doi: 10.1029/93WR02676.
Durner, W. and Or, D. (2005): Soil water potential measurement.

observed pressure heads near the surface after a longer pe-

In Anderson M.G. and McDonnell, J.J. (ed.) Encyclopedia of

riod of evaporation than used for the tensiometer measure-

Hydrological Sciences, pp.1089–1102. John Wiley & Sons,

ments. The results indicate that including water retention
data at low pressure heads can extend the applicable rangeof the model predications, at least down to approximately−104 cm.

Ltd., doi: 10.1002/0470848944.hsa077a.
Durner, W. and Iden, S.C. (2011): Extended multistep outflow
method for the accurate determination of soil hydraulic properties near water saturation. Water Resour. Res., 47: W08526,
doi: 10.1029/2011WR010632.

The benefit of using independently measured water re-

Fayer, M.H. and Simmons, C.S. (1995): Modified soil water re-tention data was further studied by including different pressure head measurement ranges in the objective function.

tention functions for all matric suctions. Water Resour. Res.,
31: 1233–1238, doi: 10.1029/95WR00173.

14

土壌の物理性Gerke, H.H. and van Genuchten, M.Th. (1993): A dual-porosity

第 125 号(2013)

model for simulating the preferential movement of water andPeters, A. and Durner, W. (2010): Reply to comment by N.
Shokri and D. Or on ”A simple model for describing hy-

solutes in structured porous media. Water Resour. Res., 29:

draulic conductivity in unsaturated porous media accounting

305–319, doi: 10.1029/92WR02339.

for film and capillary flow”. Water Resour. Res., 46: W06802,

Hamamoto, S., Perera, M.S.A. Resurreccion, A.C., Kawamoto,doi:10.1029/2010WR009181.

K., Hasegawa, S., Komatsu, T. and Moldrup. P. (2009): The

Peters, R.R. and Klavetter, E.A. (1988): A continuum model forsolute diffusion coefficient in variably compacted, unsatu-

water movement in an unsaturated fractured rock mass. Water

rated volcanic ash soils. Vadose Zone J., 8: 942–952, doi:
10.2136/vzj2008.0184.
Hopmans, J.W. and Dane, J.H. (1986): Temperature dependenceof soil water retention curves. Soil Sci. Soc. Am. J., 50: 562–

Resour. Res., 24: 416–430, doi: 10.1029/WR024i003p00416.
Priesack, E. and Durner, W. (2006): Closed-form expression
for the multi-modal unsaturated conductivity function. Vadose
Zone J., 5: 121–124, doi:10.2136/vzj2005.0066.

567, doi: 10.2136/sssaj1986.03615995005000030004x.

Saito, H., Šimůnek, J. and Mohanty, B. (2006): Numerical anal-

Hopmans, J.W., Šimůnek, J., Romano, N. and Durner, W. (2002):

yses of coupled water, vapor and heat transport in the vadose

Simultaneous determination of water transmission and reten-zone. Vadose Zone J., 5: 784–800, doi:10.2136/vzj2006.0007.

tion properties. Inverse methods. In Dane, J.H. and Topp G.C.

Sakaguchi, A., Nishimura, T. and Kato, M. (2005): The effect of(ed.) Methods of Soil Analysis, Part 4, pp. 963–1008. SSSA

entrapped air on the quasi-saturated soil hydraulic conductiv-

Book Ser. 5. SSSA, Madison WI.

ity and comparison with the unsaturated hydraulic conductiv-

Iden, S.C. and Durner, W. (2007): Free-form estimation of theity. Vadose Zone J., 4: 139–144, doi:10.2136/vzj2005.0139.

unsaturated soil hydraulic properties by inverse modeling us-

Sakai, M and Toride, N. (2007a): Optimum conditions for pre-ing global optimization. Water Resour. Res., 43: W07451, doi:

dicting unsaturated hydraulic properties using the evaporation

10.1029/2006WR005845.

method. J. Jpn. Soc. Soil Phys., 106: 33–46.

Kaveh, F. and van Genuchten, M.Th. (1992): A further look at

Sakai, M and Toride, N. (2007b): Soil water hydraulic functions

a new unsaturated hydraulic conductivity equation. Iranian J.

for a sandy soil and an aggregated soil. J. Jpn. Soc. Soil Phys.,

Agric. Sciences, 23: 24–32.

107: 63–77.

Maček, M., Smolar, J. and Petkovšek, A. (2013): Extension

Schaap, M.G. and Leij, F.J. (2000):

Improved predictionof measurement range of dew-point potentiameter and evap-of unsaturated hydraulic conductivity with the Mualem-van

oration method. Proc. 18th Intl. Conf. Soil Mechanics andGenuchten model. Soil Sci. Soc. Am. J., 64: 843–851, doi:

Geotechnical Engineering, pp. 1137–1142, Paris.

10.2136/sssaj2000.643843x.

Mallants, D., Jacques, D., Tseng, P.H., van Genuchten, M.Th.

Schelle, H., Iden, S.C., Peters, A. and Durner, W. (2010): Analysis

and Feyen, J. (1997): Comparison of three hydraulic prop-of the agreement of soil hydraulic properties obtained fromerty measurement methods. J. Hydrol., 199: 295–318, doi: 10.

multistep-outflow and evaporation methods. Vadose Zone J.,

1016/S0022-1694(96)03331-8.

9: 1080–1091, doi: 10.2136/vzj2010.0050.

Miyamoto, T., Annaka, T. and Chikusi, J. (2003): Soil aggregate

Schindler, U., Durner, W., von Unold, G., Mueller, L. andstructure effects on dielectric permittivity of an Andisol mea-

Wieland, R. (2010): The evaporation method — Extending

sured by time domain reflectometry. Vadose Zone J., 2: 90–97,

the measurement range of soil hydraulic properties using thedoi: 10.2136/vzj2003.9000.

air-entry pressure of the ceramic cup. J. Plant Nutr. Soil Sci.,

Mohanty, B.P., Bowman, R.S., Hendrickx, J.M.H. and van

173: 563–572, doi: 10.1002.jpln200900210.

Genuchten, M.Th. (1997): New piecewise-continuous hy-

Shoji, S., Nanzyo, M. and Dahlgren, R.A. (1993): Volcanic ash

draulic functions for modeling preferential flow in an intermit-

soils: Genesis, properties, and utilization, pp. 1–288. Elsevier,tent flood-irrigated field. Water Resour. Res., 33: 2049–2063,
doi: 10.1029/97WR01701.

Amsterdam, Netherlands.
Šimůnek, J., Wendroth, O. and van Genuchten, M.Th. (1998): Pa-

Mualem, Y. (1976): A new model for predicting the hydraulicrameter estimation analysis of the evaporation method for de-conductivity of unsaturated porous media. Water Resour. Res.,

termining soil hydraulic properties. Soil Sci. Soc. Am. J., 62:

12: 513–522, doi: 10.1029/WR012i003p00513.

894–905, doi:10.2136/sssaj1998.03615995006200040007x.

Nanzyo, M. (2002): Unique properties of volcanic ash soils. Glob.
Environ. Res., 6: 99–112.

Šimůnek, J., Šejna, M., Saito, H., Sakai, M. and van Genuchten,
M.Th. (2008): The Hydrus-1D software package for simu-

Othmer, H., Diekkrüger, B. and Kutilek, M. (1991): Bimodal

lating the movement of water, heat, and multiple solutes inporosity and unsaturated hydraulic conductivity. Soil Sci.,

variably saturated media. Version 4.0. HYDRUS Softw. Ser.

152: 139–150, doi: 10.1097/00010694-199109000-00001.

3. Dep. of Environ. Sci., Univ. of California, Riverside.

Peters, A. and Durner, W. (2008): A simple model for describ-

Spohrer, K., Herrmann, L., Ingwersen, J. and Stahr, K. (2006):

ing hydraulic conductivity in unsaturated porous media ac-

Applicability of uni- and bimodal retention functions for water

counting for film and capillary flow. Water Resour. Res., 44:

flow method in a tropical acrisol. Vadose Zone J., 5: 48–58,

W11417, doi: 10.1029/2008WR007136.

doi:10.2136/vzj2005.0047.

論文：Estimating the unsaturated hydraulic conductivity of Andisols using the evaporation method

15

van Genuchten, M.Th. (1980): A closed-form equation for pre-Wösten, J.H.M. and van Genuchten, M.Th. (1988): Using tex-

dicting the hydraulic conductivity of unsaturated soils. Soil

ture and other soil properties to predict the unsaturated soil

Sci. Soc. Am. J., 44: 892-898, doi: 10.2136/sssaj1980.

hydraulic properties. Soil Sci. Soc. Am. J., 52: 1762–1770.

03615995004400050002x.

Zurmühl, T. and Durner, W. (1998): Determination of parame-

van Genuchten, M.Th., Leij, F.J. and Yates, S.R. (1991): The

ters for bimodal hydraulic functions by inverse method. Soil

RETC Code for Quantifying the Hydraulic Functions ofSci. Soc. Am. J., 62: 874–880, doi: 10.2136/sssaj1998.

Unsaturated Soils, Version 1.0. EPA Report 600/2-91/065,

03615995006200040004x.

R.S. Kerr Environmental Research Laboratory, U.S. Environmental Protection Agency, Ada, OK, 85 p. (http://www.pcprogress.com/en/Default.aspx?retc).

要

旨団粒構造を持つ 2 種類の黒ボク土に対して，蒸発法を用いて bimodal van Genuchten（VG）モデルのバ
ラメータを逆解析により決定した．モデルを水分保持曲線の実測値に適合した初期値を用いると，推定パラメータの数にかかわらず収束は速い．また飽和付近から −105 cm 程度の低圧力水頭までの実測データに適合した bimodal VG モデルの水分保持曲線のパラメータ値を固定すると，透水係数の 2 個のパラメータ（Ks , ℓ）のみの最適化により不飽和透水係数 K(h) の推定が可能であった．しかし階段状の
水分保持曲線の平坦な中間圧力領域における正確な測定は難しいため，結果全体に対する最適な結果を得るために bimodal VG モデルのパラメータはすべて適合することが望ましい．乾燥領域の低圧力領域
までの水分保持曲線のデータを目的関数に含めると，モデルの適用範囲は，−104 cm 程度まで広がった．
キーワード：不飽和透水係数，水分保持曲線，団粒土，黒ボク土，蒸発法