Characterizing heterogeneity of soil water flow by dye infiltration experiments

Journal of Hydrology (2006) 328, 559– 571

available at www.sciencedirect.com

journal homepage: www.elsevier.com/locate/jhydrol

Characterizing heterogeneity of soil water flow by dye infiltration experiments Kang Wang

a,c

, Renduo Zhang

a,b,*

, Hiroshi Yasuda

c

a

State Key Laboratory of Water Resources and Hydropower Engineering Science, Wuhan University, Wuhan 430072, China School of Environmental Science and Engineering, Sun Yat-Sen (Zhongshan) University, Xinggang Xilu 135, Guangzhou 510275, China c Arid Land Research Center, Tottori University, 1390 Hamasaka, Tottori 680-0001, Japan b

Received 23 April 2005; received in revised form 5 September 2005; accepted 2 January 2006

KEYWORDS Flow heterogeneity; Measurement scale; Dye experiments; Random cascade model

Summary The inherent heterogeneity of soil water movement, such as fingering flow, is still poorly understood. To address the flow heterogeneity issue, the objectives of this study are to: (1) investigate heterogeneity of soil water flow at different measurement scales and under different boundary conditions using dye infiltration experiments; (2) characterize the heterogeneity information included in different scales using the random cascade model, and (3) predict heterogeneous soil water flow over the measurement scales. Field experiments of dye infiltration included three measurement scales (three dye source surface areas of 25 · 25, 50 · 50, and 100 · 100 cm) and three hydraulic boundary conditions (three initial water ponding depths: 1.5, 2.5, and 5.0 cm at the soil surface). The random cascade model with a lognormal distribution was used to simulate the infiltration process in soils and different methods were applied to estimate the model parameters. Results based on the experimental data and the model simulations showed that the hydraulic boundary condition and the measurement scale were important factors to affect the flow patterns in soils. To accurately describe flow transport processes at different scales, it was necessary to consider heterogeneities in the vertical and horizontal directions. As the measurement scale increased, the effect of multi-dimensional heterogeneities on the flow processes in soils became more significant. Interdependency between flow paths also increased with the scale. ª 2006 Elsevier B.V. All rights reserved.

* Corresponding author. Tel.: +86 20 31399501; fax: +86 20 84110692. E-mail address: [email protected] (R. Zhang). 0022-1694/$ - see front matter ª 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2006.01.001

560

K. Wang et al.

Introduction Because of soil heterogeneity and the interaction between soil and water, quantitatively describing soil water flow remains a challenging undertaking in the vadose zone hydrology. Highly heterogeneous transport phenomena have been observed from field experiments and the flow heterogeneity increases with the measurement scale (Hillel, 1987; Molz ¨ hrstro and Boman, 1993; O ¨m et al., 2002; Wang et al., 2003). It is desirable to predict water flow and solute transport in large scales from observations at smaller scales. Therefore, determination of scale-dependent variability of soil water flow becomes essential to quantitative prediction of the flow process in the field. Attributable to good visibility in soils and moderated costs, experiments using dye as a tracer are a useful tool for revealing flow patterns in soils and have been used in soil science for decades (Bouma and Dekker, 1978; Persson et al., 2001; Yasuda et al., 2001). Dye tracer data are able to provide more detailed visual information in the multidimensional flow domain. In recent years, significant effort has been made in developing a variety of approaches to deal with heterogeneous flow patterns, such as fingering flow in soils. These approaches mainly fall into two categories: the continuum approach and the discrete approach (Liu et al., 2003). In the continuum approach, connected macropores and soil matrix are viewed as two or more overlapped interacting continua, which are considered to coexist in a system. In such approach, the continuum mechanics formulations are used to describe flow moments in the continua. The discrete approach treats the soil water system as discrete ‘‘particles’’ or ‘‘packages’’, whose movement is controlled by simple rules to generate growth patterns with rich structure. The random cascade model (Olsson et al., 2002) and diffusion limited aggregation model (Flury and Flu ¨hler, 1995; Persson et al., 2001) belong to the discrete category. The continuum and discrete approaches have advantages and disadvantages. The traditional continuum method assumes uniformly distributed flow patterns at a subgrid scale, therefore cannot be used for representing gravity-driven flow, resulting from the soil heterogeneity. Several models have been developed to incorporate heterogeneous flow behaviors into the continuum approach. The dual-continuum approach, in which matrix and fractures are treated as two overlapping and interactive continua, is used as the baseline approach for modeling water flow and transport. These models are superior to the traditional continuum approach for flow transport in heterogeneous soils. However, a major limitation of the currently available continuum models is that

Table 1

