J. Fo~ Res. 4 : 1 0 7 - 1 1 4 (1999)

Spatial Variability of Soil Hydraulic Properties in a Forested Hillslope Hendrayanto,1 Ken'ichirou Kosugi, Taro Uchida, Sakiko Matsuda, and Takahisa Mizuyama Division of Forestry and Bio-material Sciences, Graduate School of Agriculture, Kyoto University, Sakyo-ku, Kyoto 606-8502, Japan. Spatial variability of soil hydraulic properties was measured in a forested hillslope and analyzed by applying combined water-retention-hydraulic-conductivity models (the LN and VG models) and a power function model for soil hydraulic conductivity (the Leibenzon model). Results showed that the pore-tortuosity parameter I in the LN and VG models should be treated as a fitted parameter for accurate descriptions of the unsaturated conductivity. The simultaneous optimization of both retention and conductivity curves is preferable to obtain appropriate descriptions of hydraulic properties when both retention and conductivity data are available for the parameter estimation. The Leibenzon model produced slightly poorer estimates than the LN and VG models. The exponent parameter in the Leibenzon model exhibited a spatial variation with the standard deviation of 2.38. The spatial variability of hydraulic properties was analyzed based on the spatial variation in parameters of the LN model obtained by the simultaneous optimization procedure. The parameter ~m, which has a positive correlation with median pore-radius, was generally small at the crest and upper slope locations and large at mid-slope to footslope locations. Except for the crest, surface soils had larger ~m values than the subsurface soils, suggesting a well-developed crumb structure in surface horizons of forest soils. For most soils, cr was greater than 1, indicating a relatively large width of pore-size distribution. Saturated hydraulic conductivity (K0 was generally small at crest and upper slope locations, and large at mid-slope to footslope locations. The larger K, values were attributable to larger ~, values. Key words: forest soil, hydraulic conductivity, parameter estimation, soil water retention, spatial variability

Soils are spatially variable natural bodies. Variations of parent material and vegetation across the landscape from which soils are derived affect the variability of soils even at relatively short distances (Kutflek and Nielsen, 1994). Because of this spatial variability, the transport and retention properties of a field-scale soil unit also exhibit spatial variation that influence mass transport through the subsurface zone (Beckett and Webster, 1971; Freeze, 1975; Hoeksema and Kitanidis, 1985; Warrick and Nielsen, 1980; Jury et al., 1987). The knowledge of spatial variability of soil hydraulic properties is important for accurated solving the subsurface flux of water. Hence, many studies on spatial variability of soil hydraulic properties have been conducted (Mohanty et al., 1994; Clausnitzer et al., 1992; Peck et al., 1977; Hills et al., 1989; Nielsen et al., 1973 etc.). However, most of these studies were conducted in agricultural soils; spatial variability of hydraulic properties of forest soils has rarely been measured. Forest soils contain abundant macropores (Noguchi et al., 1997; Uchida et at., 1996; Kitahara, 1994; Beven and Germann, 1982), particularly in surface layer as a result of faunal activity and high root density (Bonell, 1993). Buttle and House (1997) noted that the distribution of soil macropores is an important factor governing the variability of soil hydraulic properties. The purpose of this study is to measure and analyze spatial variation in soil hydraulic properties in a forested hillslope. Many functional models for soil water retention and hydraulic conductivity have been proposed and used to analyze water flow in forest soils. Sammori and Tsuboyama (1990) adopted the combined water retention-hydraulic conductivity model proposed by van Genuchten (1980) to analyze slope stability taking infiltration into consideration. Shinomiya et al. LCorresponding author.

(1998) recently utilized the van Genuchten model to analyze the vertical variation in unsaturated conductivity of forest soils. Kosugi (1997c) evaluated effects of pore radius distribution on vertical water movement in soil profile by using the combined water-retention-hydraulic-conductivity model based on a lognormal pore-size distribution. The simple power function model for unsaturated conductivity proposed by Leibenzon was combined with some retention models to analyze water flow at forested hillslope (e.g., Kubota et al., 1987; Tani, 1985). By applying these models to observed soil hydraulic data sets, soil hydraulic properties can be characterized with only a small number of parameters. Moreover, spatial variations of soil hydraulic properties are expressed by the spatial distributions of these parameters. However, the accuracy of these models, especially the accuracy of hydraulic conductivity models, was not examined for undisturbed forest soils. For example, the exponent parameter in the van Genuchten conductivity model is usually fixed at 0.5 which was suggested by Mualem (1976) as an average for 45 disturbed soils. Values for the exponent parameter of the Leibenzon's conductivity model were rarely determined by using measured conductivity data for forest soils. The objectives of this study are (1) to measure soil hydraulic properties in a forested hillslope, (2) to examined the performance of several models and to evaluate the methods of parameter optimization using measured hydraulic data for forest soil, and (3) to analyze the spatial variability of hydraulic properties using parameters of the applied retention and conductivity models.

Materials and Methods 1 Measurement of soil hydraulic properties Soil samples were taken in a forested hillslope in a head water catchment of the Sumiyoshi River Basin (see Fig. la).

108

J. For. Res. 4 (2) 1999:

(a)