they were mainly developed for modeling flow with the effect of soil macropores. It is difficult for such models to capture the fractal flow patterns (Liu et al., 2003). Fractal flow patterns have been often observed in variably saturated flow systems (Glass, 1993; Smith and Zhang, ¨ hrstro 2001; O ¨m et al., 2002). Persson et al. (2001) showed that observed mean power spectrum from a dye experiment in the field displayed a typical power-lower relation, an indication of flow’s fractal behavior. The discrete approach can be used to generate fractal patterns of soil water flow. Cascade-based modeling has been proved to be highly efficient for simulating heterogeneous hydrological processes (Olsson, 1995; Liu and Molz, 1997) and soil infiltration processes can be well represented by a random cascade model (Olsson et al., 2002). However, the main problem of the discrete approach is that satisfactory theories underlying the approach are still missing. Heterogeneous flow paths depend on soil properties, experimental conditions (e.g. hydraulic boundary conditions), and measurement scales (or support sizes). These issues have not been well addressed in the cascade-based modeling. To better understand flow patterns in soils and some physical insights of the random cascade model, in this study we focused on investigating flow heterogeneity information included in different measurement scales and effects of hydraulic boundary conditions on the flow heterogeneity, using dye infiltration experiments and the random cascade model. Specifically, the objectives of this paper are to: (1) investigate heterogeneity of soil water flow at different measurement scales and under different boundary conditions using dye infiltration experiments; (2) characterize the heterogeneity information included in different scales using the random cascade model, and (3) predict heterogeneous soil water flow over the measurement scales. The predicted results were compared with the experimental data.

Materials and methods Study site description Dye infiltration experiments were carried out in the Arid Land Research Center (3532 0 N, 13413 0 E), Tottori University, Japan in 2004. Two years prior to the experiments, potato was grown at the study site. After potato harvest, the site was fallow. Except some residual roots within the depth of 20 cm, there was no visible heterogeneity in the sandy soil. The natural groundwater table was deeper than 5 m. The soil was Arenosol (silicious sand, Typic Udipsamment) (Qiu et al., 1999). After photographing dye stained patterns

Soil properties in the field experiment

Depth (cm)

Sand (%)

Silt (%)

Clay (%)

Bulk density (g cm3)

Porosity (cm3 cm3)

Saturated hydraulic conductivity (m s1)

0–10 10–20 20–50 50–100

88 96.5 100 100

4 0.2 0 0

8 3.3 0 0

1.55 1.52 1.43 1.51

39.53 40.21 38.67 38.55

3.30 · 105 2.66 · 104 3.76 · 104 3.76 · 104

Characterizing heterogeneity of soil water flow by dye infiltration experiments of soil layers (see below), undisturbed soil samples of 100 cm3 were collected at depths of 0–10, 10–20, 20–50, and 50–100 cm. The samples were used to measure soil properties, including soil texture, bulk density, porosity, and saturated soil hydraulic conductivity, which are summarized in Table 1.

Experimental setup The experiment design is schematically shown in Fig. 1. In a field of 1 ha, 18 plots were set up with three inner surface areas of 25 · 25, 50 · 50, and 100 · 100 cm (or three measurement scales of dye source). For each plot, the dimension of outer frame was 200 · 200 cm. To avoid interference, the plots were separated at least 3 m apart. The outer and inner frames were inserted into the soil to a depth of 20 cm. Before the experiment, the plot surface was leveled with vegetation removed. One week prior to the experiment, we started to saturate the plots. At the first three days, the plots were irrigated for several times each day without ponding. From the fourth day, a ponding head of 5.0 cm was kept at the soil surface for about 3 h everyday. Soil within the inner and outer frames was covered with plastic film. At the beginning of the dye experiment, the inner plots were prepared by ponding with 1.5, 2.5, and 5.0 cm depths of dye solution on the top of the plastic film as three different boundary conditions, and the outer area was ponded with the same depths of water. Two replicates were performed for each of the nine treatments. Brilliant blue FCF has been commonly used for visualizing heterogeneous flow paths (Olsson et al., 2002; ¨ hrstro O ¨m et al., 2002, 2004). Therefore, the food-grade dye pigment Brilliant Blue FCF (Kiriya Chemical Co. Ltd.) was chosen as the tracer of our experiments. The tracer is readily soluble in water and the water solution gives a clearly visible blue staining to the soil. According to the previous research (Flury and Flu ¨hler, 1995), the concentration of dye solution in our experiments was adjusted to 4 kg m3. When the experiments started, the plastic film of plots was removed immediately, creating instantaneous ponding infiltration, lasting for 5–10 min. Then the plots were covered with plastic sheets to prevent evaporation and left for 24 h to complete the dye infiltration process. Finally the soil in the plots was removed layer by layer. For the first 60 cm depth, the increment between layers was 1 cm and then the increment became 5 cm to the depth where dye coverage became invisible. Dye distributions in the layers were recorded using a digital camera. The photographs were taken during day time. A white sheet was placed over each plot to diffuse the light and to avoid direction radiation (Forrer et al., 2000).

Figure 1

561

Image processing The color of each pixel was transformed to black (dye covered) and white (non-dye covered) on the basis of red, green, and blue values by the following two steps: adjusting colors and setting threshold values. Since the photographs were taken at different time during the day, they showed different tints even for the same object. Therefore, the photographs were adjusted to a common standard following the procedure by Forrer et al. (2000). The threshold algorithm by Baveye et al. (1998) was utilized to transform color images into black and white (binary) images. A trial and error method was used to obtain reasonable threshold values. Gray scale levels of colors customarily range from 0 to 255, with white and black corresponding to gray scale values of 0 and 255, respectively. First, a value between 0 and 255 was chosen as the threshold value. Second, the mean value of pixels with gray scale levels greater than the threshold value was calculated. The mean value of pixels with gray values less than and equal to the threshold value was also calculated. The average of the two mean values was used as the next threshold value and this process continued until no change of sequential threshold values. For each layer, this process was preformed independently. The pixel size of transformed black and write images was approximately 1.0 mm2.

Theory of random cascade simulation The dye infiltration depth can be described using a random cascade model as follows (Olsson et al., 2002; Schmitt, 2003): Zðx; yÞ ¼ ZA eðx; yÞ

ð1Þ

where Z(x, y) is the infiltration depth at a location with coordinates (x, y), ZA is the mean infiltration depth, and e(x, y) is the random cascade function. The construction of a random cascade begins with the mean value of dye infiltration depth over a two-dimensional (d = 2) domain. During each subdivision level of the domain, the infiltration depth at each cascade obtained from the previous step is distributed into b sub-areas (b = 2d = 4) according to a set of cascade generators W (Fig. 2). The random cascade generators are positive and independently distribute with the following property: EhWii = 1, where Eh Æ i denotes expectation. As shown in Fig. 2, at the second subdivision level, the average dye infiltration depth of the first level is divided into four sub-areas based on random cascade generators denoted by W1(x, y). At the third level, each of the above sub-areas is further subdivided into b sub-areas and the random cascade generators for the b2 sub-areas are

Diagram of the experimental setup.

562

K. Wang et al. According to the properties of random cascade generators, the subscript i at the right hand side of Eq. (3) can be any value from 1 to n. The scaling behavior of e(x, y) with k1 (or the total scale ratio) can be investigated for various q values in the form of Eheðx; yÞq i ¼ ½k1 KðqÞ

ð4Þ

n

where k1 = k and K(q) is the cumulative generating function. Eq. (4) can be rewritten as: Eheðx; yÞq i ¼ ½kn KðqÞ ¼ ½kKðqÞ n

ð5Þ

Substituting Eq. (5) into Eq. (3) yields: KðqÞ ¼ logk EhW i ðx; yÞq i

ð6Þ

The conservative property EhWi(x, y)i = 1 gives K(1) = 0 and Ehe(x, y)i = 1. Eq. (4) has the generic property of multifractal fields, obtained through multiplicative cascades. The universal multifractal model is able to reproduce the fractal dye behavior reasonably well (Schertzer and Loverjoy, 1987; Olsson et al., 2002). The random cascade generator W can be any distribution that has a power law type tail (Boufadel et al., 2000), such as the Levy distribution and the Gaussian distribution. The Levy distributed generator was documented for water flow in clay soils, with values of the multifractal index parameter (a) close to 2 (from 1.84 to 1.95) (Olsson et al., 2002). In this study, the Gaussian distribution (a = 2), corresponding to a lognormal distribution of the cascade generator, is applied because the distribution includes the minimal number of parameters. The validity to use the distribution is to be verified latter. We take the form of Wi(x, y) = eg, where g is a Gaussian random variable with the mean of G and the variance of r2: g = rg0 + G, where g0 is a centered and unitary Gaussian random variable (with zero mean and unit variance). Then we have EhW i ðx; yÞq i ¼ Eheðrg0 þGÞq i ¼ er Figure 2 Discrete random cascade model with a subdivision number b = 4.

2 q2 =2þqG

ð7Þ 2

From EhWi(x, y)i = 1 and Eq. (7), we obtain G =  r /2. The finial form of the cascade generators Wi(x, y) is as follows: r2

W i ðx; yÞ ¼ erg0  2 denoted by W2(x, y). The relationship between the second and third levels is represented by W1(x, y) · W2(x, y). This process of the subdivision is continued further down the scales. The random cascade function of the dye infiltration depth at the smallest scale l0 is related to the cascade generators as follows: eðx; yÞ ¼

n Y

W i ðx; yÞ

ð2Þ

i¼1

where n is the number of subdivision levels at the smallest scale. The cascade is developed from the length scale L0 of a study area (L0 · L0) down to l0 = L0/kn, where k is a constant scale ratio between two consecutive scales. Since all random variables are independent, the moment of order q of the e(x, y) is expressed by (Schmitt and Marsan, 2001): Eheðx; yÞq i ¼

n Y i¼1

EhW i ðx; yÞq i ¼ ½EhW i ðx; yÞq in

ð3Þ

ð8Þ

Substituting Eq. (7) into Eq. (6) gives: KðqÞ ¼ logk EhW i ðx;yÞq i ¼ lnðer

2 q2 =2qr2 =2

l Þ=lnk ¼ ðq2  qÞ ð9Þ 2

where l = r2/ln k = K(2) is a constant. Eq. (8) is used to simulate random cascades with the lognormal distribution.

Parameter estimation for the random cascade model To carry out random cascade simulations, we need to estimate two parameters: the mean infiltration depth ZA and the standard deviation r of the cascade generator. In this study, three methods were used to determine the two parameters. In the first method (Method I), the parameters were estimated from the experimental dye infiltration data. The second method was a modification of Method I (Modified Method I). In the third method (Method II), dye coverage of each layer was used to fit the parameters.