A,B, ...FI:S L Crest ~L (b) ~" ";" 870" A O-horizon-~B l...Ir;;m;;..

Uoeer

.L

Slope

Mid-Slope~L Lower J Foot

"~"

"~" Slope ~ Slope

C

_~.

_l.illlm,.,!,m]p.m.~ n "]~ll,lr.~.nll,t/, , , , , ,': ;'i~sm.~ ~41111~J~l~;k~.. IF,'.'l'~%~'~'~'ak%'~'~'~

w~n

n

51

~

%'~'d~%'¢~'~'lPl~|i;~..

A-h&~-on

,

''.':",';';'$",r'['X",!!!|st

I1-

B

o

2o

-

/ ~- . . . . . . .

40

~u

,

r

.~S ..... I /

|

• "

80

I

I~

I k.

I I

~

'" , , , , , e ¢ , . ' , ' r , . ' ~ l l i i l i , l

"".','.',' '7,%~s =;

10E i

tl

"

,~

1'°,

: Softsamplelocations |

0

2'0

40

60

I

80

100 1;0 14~0 Distance from the crest (m)

Fig. 1 (a) Contour map and (b) cross sectional view of the slope showing soil depth and No-valuesat the soil samplinglocations. The Basin is located on the western slope of Mount Nishiotafuku in Rokko Mountain range, Hyogo Prefecture, Japan. It is situated at 34 ° 46'00" N, 135 ° 15'48"E. The Basin consists of weathered granite, with very steep slopes of about 35 ° (see Fig. lb). The O-horizon is about 4-8 cm thick. The Ohorizon is thicker near the crest of the slope compared to the lower slope regions. The A-horizon is about 15-20 cm thick in the mid-slope and footslope regions and 5 cm thick near the crest. The B-horizon is more than 50 cm thick. The area is covered with a closed forest of predominately Quercus serrata, Styrax japonica, and Clethra barvinervis. Soil samples were taken from six sites distributed from the crest to the footslope (see Fig. lb) by using PVC cylinders of 19.5 cm in inner diameter, 50-60 cm long, and 10 mm thickness. The sampling technique was that proposed by Ohte et aL (1989). The PVC cylinder equipped with a sharp-edged ring was inserted vertically into the soil without disturbing the soil structure. During the process of insertion, root systems were carefully cut off from the soil layer around the sampled col-

umn. The sampling locations are referred to as A through F from the crest to the footslope, herein after. For each location, hydraulic properties of five soil layers with the same thickness of 10 cm (the 1st, 2nd, " ' , 5th layers from the soil surface) were determined. That is, 30 sets of soil water retention and hydraulic conductivity data were analyzed in this study. Figure lb shows the penetrability of soil at six sampling locations indicated by Nc-values. The Nc-value is the number of blows required for a 10 cm soil depth penetration measured using the cone-penetrometer with a cone diameter of 2.5 cm, a weight of 5.0 kg and falling distance of 50.0 cm (Osaka, 1996). The method to measure the unsaturated soil hydraulic properties was the same as the method proposed by Hendrayanto et al. (1998). Tensiometers (DIK-3150) connected to pressure transducers (Copal-PA800) and TDR (MP-917) were used to measure changes in profiles of matric potential head V and volumetric soil water content 0 during water redistribution. Unsaturated hydraulic conductivity K was computed based on the instantaneous profile method (Watson,

Hendrayanto et al.

109

1966). The steady-state method proposed by Ohte et al. (1989) was used to measure the saturated hydraulic conductivity Kx. 2 Models for soil hydraulic properties The observed hydraulic data sets for the 30 soils were analyzed by using combined water-retention-hydraulic-conductivity model proposed by Kosugi (1996) (the LN model) and van Genuchten (1980) (the VG model), and the conductivity model proposed by Leibenzon (1947). In the LN model, a lognormal pore-size distribution was combined with Mualem's model (1976) to derive a conductivity function. Because it employs physically based parameters, the LN model is effective for analyzing the hydraulic properties and soil water movement in connection with the soil pore-size distribution (Kosugi, 1997a, b, c). Based on the LN model, the water retention and hydraulic conductivity are expressed as: Se = ( 0 - Or)/(Os -- Or) = a( l n( ~lvm)l cr) (1)

K ( v ) = K ' [ Q ( l l n ( V ] ] ] t [ Q ( l l n ( IF )+6]] z L ~o" !kl//'m./.,/J l" IkO" !kl/'/mJ JJ

(2)

where, Se is the effective saturation, 0, and Or are saturated water content and residual water content, respectively, K, is saturated hydraulic conductivity, and Q is the complementary normal distribution function. The parameter ~m is the matric potential head at Se = 0.5, which is related to the median pore radius by the capillary pressure function; cr (dimensionless) is the standard deviation of log-transformed pore radius and represents a width of pore radius distribution (Kosugi, 1996). The parameter l was introduced by Burdine (1953) to evaluate soil-pore-tortuosity effects on hydraulic conductivity by using a power function of the effective saturation. The value of I was suggested to be 0.5 by Muallem (1976) as an average of 45 soils. The LN model has six parameters: &, Or, I//m, O', K,, and 1. The VG model has been extensively used for numerical modeling of soil water flow. The VG model expresses water retention as follows:

&-- 1/{ 1 +(a[ ~l)"}"

(3)

where ~zand n (n > 1) represent fitting parameters and m is related to n by m = 1 - 1/n (0 < m < 1). Combining Eq. (3) with Mualem's model (1976) yields,

~ Ki-logKi) oBJ = j=,(°' -0,) 2+ w,=~x (log

{1 +(air, I).}-"] = K(V) = K~

{1

o}n'

els to obtain the best fitting to the observed data sets. For that purpose, the non-linear least squares optimization method was used. The best fit was achieved by minimizing the residual sum of squares (RSS) between the computed and measured 0 and K values. (1) Method 1 (separate optimization) For the separate optimization procedure the 0(V) and K(V) were fitted to observed curves separately. In order to fit the water retention model to the observed 0(gt) curves, ~m and o"of the LN model (a and n of the VG model) were optimized by minimizing RSS. As a fitted parameter, & was sometimes optimized to be greater than 1. To avoid this unrealistic scenario, & was fixed at a maximum measured 0 value for each soils. The value of &. was optimized for 3 soils for which near saturation measured data were not available (the optimized & values were reasonable for these soils). The optimized values of Or were extremely small (sometimes smaller than zero) for four soils. For these soils, Or was fixed at the minimum measured 0 values. Functional K(V) curves were computed by substituting the optimized ~m and o', and ot and n into Eqs. (2) and (4), respectively. Saturated hydraulic conductivity (K~.) was fixed at measured values except for the first layer of the sampling location B; here measured K, was missing and an unsaturated hydraulic conductivity value at the highest measured V value was used as the matching point in Eqs. (2) and (4) (Luckner et al., 1989). To evaluate effects of the parameter I on the K estimation, three methods (methods la, lb, and lc) were examined. Method la used the constant I value of 0.5 as suggested by Mualem (1976). In method lb, l was optimized to minimize RSS comparing observed versus computed log K values for the entire data sets. In method lc, l was optimized for each soil to minimize RSS comparing observed versus computed log K values for each soil. (2) Method 2 (simultaneous optimization) Simultaneous optimization was used to fit water retention and hydraulic conductivity models simultaneously to observed 0(~) and K(lff) curves. Parameters, ~m, and ~ of the LN model (~x and n of the VG model) and l were optimized for each soil to minimize the objective function (OBJ) defined as follows (Yates et al., 1992):

(4)

The VG model contains six parameters: 6, Or, ~z, n, K~, and 1. Leibenzon proposed a simple power function model for hydraulic conductivity: K = K,Se # (5) where the exponent ]3 is an empirical parameter which is related to the pore-size distribution (Kutflek and Nielsen, 1994). 3 Model application to observed hydraulic data sets 1) Application of the LN and VG models Several methods were used to apply the LN and VG rood-

l

(6)

J

7 JE='eJ W

2

(7)

1 ~

T

i=~1log Ki

where 0j and t}, are the observed and calculated water content, respectively, Ki and /~i are the observed and calculated hydraulic conductivity, respectively, and J and I are the number of retention and conductivity data points, respectively. The weighting factor W corrects the difference in number of data points and difference in units for 0 and K. The simultaneous optimization used the same 0,. and Or values used in method 1.

110

J. For. Res. 4 (2) 1999: 0.8 .~ •

iI

o

i

;I

I I

"

0

i

I f

I

.2

0.6 4 o

I

i .%

at

"~

.L . . . .

~'%-

.4

0.6

•

1

I

L

1

I

I

I

I

I I

I I

I I

I I

I

I

I

I

I

t

I

" ~'~q,_

......

0.4

r

D

u

0.2 "[[ 0.8

-

'

"~ 0.6

-TI A

•

-6

l

l

l

I

1

1I

0

I I I

I I I

I I ~llo I

i

I I I

I I

I I

I ..... i I

I I

e~

I

I

1

1

1

1

I

I

I

t

_~ 0.6 "~ 0.4

~

,%

I I

me

. . . .

H. . . . . I I

4

.E -4

I I

~ 0.2 ~ 0.8

. . . . .

0.6

~

I I . . . .

I

. . . . .

4

I I

.......

0.4 0.2 0.8 0.6

.....

0

H . . . . .

~

';

-2

I

I

f

I I

I I

I

I

I

I I I

', I

I

I

I i

on o

i

.......

" ~i•

~- - - ~

~. . . . . . .

lll¥lvk

I':-.

- - -~ * - : - % ' 4 , , ~ ; ~

•'•%

~.... _" ?

T,,- ......

*

dooo

I

1

;

I

i ...... I

I

J ....... I

•

e"

¢

1 1

I

•

-4

I ;

l~

. . . .

L ...... I

a

I

o

1

I

I I

I l

I I

I

I

I

0

4 .....

4 ......

'~"

1 l

I

I . . . . . .

•

J . . . . . . .

I

I

I I

°

J

L .......

i ....

..... ~,~,,,~

I

0-10 cm

.

I

:II

D

'1

I

.

1- . . . . . . •

J . . . . . . .

,**]

. . . . . .

1

-40

.

.L . . . . . .