Characterizing heterogeneity of soil water flow by dye infiltration experiments For Method I, measured dye infiltration depths for the experimental plots were obtained from the digitized photographs. For the plots with inner surface areas of 25 · 25, 50 · 50, and 100 · 100 cm, the length scales (L0) for the cascade development were 25, 50, and 100 cm, respectively. Taking the plot with L0 = 100 cm as an example, we transformed the photographs of each layer to a file including 1024 · 1024 pixels, each of which had a value of either 1 (dye covered, black) or 0 (non-dye covered, white). As a result, we obtained the dye pattern for each layer. Then these files were put together layer by layer in the same order as that in the field experiments, resulting three-dimensional dye-stained patterns. The final dye infiltration depth at each pixel was determined from the three-dimensional dye-stained patterns and used as input information for estimation of the cascade parameters in Method I. We used the random cascade model with 11 layers (n = 11) to simulate dye infiltration. As shown in Fig. 2, at the 11th level we have 210 · 210 sub-areas, corresponding to the 1024 · 1024 pixels of the final dye infiltration distributions from the experiments. The geometry average method was used to calculate the dye infiltration depth ZA

563

from levels 10 to 1 (Jothityangkoon et al., 2000). Values of Wi(x, y) were determined from the dye distributions at each level, from which e(x, y) was evaluated. Then the function of K(q) was calculated through Eq. (4) and r was estimated using Eq. (9). It should be noticed that in Method I it was assumed that the zone above the final dye infiltration depth at each pixel was completely covered by dye, which resulted in losing heterogeneity information in the vertical direction. Therefore, we introduced an apparent dye infiltration depth at each pixel, defined by: Za ðx; yÞ ¼

k X

where k is the number of sampling layers, Dzi is the distance between layer i and layer i + 1, b(x, y) is the index of dye coverage (if a pixel in the layer is covered with dye, b(x, y) = 1, otherwise b(x, y) = 0, which can be directly determined from the three-dimensional dye distributions). Different from Method I, Modified Method I uses the average value of the apparent infiltration depths to estimate ZA.

Dye coverage 0.2

0.4

0.6

Dye coverage 0.8

1

0

0

0

20

20

Depth (cm)

Depth (cm)

0

ð10Þ

bðx; yÞDzi

i¼1

40 Plot 1 Plot 2 Plot 3 Plot 4 Plot 5 Plot 6

60

80

0.2

0.4

0.6

0.8

1

40 Plot 7 Plot 8 Plot 9 Plot 10 Plot 11 Plot 12

60

80

100

100

a

b

Dye Coverage 0

0.2

0.4

0.6

0.8

1

0

Depth (cm)

20

40

60

80

100

Plot 13 Plot 14 Plot 15 Plot 16 Plot 17 Plot 18

120

c

Figure 3

Dye coverage versus depth for the plots with ponding water heads of (a) 1.5 cm, (b) 2.5 cm, and (c) 5.0 cm.

564 For Method II, the dye coverage of each layer was used to estimate ZA and r. With a chosen set of ZA and r, Wi(x, y) were calculated using Eq. (8) and then dye infil-

K. Wang et al. tration depth distributions were generated using Eqs. (1) and (2). The generated dye coverage of each layer was compared with the experimental results and correlation

Figure 4 Black and white images of flow patterns for plots with ponding water head of (a) 1.5 cm, (b) 2.5 cm, and (c) 5.0 cm. For each group of the three plots, the plot dye source areas are in the order 100 · 100, 50 · 50, and 25 · 25 cm. In the figures, z is the vertical distance starting from the soil surface.

Characterizing heterogeneity of soil water flow by dye infiltration experiments

565

Figure 4 (continued)

coefficients between the generated results and experimental data were calculated. By repeating the procedure with different sets of ZA and r, the values of ZA and r resulting in the highest correlation coefficient were considered the final estimated parameters. Random cascade realizations were different even if the same parameters were used. To determine how many random cascade realizations were required to arrive at stable results, 200

Table 2

Measurement scale (cm2)

1.5

100 · 100 50 · 50 25 · 25 100 · 100 50 · 50 25 · 25

5.0

Correlation analysis Fig. 2 shows the paths to construct e(x, y) at two points, (x1, y) and (x2, y). The paths of the two points were in

Flow heterogeneity information estimated from dye infiltration experiments

Ponding head (cm)

2.5

realizations were performed for each set of ZA and r values. Results showed that 40 realizations were sufficient to achieve good stability.

100 · 100 50 · 50 25 · 25

Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

Maximum depth (cm)

50% Dye coverage depth (cm)

Maximum depth/50% depth

75.0 77.0 56.0 34.0 53.0 32.0 92.0 114.0 86.0 77.0 52.0 58.0 117.0 124.0 108.0 122.0 124.0 119.0

16.0 11.0 14.0 9.0 16.0 16.0 32.0 25.0 18.0 17.0 11.0 4.0 56.0 45.0 39.0 55.0 44.0 45.0

4.69 7.00 4.00 3.78 3.31 2.00 2.88 4.56 4.78 4.53 4.73 14.50 2.09 2.76 2.77 2.22 2.82 2.64

566

K. Wang et al.

the same sub-areas from the cascade level 1 to n  m. At level n  m + 1, the two paths begin to separate, as indicated by the gray areas in n  m + 1 and n levels in Fig. 2. In the construction of e(x1, y) and e(x2, y), the values of Wi are identical for i from 1 to n  m and different from n  m + 1 to n. We use r = km to estimate the distance of two points in the x- (or y-) direction. The definition of cascade model gives: * Eheðx; yÞeðx þ r; yÞi ¼ E

n Y

W i ðx; yÞ 

i¼1

n Y

+ W j ðx þ r; yÞ

j¼1

¼ ðEhWðx; yÞ2 iÞnm ðEhWðx; yÞiÞ2m

ð11Þ

Introducing K(q) gives Kð2Þ

Eheðx; yÞeðx þ r; yÞi ¼ k1 rKð2Þ

ð12Þ

Eq. (12) can be used to calculate the correlation between two points in the x- (or y-) direction.

0.25

Plot 1

Results and discussion Dye stained patterns Fig. 3(a)–(c) show the dye coverage versus soil depth for the plots with the initial ponding water heads of 1.5 cm (Plots 1–6), 2.5 cm (Plots 7–12), and 5.0 cm (Plots 13– 18), respectively. Plots 1, 2, 7, 8, 13, and 14 had the inner dye surface areas of 100 · 100 cm, i.e., the largest measurement scale of dye infiltration. Plots 3, 4, 9, 10, 15, and 16 had the inner dye surface areas of 50 · 50 cm. Plots 5, 6, 11, 12, 17, and 18 had the inner dye surface areas of 25 · 25 cm, i.e., the smallest measurement scale of dye infiltration. Fig. 4(a)–(c), respectively, present black and white images of flow patterns for plots with 1.5 cm ponding water head (Plots 1, 3, and 5), 2.5 cm ponding head plots (Plots 8, 9, and 11), and 5.0 cm ponding head plots (Plots 13, 16, and 17). The figures show great variations from section to section and from plot to plot. The dye patterns exhibit fingering flow for all the measurement scales and the dye coverage decreases rapidly with depth. A complete

0.2

R 2 = 0.9889

0.2

0.15 R 2 = 0.9905

0.1

R 2 = 0.9940

0.05

Method I

0.05

Modified Method I

0.1

K(q)

K(q)

0.15

R2 = 0.9914

Method I

Plot 7

Modified Method I

0 -1 -0.05

1

3

5

0 -1 -0.05

7

2

(q -q )/2 0.25

0.35 Plot 9

1

3

5

7

2

(q -q )/2 R2 = 0.9868

Plot 11

R2 = 0.9830

0.25 R2 = 0.9899

Modified Method I

K(q)

K(q)

0.15

Method I

0.15

R2 = 0.9848

0.05

Method I

0.05 -1 -0.05

1

3

5

7

Modified Method I

-1 -0.05

1

2

3

5

7

2

(q -q )/2

(q -q )/2 0.16

K(q)

0.12

Method I

R2 = 0.9993

Modified Method I R2 = 0.9981

0.08 0.04

Plot 14 0 -1 -0.04

1

3

5

7

2

(q -q )/2

Figure 5 Relationships between K(q) (Eq. (6)) and (q2  q)/2.0 for some plots, whose slopes give standard deviation values of the random cascade model.

Characterizing heterogeneity of soil water flow by dye infiltration experiments dye coverage existed only within a few centimeters below the soil surface. Table 2 lists the maximum depth of dye infiltration and the depth where 50% of soil was covered with dye (50% dye coverage depth) in the plots. As shown in the table, for the plots with the initial water heads of 1.5 and 2.5 cm, the maximum infiltration depths of the 100 · 100 cm plots are greater than those of 50 · 50 and 25 · 25 cm plots. However, there seems no significant difference of both maximum depth and 50% dye coverage depth among the three measurement scales with the initial water head of 5.0 cm. The ratios of the maximum depth and the 50% dye coverage depth decreased with the initial water head for all the measurement scales. The results indicated that the hydraulic boundary condition was one of the factors contributing to the heterogeneity of soil water flow.

Random cascade simulation For Method I, the observed dye infiltration data were used to calculate Wi(x, y). The ZA value was obtained from the geometric average of the final dye penetration depths at all pixels. The value of e(x, y)q at each pixel was calculated using Eq. (6), taking q from 0 to 4 with an interval of 0.5 (Schmitt, 2003). Finally Eq. (9) was used to fit l, from which r was calculated. Fig. 5 shows some examples of the relationship between calculated K(q) versus (q2  q)/2, in which the slope was used to determinate r. Values of ZA and r determined from Method I are listed in Table 3 labeled as ð1Þ ZA and r(1). Modified Method I was applied in the same way as Method I but involving Eq. (10) and the resulting ð1Þ parameters are also listed in Table 3 as ZA and r(1)*. With the parameters from the methods, we generated dye infiltration depth distributions for the 18 plots. A comparison of observed and simulated dye coverage with Method I, Mod-

Table 3

Measurement scale (cm2)

100 · 100 50 · 50 25 · 25 100 · 100 50 · 50 25 · 25

5.0

Correlation analysis In Modified Method I, was determined from the dye distribution data. Therefore, we regarded the calculated correlations using e(x, y) and e(x + r, y) as the ‘‘observed’’ results. Since the correlations in the x- and y-directions are the same (Eq. (12)), the average correlation in the x- and y-directions was calculated based on the observed data

Method I ð1Þ

2.5