I

1

l---

.

I

0

-20

~dr "L . . . . . i

.

.~," I o °

I . ...... •

1

F _'~_~.~...L

0.2

_~dF~'. • ee0e •o~

a~'*~k.._

~e

I

-~. . . . .

,°

_ _ _ _.~.,.~,,,~,,.~_.~

. . . . . .

I

~- 0-~' ....

l

[

E

v~lul

0.4

I

I

i

I

i

I

,~. . . . . .

I

.~-

l

I

o

I

i

i

-* 0

1

,

oo

o

d

~.......

"1". . . . . .

C

I

I

od!¢1oo

I I

I I . . . . .

°

?...* *°°°°~ ° " " . .

I I I

° ,e

M

ee

* * * *

.

,"...--F__--~=ll.~l.o].,.o I i

~

I n

s I

. . . . . . .

[ I

• I I _ , l % ~ _ . a l , ~ ~ lJ

1Pk IJ'

--I

I

I I

. . . . . .

.

i~ill_•

.......

I I

I

0.8

.

I

n

I 1 I

......

.~ 0.2

L . . . . . .

.4

I I

I t

.

,

; I

I I

.

i

I I

I I

.

I

•

~ 0.4

--

-2 -~.g3....;-- ,~. . . . . . .....

I I

I I

_1 . . . . . I

0

~

i

t

-6

I I

1

•

i

I I

I

0.4 0.2 0.8

]

I I

"

*er,

L . . . . . . .

•l

~PIiB

oi'~,o•

1

1

I

1

1

i

I

I

~

""

,

-4

i

J

I l

~-1 ......

-60

-80

1 1 -100 0 -20 -40 Pressure head 1/I(cm) 10.20 cm • 20.30 cm - 30-40 cm •

1 -60

-80

40-50 cm

Fig. 2 Water retention 0(V) and hydraulic conductivity K(V) curves observed for the soil sampling locations A through F. The saturated conductivity K~ was fixed at the observed values. 2)

A p p l i c a t i o n o f L e i b e n z o n m o d e l ( m e t h o d 3)

The Leibenzon model was applied to the observed conductivity data sets by using the following two methods. In method 3a, ~ in Eq. (5) was optimized to minimize the RSS comparing measured v e r s u s computed log K values for the entire data sets (that is, a single ,/~ value was used for the entire data set for 30 soils), whereas in method 3b, ~ was optimized for each soil. Both methods used the same 0~ and Or values determined by the LN model in method l, and K, was fixed at the measured values.

Results a n d D i s c u s s i o n 1 Spatial variation of m e a s u r e d soil hydraulic properties

Figure 2 shows the measured water retention, 0(~/), and hydraulic conductivity, K(~), curves for all 30 soils. Water retention curves for the locations B through F have large changes in 0 in the range of 0 > V > -- 30 cm. Such changes are typical for forest soils (Kosugi, 1994). These water retention curves indicate existence of macropores. Water retention

curves for surface soils (0-30 cm deep) tend to have greater changes in 0 in the range of 0 > ~ / > -- 30 cm than subsurface soils (30-50 cm deep). Water retention curves for 0-10 and 30-50 cm deep of the location F have relatively small 0 values. At the crest (location A), retention characteristics are dif-

Hendrayanto et al.

//

0

0.8

(b)Methodla ~ ¢ ~

6

0.6

4

0.4

2

}.6

0

0.2

111

-2 RSS 0 1787

0.0 0.0

0.2

0.4

0.6

-10

/

-10

0.8

RSS= 2090.65 -8

-6

-4

Oobs 0

(c) Method~

0

-6

IOgXob~

-8

-2

-2

°~-~

-8

~

-10

-10

-S

-6

-4

-2

0

tOggobs

-4

-2

0

Methodla Method lb Method lc

2

4

6

8

Fig. 4 Relationshipbetween I values for the LN (I-LN)and VG (/-VG) models used in methods la through lc.

-8

RSS= 1684.39

-6

~

/-LN

_o.6

-6

o / ~ o

-8

(d) Method1¢ a , ~ o ~ ' /

-2

/

-4

-10 -10

-8

-6

-4

-2

0

IogKobs

3 Relationships between observed (0ot,.0and estimated (0e.,.3water contents, and observed (Koh.O and estimated (Ke~t) hydraulic conductivities resulting from the separate parameter optimization (method 1) using the LN model. Notes: Method la used the constant l of 0.5. Method lb used the constant I of - 0.323. In method lc, I was optimized for each soil.

Table 1 Residual sum of squares (RSS) for 0 and log K derived from separate (method 1) and simultaneous (method 2) parameter optimizations of the LN and VG models. Model

Fig.

ferent from those at slope regions; 0(N) curves decrease continuously in 0 over the entire range of vespecially for the upper layers (0-30 cm). At every location, K, is larger for surface soils than for subsurface soils. Surface soils have larger decreases in K with decreasing in gcompared with subsurface soils. As a result, K values in dry region are smaller for surface soils than those for subsurface soils. This trend is most typical for the upper and middle portions of the slope (locations B, C, and D). Vertical variations in K(V) curves for locations A and F are relatively small in spite of the relatively large vertical variations in 0(N) curves for these locations. The soils at the foot slope (location F) have larger K,. values and greater changes in K than the soils at the crest (location A). 2 Applications of the soil water retention and hydraulic conductivity models 1) The LN and VG models Figure 3a shows the relationship between measured and estimated 0 using method 1 for the LN model. The data fall closely around the 1:1 line which denotes where the measured and estimated values are equal. This means that the retention model expressed as Eq.(1) can express the observed water retention curves adequately. Figures 3b, c, and c show the relationships between measured and calculated log K from the LN model by using methods la, lb, and lc, respectively. The points in Fig. 3b, where K values were calculated using l : 0.5 (Method la), were scattered widely around the 1:1 line, which indicates that the predictive model does not adequately describe the observed K values. Many data points were below the 1:1 line, indicating that the model generally under-

LN VG

RSS (0) RSS (log K) Method I Method2 Method la Method Ib Method lc Method2 0.1787 0.5312 2090.7 1684.4 369.5 107.6 0.1949 0.5704 1625.9 1553.2 335.0 110.4

estimates K(V). Similar results were shown by Yates et al. (1992). Optimizing the I value for all soil hydraulic conductivity data (method lb) resulted in an I value smaller than 0.5 (1 = - 0.323) because the l value of 0.5 tended to underestimate the K(V) data. However, the RSS was still large and the data points were scattered far from the 1:1 line (Fig. 3c). The results of methods la and lb indicate that predicting K(V) values from retention data only could produce large errors in numerical modeling of soil water flow. Optimizing I value for hydraulic conductivity data of each soil (method l c) greatly improved the estimation of K (Fig. 3d). The RSS was reduced significantly, and the data points were close to 1:1 line compared with Figs. 3b and c. Many of the data points fell within approximately one order of magnitude. Consequently, I should be treated as a fitted parameter for accurate descriptions of K(V) curves. Table 1 summarizes the RSS values resulting from fitting the LN and VG models to the observed O(V) and K(V) data. The results derived by the VG model were similar to those derived by the LN model. This corresponds to the results of Kosugi (1996), which showed that both LN and VG models generally produce similar fitting results. Figure 4 shows the relationship between I values for the LN (/-LN) and VG (/-VG) models used in methods la through lc. The optimized 1 values by method lb ( - 0.323 and - 0.009 for the LN and VG models, respectively) were smaller than the I values of 0.5 used in method 1a. Points for the optimized 1 values by method lc (open circles in Fig. 4) are close to the 1:1 line, indicating that the/-LN and/-VG for each soil are similar to each other. The/-LN values for method lc ranged from -- 5.63 to 5.00 with a standard deviation of 2.20, and IVG ranged from - 6.18 to 5.82 with a standard deviation of 2.72. Such large variations in 1 values were also demonstrated by Wtsten and van Genuchten (1988), Schuh and

112

J. For. Res. 4 (2) 1999:

0.8

0 (b) -2

0.6

,%.°.s~ * ,,*' ° ~

OA

16 14 12

................................................

'i-'"o"

................................................

÷ ........

lo 8

.................................

...................................

O~ L

0.2 0.4

0.2

0.6

0.8

-10 -10 -8 -6

-4 -2

0

/

/

-10 -10 -8

~

°

Co)Method3b ~,.,~N~'~,s~.,oaor ,~..4 • -2" _o-6.

° l

° v o

-8 • / / ° * v - o RS°S 412.28 -10 -10 -8 -6 -4 -2 0

RSS= 1927.14

-6 -4

-2

logKobs

. . . . .

.~. ......

i

IogKobs

Fig. 6 Relationshipsbetween observed (Kol,.~) and estimated (K~st) hydraulic conductivitiesobtainedby using the Leibenzonmodel. Notes: Method 3a used the constant/3= 4.13. In method 3b,/3 was optimized for each soil. Cline (1990), and Yates et al. (1992). Schuh and Cline (1990) showed that there was no predictive relationship between 1and the parameters of the VG model. The same was true for the forest soils analyzed in this study; no clear relationships between l and the parameters of the LN and VG models were found. Figures 5a and b show the relationships between the measured and estimated 0 and K values, respectively, resulted from the simultaneous optimization (method 2) using the LN model. Although the RSS of log K obtained by method 2 was higher than that by method 1 (Table 1), the data points are still close to 1:1 line (Fig. 5a). The RSS of log K obtained by method 2 was much smaller than that of method lc (Table 1), and most points in Fig. 5b fell within one order of magnitude of the 1:1 line. Similar results were obtained for the VG model (see Table 1). In conclusion, the simultaneous optimization procedure produces the best descriptions of hydraulic properties using the LN and VG models. Therefore, when both retention and unsaturated conductivity data are available for the parameter estimation, simultaneous optimization procedure should be used. 2) The Leibenzon model Figure 6a shows the relationship between measured and estimated log K obtained by using method 3a, where/6 of the

........

~ .........

.

.

!..,..~.~

.

............................

.

.

-10 -20 -30 -40 -50 -60 Ip'=(cm)

IogKobs

-2 (a)M e t h o d ~

o~

0

0