ified Method I, and Method II is presented in Fig. 6 for some examples. The best fitting parameters of Method II are also ð2Þ listed in Table 3 labeled as ZA and r(2). As mentioned before, the zone above the final dye penetration depths was not always covered by dye. By assuming that the zone above the final dye penetration depths was covered by dye, Method I ignored the vertical heterogeneity. Therefore, Method I overestimated the mean infiltration depth and underestimated the r value, as a result, overestimated the dye coverage for all plots (Fig. 6). Modified Method I overcame this problem to some extent and the dye coverage curves using the method were quite close to the observed data in some plots. As shown in Fig. 5, the slopes of the Modified Method I (solid lines) are greater than those of Method I (gray lines), which means that more heterogeneity information (or higher r values) was included using Modified Method I than using Method I. The best fitting r values using Method II, which included the heterogeneities in the horizontal and vertical directions, were much higher than those obtained with Method I and Modified Method I (Table 3). The results indicated that it is necessary to consider heterogeneities in the horizontal and vertical directions to accurately describe flow patterns in soils. The excellent linear relationships in Fig. 5 verified that the assumption of lognormal distribution for the cascade generators was reasonable.

Parameters of the cascade model estimated with Method I, Modified Method I and Method II

Ponding head (cm)

1.5

567

100 · 100 50 · 50 25 · 25

Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot Plot

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18

Modified Method I ð1Þ

r(2)/r(1)

Method II ð2Þ

ZA

r(1)

ZA

r(1)*

ZA

r(2)

30.4 31.2 23.2 14.3 23.1 18.0 52.3 44.9 47.7 36.9 17.2 14.2 78.9 61.7 62.3 66.7 67.5 80.8

0.153 0.170 0.158 0.141 0.117 0.092 0.084 0.125 0.096 0.116 0.135 0.127 0.080 0.086 0.097 0.087 0.103 0.090

27.1 25.2 19.9 11.6 13.6 15.0 37.3 37.1 36.5 31.7 14.6 13.7 54.1 43.7 30.1 38.2 37.2 45.9

0.163 0.199 0.173 0.166 0.150 0.122 0.835 0.303 0.897 0.291 0.346 0.306 0.127 0.131 0.122 0.117 0.130 0.143

21.1 15.4 17.2 10.6 18.6 15.0 34.9 42.8 43.0 33.3 12.0 12.0 74.2 52.8 42.9 51.5 48.7 47.5

0.801 1.192 0.727 0.531 0.307 0.352 0.955 1.434 0.922 0.801 0.590 0.401 0.415 0.355 0.335 0.329 0.330 0.348

5.25 7.02 4.59 3.76 2.62 3.84 11.4 11.5 9.60 6.91 4.40 3.20 5.18 4.13 3.43 3.77 4.16 3.87

568

K. Wang et al. Dye coverage 0

0.2

0.4

Dye coverage

0.6

0.8

0

1

0.2

0.4

0.6

0.8

1

0

0

Plot 1

Plot 3

Depth (cm)

Depth (cm)

20

40

Observed

60

Observed

40

Method I

80

20

Method I

Modified Method I

Modified Method I

Method II

Method II

60

100

Dye coverage 0

0.2

0.4

Dye coverage

0.6

0.8

1

0

0

0.2

0.4

0.6

0.8

1

0

Plot 7

Plot 6 20

Depth (cm)

Depth (cm)

10

20 Observed

40 60 Observed

80

Method II

Method I

30

Modified Method I

Method I

100

Modified Method I

Method II

120

40

Dye coverage 0

0.2

0.4

0.6

Dye coverage 0.8

1

0

0

40

20

60 Observed

Depth (cm)

Depth (cm)

10

0.6

0.8

1

30 40

Method II

100

0.4

Plot 11

Plot 9 20

80

0.2

0

Method I

Method II

50

Method I Modified Method I

Modified Method I

120

Observed

60

Figure 6 Comparison of measured dye coverage and calculated coverage using random cascade model with parameters estimated with Method I, Modified Method I and those from best fitted parameters using Method II.

and used for comparison. To show lag correlations in the three scales, r obtained from Modified Method I was used in Eq. (12) for the calculation. The mean values of the final and apparent penetration depth distributions were used to calculate the observed e(x, y) correlation. A comparison of observed and cascade estimated correlations of e(x, y) at Plots 8 (100 · 100 cm), 10 (50 · 50 cm), and 11

(25 · 25 cm) is presented in Fig. 7. For all plots, the correlations decreased with lag distance. However, the range at which samples became independent increased with the experimental scale. The calculated correlations by Eq. (12) were higher than the observed ones, implying loss of heterogeneity information when using Modified Method I. The differences between the cascade model estimated

Characterizing heterogeneity of soil water flow by dye infiltration experiments Dye coverage 0

0.2

0.4

0.6

Dye coverage 0.8

1

0

0 20

0.2

0.4

0.6

0.8

1

0

Plot 16

Plot 14

20

40

40

Depth (cm)

Depth (cm)

569

60 80

Observed Method I

100

60 80 Observed

100

Method I

Modified Method I

120

Modified Method I

120

Method II

Method II

140

140

Dye coverage 0

0.2

0.4

0.6

0.8

1

0

Plot 17 20

Depth (cm)

40 60 80

Observed

100

Method I Modified Method I

120

Method II

140

Figure 6 (continued)