Fig. 5 Relationshipsbetween observed(O,,t,x) and estimated(0~.,r)water contents, and observed (Koh0 and estimated (K~t) hydraulic conductivities resultingfrom the simultaneousparameteroptimization(method 2) using the LN model.

-6

2

RSS= 107.57

Oobs

o

:!O

.........

4

-8

0.0 0.0

ii: .......... i "..i

Fig. 7 Relationshipbetween ~,, of the LN model and/3of the Leibenzon model. Leibenzon model (Eq. (5)) was optimized to all conductivity data (the optimized ,8 value was 4.13). Data were scattered far from the 1:1 line and RSS was large. The differences between observed and estimated K values were more than two orders of magnitudes for many data sets. Optimizing the ,8 value for hydraulic conductivity data of each soil (method 3b) greatly improved the estimation of soil hydraulic conductivity (Fig. 6b). The RSS of log K was reduced significantly and the data points were much closer to the 1:1 line. Most of the data fall within approximately one order of magnitude. The RSS for method 3b was slightly larger than RSS for method lc (see Fig. 3d and Table 1). Figure 7 shows the relationship between the optimized ,8 values for the 30 soils and Vm of Eq. (1) determined by method 1. Although Averjanov and Irmay suggested using constant ,8 values of 3.5 and 3.0, respectively Mualem (1978), Mualem (1978) demonstrated that ,8 varied from 1 to 24.5 and ]3 was larger for fine-textured soils than for coarse-textured soils. In Fig. 7, the values of ,8 for the 30 forest soils ranged from 1.64 to 13.80 with the mean value of 5.4 and the standard deviation of 2.38. Moreover, it is shown that ,8 tends to increase as ~,, decreases. That is, the soil with smaller median pore-radius tend to have a larger ,8 value. 3 Analyses of the spatial distribution of soil hydraulic properties using the LN model Based on the results of model applicability discussed in the previous section, the simultaneous optimization procedure employing the LN and VG models should be used to provide accurate description of soil hydraulic properties at the forested hillslope. Comparing the two models, only the LN model uses the physically-based parameters related to the statistics of soil pore-size distributions. Therefore, the parameters of the LN model obtained by method 2 were used in the following analyses of the spatial variability of forest soil hydraulic properties. Table 2 summarized the mean and standard deviation of the hydraulic parameters. Statistics of log ( - Nm) and log K,. were shown because these parameters were lognormally distributed rather than normally distributed based on the Kolmogorov-Smirnov test. The test suggested that the normal distribution adequately fit the 0, Or, 0 , - Or, and O"distributions. The coefficient of variation (CV), which is defined as

Hendrayanto et al.

113

Table 2 The mean and standarddeviationof simultaneously-optimized parameters (method2) of the LN model. LN model parameters 0~ Or

log[

O~ - - Or

I

in cm) ol log K~.(K~.in cm/s)

Mean 0.5928 0.3141 0.2787 1.16 1.23 0.22 - 1.11

Standard deviation 0.1200 0.0732 0.1155 0.29 0.38 1.38 0.66

la) Or Or

O.I 0.2

4) 20 40 30

60

90

(b) Iog(-~,.) (~,,, in cm)

! 20 L •i 0.75

0

; , 1.00 1.25 1.511

20 40 60

3O

the ratio of the standard deviation to the mean value, varied from 0.2 to 0.25 for the parameters 0,, Or, log (-- v/m),and a. The CV for the effective porosity, 0 ~ - Or, was about 0.4. The maximum and minimum values for l were 4.79 and - - 1.59, respectively, and maximum and minimum log K~ values were - 0.38 (K~.= 4.2 × 10 I cm/s) and - 3.25 (K, = 5.61 × 10 - 4 cm]s), respectively. The spatial distributions of 0,, - Or, l o g ( - v/m), O', log K,., and I are presented in Fig. 8. The effective porosity, 0,. - Or is large near the crest (the location A) and small in the lower portion of the slope (E and F). Effective porosity generally decreases as soil depth increases. Values v/m are generally small (that is, l o g ( - V/m)is large) at the crest and upper slope segments (A, B, and C) and large at mid-slope to footslope locations (D, E, and F). These differences coincide well with results derived by Kosugi (1997a) who showed v/m tended to be small for drier forest soils. Except at the crest (A), V/,~ tends to decrease as soil depth increases. The value of v/m is closely related to the soil structure (Kosugi, 1997a). Large v/m values at the locations D, E, and F and at surface layers at locations B and C, indicate the existence of a well-developed crumb structure in forest soils. Width of pore radius distribution (or) generally small at the crest and upper slope (A and B). Generally, cr is greater than 1, suggesting relatively a large width of pore radius distribution which is typical for forest soils (Kosugi, 1997a). Saturated hydraulic conductivity (K~) is generally small at the crest and upper slope segments (A, B, and C) and large at mid-slope to footslope locations (D, E, and F). At every location, K,. decreases as the soil depth increases. The value of K.~.is expected to be larger as the effective porosity becomes larger. Moreover, Kosugi (1997b) showed that K,. is inversely proportional to the square of v/re. The value of K, is small at sub-surface layers at the location B and C because of small 0~ -- Or and v/m values. Small v/m values resulted in small K~ values at the crest of the slope in spite of the high 0, - Or values. Large v/m values resulted in large K, values at mid-slope to footslope locations (D, E, and F). The pore tortuosity parameter (1) for samples from the crest and upper slope segments (locations A-C) is large in the soil surface and decreases with depth. At the mid-slope to lower slope locations (locations D, E), l tends to increase with increasing soil depth. At the footslope, 1 is uniform. The I values for the surface soils (0-30 cm deep) are more uniform than those for the subsurface soils (30-50 cm deep), indicating more uniform soil-pore tortuosity for the surface soils.

11.4

0.3

90

(c~ t~ 0

120 0

1'.0

.1.5

2O 40 30

60

90

120

Id) log/~ (K~ in cm/s) -2,5 .2.0 .I.5 .I.0 -0,5

0 20 40 3O

60

90

(e) l 0

-4.0

120 -2.0

0.0

2.0

4.0

20 4O 3O

60

90

120

Distance f r o m the crest ( m ) A

B

C

D

E

F

Soil sample locations

Fig. 8 Spatial distributions of (a) effective porosity 0~.- Or, (b) log ( - Vm) of the LN model, (c) o"of the LN model, (d) log K,., and (e) pore-tortuosityparameter I.

Conclusions

Soil hydraulic properties showed a considerable spatial variation at the forest hillslope in a headwater catchment of the Sumiyoshi River Basin. Applying the combined water-retention-hydraulic-conductivity models proposed by Kosugi (1996) and van Genuchten (1980) (the LN and VG models, respectively), it was shown that the pore-tortuosity parameter 1 introduced by Mualem (1976) should be treated as a fitted parameter for accurate descriptions of hydraulic conductivity curves. Moreover, it was concluded that the simultaneous optimization of both water retention and hydraulic conductivity curves is preferable to obtain appropriate descriptions of hydraulic properties when both retention and conductivity data are available for the parameter estimation. Estimations of hydraulic conductivity using Leibenzon's model showed that the exponent parameter had a spatial variation with the standard deviation of 2.38. The Leibenzon's model produced slightly poorer estimations than the LN and VG models. The parameters of the LN model obtained by the simultaneous optimization procedure were used for analyses of spatial variation in soil hydraulic properties. The parameter v/m, which has a positive correlation with median pore-radius, was generally small at crest and upper slope and large at mid-slope to footslope locations. Except for the crest, surface

114 soils had larger grin values than the subsurface soils, suggesting the existence of developed crumb-structure forest soils. For most soils, cr was greater than 1, indicating relatively large width of pore-size distribution which is typical for forest soils. Saturated hydraulic conductivity, Ks, was generally small at crest and upper slope and large at mid-slope to footslope locations. Although the soils at crest had large effective porosity values, K~. values for these soils were small because of small gtm values. The larger K~. values at mid-slope to footslope locations were attributable to the larger ~m values. The considerable spatial variability in soil hydraulic properties presented in this study is important to be considered in water-flux studies based on Darcy-Buckingham equation. However, Darcy-Buckingham equation does not account for water-flux in connected macropores or pipe flow which reportedly can contribute to storm discharge up to 99.5% (Noguchi et al., 1997). To better understand hillslope scale discharge phenomena, a combined micropore-macropore system should be analyzed. L i t e r a t u r e cited Beckett, P.H.T. and Webster, R. (1971) Soil variability: a review. Soils Fert. 34: 1-15. Beven, K. and Germann, P. (1982) Macropores and water flow in soils. Water Resour. Res. 18:1311-1325. Bonell, M. (1993) Progress in understanding runoff generation dynamics in forests. J. Hydrol. 150: 217-275. Burdine, N.T. (1953) Relative permeability calculation from size distribution data. Trans. Am. Inst. Min. Metall. Pet. Eng. 198: 71-78. Buttle, J.M. and House, D.A. (1997) Spatial variability of saturated hydraulic conductivity in shallow macroporous soils in a forested basin. J. Hydrol. 203: 127-142. Clausnitzer, V., Hopman, J.W., and Nielsen, D.R. (1992) Simultaneous scaling of soil water retention and hydraulic conductivity curves. Water Resour. Res. 28: 19-31. Freeze, R.A. (1975) A stochastic conceptual analysis of one-dimensional groundwater flow in nonuniform homogeneous media. Water Resour. Res. 11:725-741 Hendrayanto, Kosugi, K., and Mizuyama, T. (1998) Field determination of unsaturated hydraulic conductivity of forest soils. J. For. Res. 3: 11-17. Hills, R.G., Hudson, D.B., and Wierenga, P.J. (1989) Spatial variability at the Las Cruces Trench site. In Proceeding of the international workshop on indirect methods for estimating the hydraulic properties of unsaturated soils, van Genuchten, M. Th. and Leij F.J. (eds.), 718 pp, Riverside, California, 529-538. Hoeksema, R.J. and Kitandis, P.K. (1985) Analysis of spatial structure of properties of selected aquifers. Water Resour. Res. 21: 563-572. Jury, A.W., Russo, D., Sposito, G., and Elabd, H. (1987) The spatial variability of water and solute transport properties in unsaturated soil. I. Analysis of property variation and spatial structure with statistical models. Hilgardia 55(4): 1-32. Kitahara, H. (1994) A study on the characteristics of soil pipes influencing water movement in forested slopes. Bull. For. For. Prod. Res. Inst. 367:63-115. Kosugi, K. (1996) Lognormal distribution model for unsaturated soil hydraulic properties. Water Resour. Res. 32: 2697-2703. Kosugi, K. (1997a) A New model to analyze water retention characteristics of forest soils based on soil pore radius distribution. J. For. Res. 2: 1-8. Kosugi, K. (1997b) New diagrams to evaluate soil pore radius distribution and saturated hydraulic conductivity of forest soil. J. For.