and observed correlations were larger in the large plot (Plot 8) than in the smaller plots (Plots 10 and 11).

important for flow transport in soils. The same result was also supported by the above correlation analysis.

Analyses of flow heterogeneity information

Predictions of flow heterogeneity over scales

As discussed above, the standard deviation r can be regarded as an index to evaluate the flow heterogeneity in soils. To characterize the heterogeneity information, we need to consider both horizontal and vertical heterogeneities. The calculation of r(1) was based on the assumption that the zone above the final dye infiltration depth was completely covered by dye, which ignored the heterogeneity along the vertical direction, whereas r(2) was a comprehensive index to characterize both horizontal and vertical heterogeneities and gave more accurate description of flow behavior (Fig. 6). Therefore, we used r(2)/r(1) to compare the variation in the flow domain. For the treatments, the values of r(2)/r(1) ranged approximately from 3 to 12 (Table 3). Compared with Method II, Method I smoothed the flow fields and resulted poorer simulations of heterogeneous flow patterns. Furthermore, the values of r(2)/r(1) increase with the measurement scale, implying that as the scale increases, the multi-dimensional heterogeneities become more

To characterize the flow heterogeneity information at different measurement scales, we plotted values r(2) (average r(2) values of two plots with the same treatment) of the 100 cm plots (100 · 100 cm) versus those of the 50 cm plots (Fig. 8). For each pair of r(2) values, the initial ponding head was the same. Similarly, values of r(2) of the 50 cm plots versus those of the 25 cm plots are also presented in Fig. 8, which shows a linear tread. Based on the linear relationship, we can estimate r(2) values at different measurement scales. For example, we estimated r(2) values of the 100 cm plots from r(2) values of the 50 cm plots or of the ð2Þ 25 cm plots. Although the values of ZA were dependent on both the initial ponding head and the measurement scale, the initial ponding head has a more profound effect ð2Þ on ZA (Table 3). For the same initial ponding head, we used ð2Þ a linear interpolation method to estimate ZA values at difð2Þ ferent measurement scales. For example, average ZA values of the six 25 cm size plots and the six 50 cm size plots were 25.6 and 33.1 cm, respectively. Assuming that average

570

K. Wang et al. 1

1

Plot 10

Plot 8 0.9

Correlation

Correlation

0.9

0.8

0.7

0.8

0.7 Observed

Observed

Simulated

Simulated 0.6

0.6 0

20

40

60

80

100

0

10

20

30

40

50

Lag distance (cm)

Lag distance (cm) 1

Plot 11

Correlation

0.9

0.8

0.7 Observed Simulated

0.6 0

5

10

15

20

25

Lag distance (cm)

(2)

1.4

Dye coverage 2

R = 0.8894

1.2

0

0.2

0.4

0.6

0.8

1

0

1 20

0.8

Depth (cm)

values of 100 and 50 cm plots

Figure 7 Correlation analyses for plots 8, 10, and 11, with the dye source areas of 100 · 100, 50 · 50, and 25 · 25 cm, respectively.

0.6 1.5cm 2.5cm 5.0cm 1.5cm 2.5cm 5.0cm

0.4 0.2 0

0

0.2 (2)

0.4

ponding head(100 vs.50cm) ponding head(100 vs.50cm) ponding head(100 vs.50cm) ponding head(50 vs.25cm) ponding head(50 vs.25cm) ponding head(50 vs.25cm)

0.6

0.8

1

values of 50 and 25 cm plots

Figure 8 Relationship between the standard deviation of the random cascade model and the measurement scale.

40 60 80 100 120

Observed (0.015m ponding head) Predicted (0.015m ponding head) Observed (0.025m ponding head) Predicted (0.025m ponding head) Observed (0.05m ponding head) Predicted (0.05m ponding head)

Figure 9 Prediction of the dye coverage distributions of the 100 cm plots using the random cascade model.

ð2Þ

ZA values of the 50 cm size plots and the 100 cm size plots followed the same ratio (33.1/25.6 = 1.29) and based on the ð2Þ average ZA values of the 50 cm size plots at the different ð2Þ initial ponding heads, we obtained Z A values of 17.95, 49.22, and 60.89 cm for the 100 cm size plots with the initial ponding heads of 1.5, 2.5, and 5 cm, respectively. Using the

random cascade model with the estimated values of r(2) and ð2Þ ZA , we predicted the dye coverage distributions of the 100 cm plots. In general, the predicted results matched with the measured data (the averaged values of the same treatments) quite well (Fig. 9).

Characterizing heterogeneity of soil water flow by dye infiltration experiments

Summary To characterize heterogeneous flow in soils, dye infiltration experiments were carried out in the field and the random cascade model was used to simulate the flow process. The experiments included three measurement scales and three hydraulic boundary conditions (the initial ponding heads at the soil surface). The ratios of the maximum depth and the 50% dye coverage depth decreased with the ponding head for all the measurement scales, which indicated that the hydraulic boundary condition is one of the important factors to affect the patterns of soil water flow. With appropriate methods to estimate the parameters, the random cascade model well described the dye infiltration processes. Based on the experimental data and the random cascade model, different methods were applied to calculate the horizontal and total variations, which characterized the flow heterogeneities in the domain. For the treatments, the ratios of the total flow variation to the horizontal flow variation ranged approximately from 3 to 12. The ratios increased with the measurement scales, indicating that the effect of multi-dimensional heterogeneities on flow transport in soils became more significant as the scale increased. Results from a correlation analysis showed that interdependency between flow paths increased with the experimental scale. A linear relationship was found between the standard deviation of the random cascade model and the measurement scale. By assuming a linear relationship between the mean infiltration depth and the measurement scale, the random cascade model was successfully applied to predict flow heterogeneity over scales.