J. For. Res. 4 (2) 1999 Res. 2: 95-101. Kosugi, K. (1997c) Effect of pore radius distribution of forest soils on vertical water movement in soil profile. J. Jpn. Soc. Hydrol. Water Resour. 10: 226-237. Kubota, J., Fukushima, Y., and Suzuki, M. (1987) Observation and modeling of the run-off process on a hill slope. J. Jpn. For. Soc. 69: 258-269. (in Japanese with English summary) Kutflek, M. and Nielsen, D.R. (1994) Soil hydrology. 370pp, Catena Verlag, Cremlingen, Germany. Luckner, L., van Genuchten, M. Th., and Nielen, D.R. (1989) A consistent set of parametric models for the two-phase flow of immiscible fluids in the subsurface. Water Resour. Res. 25: 2187-2193. Mohanty, B.P., Ankeny, M.D., Horton, R., and Kanwar, R.S. (1994) Spatial analysis of hydraulic conductivity measured using disc infiltrometer. Water Resour. Res. 30: 248%2498. Mualem, Y. (1976) A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12:513-522. Mualem, Y. (1978) Hydraulic conductivity of unsaturated porous media: generalized macroscopic approach. Water Resour. Res. 14: 325-334. Nielsen D.R., Biggar, J.W., and Erh, K.T. (1973) Spatial variability of field-measured soil-water properties. Hilgardia 42:215-259. Nogucbi, S., Tsuboyama, Y., Sidle, R.C., and Hosoda, I. (1997) Spatially distributed morphological characteristics of macropores in forest soils of Hitachi Ohta Experimental Watershed, Japan. J. For. Res. 2: 207-215. Ohte, N., Suzuki, M., and Kubota, J. (1989) Hydraulic properties of forest soils (I) The vertical distribution of saturated-unsaturated hydraulic conductivity. J. Jpn. For. Soc. 71: 137-147. (in Japanese with English summary) Osaka, A. (1996) Land form characteristics measurements. In Hydrogeomorphology: The interaction of hydrologic and geologic process. Onda, Y., Okunishi, K., Iida, T., and Tsujimura, M. (eds.), 267pp, Kokonssoi, Tokyo. (in Japanese) Peck, A.J., Luxmoore, Janice, R.J., and Stolzy, L. (1977) Effect of spatial variability of soil hydraulic properties in water budget modeling. Water Resour. Res. 13: 348-354. Sammori, T. and Tsuboyama, Y. (1990) Study on method of slope stability considering infiltration phenomenon. J. Jpn. Soc. Erosion Control Eng. 43: 14-21. (in Japanese with English summary) Schuh, W.M., and Cline, R.L. (1990) Effect of soil properties on unsaturated hydraulic conductivity pore interaction factors. Soil Sci. Soc. Am. J. 54: 1509-1519. Shinomiya, Y., Kobiyama, M, and Kubota, J. (1998) Influence of soilpore connection properties and soil pore distribution properties on the vertical variation of unsaturated hydraulic properties of forest slopes. J. Jpn. For. Soc. 80:105-111. (in Japanese with English summary) Tani, M. (1985) Analysis of one-dimensional, vertical, unsaturated flow in consideration of runoff properties of a mountainous watershed. J. Jpn. For. Soc. 67: 449-460. (in Japanese with English summary) Uchida, T., Kosugi, K, Ohte, N., and Mizuyama, T. (1996) The influence of pipe flow on slope stability. J. Jpn. Soc. Hydrol_ Water Resour. 9: 330-339. Van Genuchten M.Th. (1980) A closed-form equation for predicting the hydraulic conductivity of unsaturated sols. Soil Sci. Am. J. 44: 892-898. Warrick, A.W. and Nielsen, D.R. (1980) Spatial variability of soil physical properties in the field. In Application of soil physics. Hillel, D. (ed.), 385pp, Academic Press, New York, 319-344. Watson, K.K. (1966) An instantaneous profile method for determining the hydraulic conductivity of unsaturated porous materials. Water Resour. Res. 2: 70%715. W6sten, J.H.M. and van Genuchten, M.Th. (1988) Using textures and other soil properties to predict the unsaturated soil hydraulic functions. Soil Sci. Soc. Am. J. 52: 1762-1770. Yates, S.R., van Genuchten, M.Th., and Lij, F.J. (1992) Analysis of measured, predicted, and estimated hydraulic conductivity using the RETC computer program. Soil Sci. Soc. Am. J. 56: 347-354. (Accepted December 25, 1998)