Acknowledgement This research was financially supported in part by grants of the National Science Foundation of China (Nos. 50579079, 50578127 and 50528910).

References Baveye, P., Boast, C.W., Ogawa, S., Parlange, J.-Y., Steenhuis, T., 1998. Influence of image resolution and threshold on the apparent mass fractal characteristics of preferential flow patterns in field soils. Water Resour. Res. 34, 2783–2796. Boufadel, M.C., Liu, S., Molz, F.J., Lavallee, D., 2000. Multifractal scaling of the intrinsic permeability. Water Resour. Res. 36, 3211–3222. Bouma, J., Dekker, L., 1978. A case study on infiltration into dry clay soil: I. Morphological observations. Geoderma 20, 27–40. Flury, M., Flu ¨hler, H., 1995. Modeling solute leaching in soils by diffusion-limited aggregation: Basic concepts and applications to conservative solutes. Water Resour. Res. 31, 2443–2452.

571

Forrer, I., Parrita, A., Kasteel, R., Fluhler, H., Luca, D., 2000. Quantifying dye tracers in soil profiles by image processing. Eur. J. Soil. Sci. 51, 313–322. Glass, R.J., 1993. Modeling gravity-driven fingering using modified percolation theory. In: Proceedings of the Fourth Annual International Conference on High Level Radioactive Waste Conference, Las Vegas, Nevada. Hillel, D., 1987. Unstable flow in layered soils: A review. Hydrol. Proc. 1, 143–147. Jothityangkoon, C., Sivapalan, M., Viney, N.R., 2000. Tests of a space-time model of daily rainfall in southwestern Australia based on non-homogeneous random cascades. Water Resour. Res. 36, 267–284. Liu, H.H., Molz, F.J., 1997. Multifractal analyses of hydraulic conductivity distributions. Water Resour. Res. 33, 2483–2488. Liu, H.H., Zhang, G., Bodvarsson, G.S., 2003. The active fracture model: Its relation to fractal flow patterns and an evaluation using field observations. Vadose Zone J. 2, 259–269. Molz, F.J., Boman, G., 1993. A fractal-based stochastic interpolation scheme in subsurface hydrology. Water Resour. Res. 29, 3769–3774. ¨ hrstro O ¨m, P., Persson, M., Albergel, J., Zante, P., Nasri, S., Berndtsson, R., Olsson, J., 2002. Field-scale variation of preferential flow as indicated from dye coverage. J. Hydrol. 257, 164–173. ¨ hrstro O ¨m, P., Hamed, Y., Persson, M., Berndtsson, R., 2004. Characterizing unsaturated solute transport by simultaneous use of dye and bromide. J. Hydrol. 289, 23–35. Olsson, J., 1995. Limits and characteristics of the multifractal behaviors of a high resolution rainfall time series. Nonlinear Proc. Geophys. 2, 23–29. Olsson, J., Persson, M., Albergel, J., Berndtsson, R., Zante, P., ¨ hrstro O ¨m, P., Nasri, S., 2002. Multiscaling analysis and random cascade modeling of dye infiltration. Water Resour. Res. 38, 1263. doi:10.1029/2001WR00080. Persson, M., Yasuda, H., Albergel, H., Berndtsson, R., Zante, P., ¨ hrstro Nasri, S., O ¨m, P., 2001. Modeling plot scale dye penetration by a diffusion limited aggregation (DLA) model. J. Hydrol. 250, 98–105. Qiu, G.Y., Momii, K., Yano, T., Lascano, R.J., 1999. Experimental verification of a mechanistic model to partition evapotranspiration into soil water and plant evaporation. Agric. For. Meteorol. 93, 79–93. Schertzer, D., Loverjoy, S., 1987. Physical modeling and analysis of rain and clouds by anisotropic scaling multiplicative processes. J. Geophys. Res. 92, 9693–9714. Schmitt, F.G., 2003. A causal multifractal stochastic equation and its statistical properties. Eur. Phys. J. B 34, 85–98. Schmitt, F.G., Marsan, D., 2001. Stochastic equations generating continuous multiplicative cascades. Eur. Phys. J. B 20, 3–6. Smith, J.E., Zhang, Z.F., 2001. Determining effective interfacial tension and predicting finger spacing for DNAPL penetration into water-saturated porous media. J. Contam. Hydrol. 48, 167–183. Wang, Z., Wu, L., Harter, T., Lu, J., Jury, W.A., 2003. A field study of unstable preferential flow during soil water redistribution. Water Resour. Res. 39, 1075. doi:10.1029/2001WR000903. Yasuda, H., Berndtsson, R., Persson, H., Bahri, A., Takuma, K., 2001. Charactering preferential transport during flood irrigation of heavy clay soil using the dye Vitasyn Blau. Geoderma 100, 49– 66.