WATER RESOURCES RESEARCH, VOL. 34, NO. 5, PAGES 1251–1263, MAY 1998 Stochastic analysis of the relationship between topography and the spatial distribution of soil moisture Pat J.-F. Yeh and Elfatih A. B. Eltahir Ralph M. Parsons Laboratory, Department of Civil and Environmental Engineering Massachusetts Institute of Technology, Cambridge Abstract. This paper deals with the issue of the spatial horizontal variability of soil moisture in the root zone of a shallow soil at the large scale. The problem of water flow in the unsaturated zone is formulated so that topography appears explicitly as a forcing for the movement and redistribution of soil moisture. This formulation emphasizes the role of the lateral redistribution of water that is induced by topography. A stochastic theory is developed to relate the statistical distribution of soil moisture to that of elevation. This approach will ultimately facilitate the use of the readily available data sets describing topography for the purpose of defining the large-scale distribution of soil moisture. The steady state horizontal distribution of soil moisture under homogeneous bare soil conditions is regulated by three distinct factors: topography, climate, and soil properties. First, topography, forces a distribution of soil moisture that tends to mimic the elevation field at large scales. The other two factors are the vertical divergence of water in response to the climate forcing (evaporation) and the capillary resistance to water movement. The climate forcing tends to smooth the spatial distribution of soil moisture. However, the capillary forces exerted by the soil matrix tend to resist displacement of water and hence exert adverse effects against the topography and climate forcings. The variance of the soil moisture distribution increases with the variance of the elevation field and decreases with the correlation scale of the elevation field and the magnitude of the climate forcing. The impact of capillary forces on the vertical fluxes of water is more significant than their impact on the topographically induced horizontal fluxes, owing to the larger hydraulic gradient in the vertical direction resulting from the disparity in scale between the vertical and horizontal directions. 1. Introduction and Motivation The accurate specification of soil moisture distribution is a critical step in two important and related research areas: (1) representation of land and surface processes in climate models [Entekhabi and Eagleson, 1989] and (2) the development of large-scale models of vegetation dynamics [Solomon and Shugart, 1993]. Under certain conditions, soil moisture variability may control land surface processes such as evaporation, runoff, and vegetation growth. The relations between soil moisture and land surface processes are, in general, nonlinear, which suggest that the scaling up of these processes from small scales to large scales has to be carried out with careful consideration to the spatial variability in soil moisture distribution over large areas. Topography, soil type, vegetation, and climate are the key physical factors that control soil moisture distribution over large scales. Although significant information is readily available about topography (for example, from digital elevation maps (DEMs) at a resolution of about 30 m), soil type (from soil survey maps), vegetation (from global land cover data sets), and rainfall (from surface stations and radar networks), only very little information is usually available about the largescale distribution of soil moisture. Some of our ignorance is due to the difficulty in obtaining synoptic measurements of soil Copyright 1998 by the American Geophysical Union. Paper number 98WR00093. 0043-1397/98/98WR-00093$09.00 moisture over large areas using standard equipment such as neutron probes. As a potential solution, microwave remote sensing techniques have been developed for measuring soil moisture from space [Jackson and Schmugge, 1989]. Although these techniques are promising because they offer the potential for coverage of large areas, they would provide wetness information about only the upper few centimeters of the soil profile. The focus of this paper is soil moisture averaged over the total depth of the upper soil layer. The emphasis is on developing a general theory that would ultimately facilitate the use of the significant amounts of data on topography, together with information about rainfall and soil type, for the purpose of defining the large-scale distribution of soil moisture. Although physical intuition suggests that for a given soil type and climate, areas of topographic convergence (divergence) are favorable locations for relatively wet (dry) soil moisture conditions, there is no rigorous theory relating this degree of soil wetness to the corresponding topographic conditions. The primary objective of this paper is to develop a theoretical framework for relating the spatial distribution of soil moisture to topography. We are interested in the spatial distribution of soil moisture (measured in terms of the level of soil saturation relative to the storage capacity of the soil layer), at the following scales: (1) horizontal spatial scale of ⬃10 km and (2) vertical spatial scale (soil depth) of the upper 1 m. The 10-km horizontal scale is considered as a first step before addressing the problem of soil moisture distribution at larger scales (⬃100 km) compatible with the typical resolution of a climate model. 1251 1252 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 1. Methodology. However, the techniques developed in this paper are applicable to the latter problem. The classical treatment of water flow in the unsaturated zone [e.g., Philip, 1957] considers the one-dimensional vertical infiltration into the horizontal surface of an isotropic and homogeneous soil. By definition this classical theory deals with the vertical flow of water. In a more recent study, Giorgini et al. [1984] considered the more general case of two-dimensional flow into sloping surface of anisotropic homogeneous soil. Their theory describes the situations in which lateral flow would be significant. Several recent studies suggested that lateral flow of water is a significant process in the unsaturated zone. McCord and Stephens [1987] observed significant lateral unsaturated flow in a hillslope located in New Mexico. They concluded that [McCord and Stephens, p. 225] “there is strong lateral component to unsaturated flow on a hillslope, even in the absence of apparent sublayers of much lower permeability.” On the basis of theoretical considerations and observations, Zaslavasky and Sinai [1981] indicated that the slope of the soil surface is a significant cause of lateral flow in unsaturated zone. Overall, the results of these studies suggested that the primary causes for lateral flow are the anisotropy of the medium and the gradient of the surface (topography). Hewlett and Hibbert [1963] investigated the flow of water in a sloping soil mass using an experimental approach and found a significant lateral flow component. Many of the observations about lateral flow come from runoff generation studies. Genereux and Hemond [1990] estimated that about 70% of the water flowing into the reach of a small stream in central Massachusetts is flowing from the unsaturated zone. In this paper we emphasize the role of lateral flow in the redistribution of soil moisture under topographic forcing. A simple and popular approach for relating the spatial distribution of soil moisture to topography was suggested by Beven and Kirkby [1979] and is known as the TOPMODEL. The basic idea in this model assumes that the groundwater table intersects with the surface at those locations where the capacity of the saturated soil profile to transport water is smaller than the flux of water. The relation to topography comes from the assumption that the hydraulic gradient for this saturated flow is equivalent to the elevation gradient at the surface. Hence this model is based on considerations of flow in the saturated zone. It neglects lateral flow of water in the unsaturated zone. Using the method of characteristics, Hurley and Pantelis [1985] extended the theory of kinematic subsur- face stormflow proposed by Beven [1981] to include the downslope lateral flow in a porous layer overlying on an impervious bedrock. Their numerical model can be used to calculate the recharge and subsequent drainage following a rainfall event. Using an analytical approach, Stagnitti et al. [1986] proposed a simple hillslope hydrological model for the prediction of the drainage from a sloping shallow soil. The experimental data collected by Hewlett and Hibbert [1963] were used to test their model. The focus of this paper is the root zone of a shallow unsaturated soil where the groundwater table is deep below the surface. The interactions between the root zone and the location of the shallow groundwater table add more complexity to the problem and will not be addressed in this paper. These interactions will be considered as a subject for future research. This paper is organized in six sections. A general methodology for studying the relationship between topography and the distribution of soil moisture is described in section 2. The problem of unsaturated flow of water in the top soil layer is formulated and presented in section 3. A stochastic analysis of the relationship between topography and soil moisture distribution is presented in section 4. The last two sections of the paper include the discussions and conclusions. 2. Methodology This section describes a framework that is proposed for relating the distribution of soil moisture to topography, soil type, and climate. This methodology is summarized in Figure 1. The input is statistical information about topography, soil type, and climate. This information serves as forcing for the equation describing flow of water in the unsaturated zone. The output is statistical characterization of soil moisture distribution at the large scale. This methodology follows a statisticalanalytical approach. Analytical techniques are applied to address the soil moisture distribution problem. The same problem can be addressed by applying numerical techniques using distributed models of catchment processes. These two different approaches are complementary; we chose the analytical approach in order to develop a theory that is general enough to provide some insight regarding the role of topography in determining the spatial distribution of soil moisture. Statistical representations are used to describe the natural heterogeneity of environmental variables. The same approach has been applied successfully to stochastic theories of ground- YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION water hydrology [Gelhar, 1993]. In this paper, similar techniques will be applied to derive the statistical distribution of soil moisture in space. The product of the proposed analysis will be to relate the statistical distribution of soil moisture to the statistical distribution of elevation. This product is adequate for the purposes defined in the introduction of this paper. The horizontal scale considered in this paper is of the order of kilometers, and the vertical scale considered is of the order of meters. This comparison of horizontal and vertical scales suggests that the aspect ratio between the horizontal scale and the vertical scale is of the order of 103. Hence for simplicity we will address the problem in two dimensions, considering the two horizontal coordinates. The treatment of the problem as a two-dimensional horizontal problem brings a significant advantage. In a two-dimensional description the horizontal derivatives of elevation (the contribution of gravity to hydraulic gradient) will no longer be zero, and hence the role of topography becomes explicit. This point will become clear following the formulation of the problem, which will be developed in the next section. The vertical fluxes of water (evaporation and percolation) can be parameterized in terms of soil moisture state as a forcing term in the two-dimensional equation. In the following section the unsaturated flow problem will be formulated following this general methodology. 3. Formulation of the Problem It is widely accepted that the flow of water in the unsaturated zone may be described by q i ⫽ ⫺K共h兲 ⭸h ⭸ xi (1) where q i is specific flow in the ith direction, K is unsaturated conductivity, and x i is distance. The variable h denotes hydraulic head defined by h ⫽ ⫺ ␺ ⫹ z, where ␺ is capillary pressure (a negative sign indicates hydraulic suction) and z is elevation (positive upward). Equation (1) is a nonlinear equation, which complicates treatment of the unsaturated flow problem. Several empirical relations have been suggested in literature to relating K with ␺ and to relate both variables with ␪, defined as available water content per unit total volume. Following Gardner [1958] and Mantoglou and Gelhar [1987], we assume that K ⫽ K se ⫺␣␺ (2) ␪ ⫽ ␪s ⫺ C␺ where K s is the saturated conductivity, ␪ s is the saturated soil moisture content, and C is the specific moisture capacity, which is the slope of the soil-water retention curve (i.e., the rate of change of soil moisture content with respect to capillary pressure head). Here ␣ represents the relative rate of decrease of hydraulic conductivity with increasing capillary pressure head and is associated with the width of soil pore size distribution, and ␣⫺1 is the thickness of capillary fringe, a measure of the relative importance of capillary force to gravity force for soil moisture movement in a specific soil. Fine-textured soils in which capillary force tends to dominate have greater thickness of capillary fringe than coarse-textured soils, in which gravity effects manifest themselves most readily [Philip, 1969]. According to published literature (for a review, see Pullan [1990]), the range of ␣⫺1 covers 0.01–10 m. However, 0.2–5 m seems to be the typical range of values for ␣⫺1 [Philip, 1969]. 1253 The substitution of (2) into (1) results in a significant simplification qi ⫽ ⫺ ⭸z 1 ⭸K ⫺K ␣ ⭸ xi ⭸ xi (3) Equation (3) is a linear equation on K which represents a significant advantage when compared with (1). Topography is represented by the elevation gradient, which appears explicitly in (3). (Note that since the problem is considered in two dimensions, the partial derivatives of z are no longer zero.) Conservation of water mass requires D ⭸␪ ⭸q i ⫺s⫹R ⫽ ⫺D ⭸t ⭸ xi (4) where D is the depth of root zone (assumed to be constant) and s is the sink of water mass due to vertical fluxes. The variable R is effective rainfall, defined as the fraction of rainfall that infiltrates into the soil. Although the value of this fraction depends on the soil moisture conditions at the surface, for simplicity we assume that effective rainfall is a constant fraction of surface rainfall. In most natural settings the top layer of the soil (1–2 m) has significantly larger hydraulic conductivity than the underlying layers [Beven, 1982]. Here we assume that the underlying layer is completely impermeable. Hence the sink of moisture due to vertical divergence, s, represents the evaporation flux upward to the atmosphere. Vertical fluxes of transpiration and percolation are not considered in this study. This evaporation flux can be parameterized as s ⫽ ␤ K, assuming that the evaporation flux is proportional to the unsaturated hydraulic conductivity of the soil. Note that physically, ␤ is equivalent to the vertical hydraulic gradient near the ground surface. From the Gardner [1958] relationship in (2), the evaporation flux s for ␪ ⱕ ␪ s , can be written as s ⫽ ␤ K ⫽ ␤ K se 关␣共␪⫺␪s兲兴/C (5) which increases exponentially to approach its maximum of ␤ K s as ␪ increases from the residual moisture content to the saturation. Note that the exponential dependence in (5) resembles the widely recognized pattern of the observed relationships between evaporation and soil moisture level of a bare soil [see Lowry, 1959, Figure 1; Rodriguez-Iturbe et al., 1991, Figure 2]. This fact provides an observational basis for the proposed parameterization as a rough description of evaporation flux. For a specific soil type the magnitude of the gradient ␤ depends on climate and land cover. However, this parameterization may not be appropriate for a deep soil layer where variability in the vertical soil moisture distribution is significant. For deep soils the rate of evaporation decays with depth, and the vertical flow induced by gravity becomes significant. This vertical flow exerts additional influences on both the vertical climate flux and horizontal soil moisture redistribution, which needs the consideration of a three-dimensional description. Therefore a shallow, unsaturated soil overlying an impervious soil layer or bedrock is the condition considered in this paper. The more complicated treatment of three-dimensional flow in deep soils is left as a topic for future research. Typical values of ␤ can be estimated from typical values of evaporation. The range of evaporation is assumed to be 50 – 100 cm yr⫺1; typical values for saturated hydraulic conductivity K s are 8 ⫻ 10⫺3 cm s⫺1 for sand and 3 ⫻ 10⫺5 cm s⫺1 for clay 1254 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 2. (a) Elevation distribution and (b) correlation function for a 3.2-km transect (5-m resolution) from the Mahantango Creek drainage basin, Pennsylvania. The correlation function is fitted by two different theoretical correlation functions: exponential function in (22) and hole-type function in (23). The dotted line in the figure marks the correlation scale corresponding to the correlation of e ⫺1 . [Bras, 1990, p. 352]. Therefore the range of ␤, which is an indicator of the strength of evaporation, is approximated as 10⫺4–10⫺1 in this paper. Substituting for q and s in (4) results in a general formulation of the problem D ⭸␪ ⭸ ⫽D ⭸t ⭸ xi 冉 冊 1 ⭸K ⭸z ⫺ ␤K ⫹ R ⫹K ␣ ⭸ xi ⭸ xi (6) The basic variable in this formulation is the unsaturated hydraulic conductivity K, which serves as a surrogate for soil moisture content. Even if soil properties are constant in space, K would vary in space depending on soil moisture variability. All the variables and parameters considered in this study are effective variables and parameters representing a block of soil with a vertical dimension of about 1 m. The root zone in a natural environment is composed of heterogeneous soil, roots, and possibly macropores. In this study it will be assumed for simplicity that classical theories describing flow of water in the unsaturated zone (such as (1) and (2)) are valid for characterizing the behavior of the root zone. Similar assumptions are made in most field applications that deal with groundwater problems. In the following, the steady state case will be considered. It is reasonable to assume steady state when considering long timescales. Equation (4) can be averaged in time to obtain D ⭸E关 ␪ 兴 ⭸E关q i兴 ⫽ ⫺D ⫺ E关s兴 ⫹ E关R兴 ⭸t ⭸ xi (7) where E[ ] denotes the expected value or time average. Since large-scale soil moisture distribution is of interest, and since evaporation is the only process responsible for the vertical divergence of soil moisture (i.e. monotonic drying), the hysteresis effect is disregarded in this paper. At steady state, the left-hand side of (7) is negligible. Substituting for q and s into (7) results in D ⭸ ⭸ xi 冉 冊 1 ⭸E关K兴 ⭸z ⫹ E关K兴 ⫺ ␤ E关K兴 ⫹ E关R兴 ⫽ 0 ␣ ⭸ xi ⭸ xi (8) Equation (8) is an ordinary differential equation on E[K]. Time-averaged effective rainfall, E[R], is almost uniform over the scale of interest (i.e., order of kilometers). Hence the main forcing in (8) is spatial variability in topography, which is represented by the gradients of elevation. This equation relates unsaturated hydraulic conductivity K, which is a surrogate for soil moisture content, to topography, described by elevation z. These two variables will be considered as random processes in the stochastic analysis of the relationship between soil moisture and topography presented in the next section. 4. Stochastic Analysis of the Relationship Between Topography and Soil Moisture Observations of elevation fields indicate significant variability for a range of spatial scales. S. Lancaster (unpublished report, 1993) analyzed several transects of elevation extracting from DEMs for different regions of the continental North America. The resolution of DEMs is 30 m, and the length of these transects is about 30 km. After appropriate detrending, the correlation scales of these elevation distributions are found to lie within the range from several hundred meters to several kilometers (i.e., 102–104 m). Note that the estimated correlation scale is closely dependent on the resolution of data: finerresolution data have the impact of shortening the correlation scale estimated from a coarser one. We also have analyzed the statistical distribution of elevation which provides an adequate representation of observed topography. Figure 2 shows one of the transects of elevation distribution and the corresponding autocorrelation function from the Mahantango Creek drainage basin, Pennsylvania. The resolution of the pixels is 5 m. The standard deviation and correlation scale of this elevation transect are roughly 20 m and 500 m, respectively. Then the correlation function is fitted by two frequently used theoretical models of correlation [Gelhar, 1993, p. 44]: exponential and hole-type function (which are given later in (22) and (23)). Following Bakr et al. [1978], the correlation scale of hole-type correlation function is taken as 2.5 times of that of the exponential correlation function for a consistent fit of a given data. The assumption that the elevation field, which appears in (8), is a random field transforms the equation into a stochastic differential equation. In this paper we investigate the solutions of this stochastic differential equation when forced by statistical descriptions of the elevation fields. The mathematical technique that will be used in solving this equation is the perturbation method. Application of this method results in equations relating small perturbations in unsaturated hydraulic conductivity to perturbations in elevation field. The solutions of such equations will be explored using spectral representations of stationary random fields. The final result will be to relate the statistical distribution of unsaturated conductivity to the statis- YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION tical distribution of elevation. Following this step, the statistical distribution of ␪ will be estimated from the statistical distribution of K through (2). 4.1. ⭸ ⭸ xi 冉 䡠 冊 冋 ៮ ⭸ z៮ 1 ⭸K ⭸K⬘ ⭸ z⬘ ៮ ៮ ⫹ R ⫹ DE ⫹K ⫺ ␤K ␣ ⭸ xi ⭸ xi ⭸ xi ⭸ xi 册 ៮ ⫽ R/ ␤ K That is, the large-scale, steady state average unsaturated hydraulic conductivity is proportional to effective rainfall and inversely proportional to vertical hydraulic gradient (i.e., evaporation). Subtracting the mean equation (9) from (8), the perturbation equation is given by 冉 冊 冕 冕 K⬘ ⫽ In order to relate soil moisture content (hereinafter referred as SM for brevity) and topography, we analyze the relationship between SM and K. Combining Gardner’s relationships from (2), we get ␪ ⫽ ␪s ⫹ (12) ⬁ e i共kxx⫹kyy兲 dZ K共k x, k y兲 ␪⬘ ⫽ S K共k x, k y兲 ␣ 2K៮ 2共k 2x ⫹ k 2y 兲 2 S 共k , k 兲 共k ⫹ k ⫹ ␣␤ /D兲 2 ⫹ ␣ 2关k x共⭸ z៮ /⭸ x兲 ⫹ k y共⭸ z៮ /⭸ y兲兴 2 z x y 2 x 2 y (13) Therefore the variance of K, ␴ 2K , is related to that of z, ␴ 2z , by the following relation: ␴ 2K ⫽ 冕冕 ⫺⬁3⬁ where ៮ 2␻ 20␴ 2z S K共k x, k y兲 dk x dk y ⫽ ␣ 2K C ␣ K៮ K⬘ (14) 兩 ␨ 兩⬍⬍1 (16) (17) Using the spectral representation of ␪⬘ 冕 ⬁ e i共kxx⫹kyy兲 dZ ␪共k x, k y兲 (18) ⫺⬁ (17) leads to the following relations between the spectrum and variance of SM, S ␪ (k x , k y ) and ␴ 2␪ , and those of the unsaturated hydraulic conductivity: S ␪共k x, k y兲 ⫽ C2 S K共k x, k y兲 ␣ 2K៮ 2 C2 ␴ ␪2 ⫽ 2 ៮ 2 ␴ 2K ␣K (19) The goal is to relate the variance in SM to the variance of elevation. This is achieved by combining (14) and (19) to obtain ␴ ␪2 ⫽ C 2␻ 20␴ 2z where k x , k y are wave numbers in x and y directions, and dZ is the spectral amplitude of a random process. By substituting the spectral relationships (12) into (11), we find that the spectral density function of K, S K (k x , k y ), is related to that of z, S z (k x , k y ), through (15) is used with its first-order term only to deduce that ⫺⬁ ⫽ ␨2 ␨3 ⫹ ⫺ ⫹··· 2 3 ln 共1 ⫹ ␨ 兲 ⫽ ␨ ⫺ (11) e i共kxx⫹kyy兲 dZ z共k x, k y兲 C K ln ␣ Ks ៮ ⫹ K⬘, and K⬘ ⬍⬍ K ៮ , the Assuming that ␪ ⫽ ␪៮ ⫹ ␪ ⬘, K ⫽ K Maclaurin series ⬁ ⫺⬁ 2 S z共k x, k y兲 dk x dk y ␴ 2z ␪⬘ ⫽ In the derivation of (11) it has been assumed that the second-order term, (⭸K⬘/⭸ x i )/(⭸ z⬘/⭸ x i ), which includes the product of perturbation terms, is equal to its expected value in brackets in the mean equation (9). Equation (11) will be solved using spectral techniques by defining the spectral representations of the two dimensional stationary random fields K⬘ and z⬘ as follows [Gelhar, 1993]: z⬘ ⫽ 2 y (9) (10) ⭸ z⬘ ⭸ z៮ 1 ⭸K⬘ ៮ ⫺ ␤ K⬘ ⬇ 0 ⫹K ⫹ K⬘ ␣ ⭸ xi ⭸ xi ⭸ xi 共k 2x ⫹ k 2y 兲 2 共k ⫹ k ⫹ ␣␤ /D兲 ⫹ ␣ 2关k x共⭸ z៮ /⭸ x兲 ⫹ k y共⭸ z៮ /⭸ y兲兴 2 2 x ⫽0 i ⫽ 1, 2 ៮ ៮ can be Assuming constant large-scale averages K and z៮ , K derived if the second-order product term in the bracket of (9) is neglected: ⭸ D ⭸ xi 冕冕 ⫺⬁3⬁ Preliminary Analysis We first write (8) in two dimensions and drop the implied expected value notation for simplicity. Each of K and z is decomposed into two terms, a spatial mean and a perturbation ៮ , K⬘, z៮ , z⬘). Then by taking the average of (8), we term (K obtain the following partial differential equation: D ␻ 20 ⫽ 1255 (20) where ␻20 is defined in (14). 4.2. One-Dimensional Analysis The one-dimensional unsaturated flow in a varying elevation field is considered. First, by using the one-dimensional version of the spectral relationships, (13) and (19), and assume that the trend of elevation does not exist, we derive S ␪共k兲 ⫽ C 2k 4 S 共k兲 共k ⫹ ␣␤ /D兲 2 z 2 (21) In order to evaluate the variance and covariance function of SM, the exponential covariance and hole-type covariance functions plotted in Figure 2b are selected to describe the elevation field. They will be compared to investigate the impact of the functional form of input statistics on the analytical results. The covariance-spectrum pair for the exponential covariance function is 1256 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 3. One- and two-dimensional normalized SM variances ␴ 2␪ / ␴ 2z C 2 in (25), (27), and (34), resulting from the input of one-dimensional exponential and hole-type covariance functions of elevation and from the two-dimensional Bessel-type covariance function of elevation, respectively. (The correlation scale ␩ of holetype covariance function is taken as 2.5 times of that of the exponential covariance function, according to Bakr et al. [1978]). R z共 ␰ 兲 ⫽ ␴ 2z exp 冉 冊 ⫺ in (23) can be substituted into (21), and the covariance function R ␪ ( ␰ ) is evaluated as (see Appendix A) 兩␰兩 ␭ (22) 2 ␴ 2z ␭ S z共k兲 ⫽ ␲ 共1 ⫹ ␭ 2k 2兲 R ␪共 ␰ ␩兲 ⫽ ␴ 2z C 2 and the covariance-spectrum pair for the hole-type covariance function is R z共 ␰ 兲 ⫽ ␴ 2z 共1 ⫺ 兩 ␰ 兩/ ␩ 兲 exp 共⫺兩 ␰ 兩/ ␩ 兲 2 z (23) 3 2 2␴ ␩ k S z共k兲 ⫽ ␲ 共1 ⫹ k 2␩ 2兲 2 ⫹ 冋 1/ 2 ␭ g exp 共⫺g ␭1/ 2␰ ␭兲共1 ⫹ g ␭1/ 2␰ ␭兲 2共1 ⫺ g ␭兲 g ␭1/ 2共 g ␭ ⫺ 2兲 1 exp 共⫺g 1/ 2␰ ␭兲 ⫹ exp 共⫺␰ ␭兲 共 g ␭ ⫺ 1兲 2 共 g ␭ ⫺ 1兲 2 g ␭ ⫽ ␣␤␭ 2/D 册 (24) ␰ ␭ ⫽ 兩 ␰ 兩/ ␭ where g ␭ is a dimensionless parameter and ␰␭ is the lag normalized by the correlation scale ␭. The variance of SM is obtained when ␰␭ ⫽ 0 ␴ ␪2 ⫽ 共 g ␭3/ 2 ⫺ 3g ␭1/ 2 ⫹ 2兲 2 2 ␴zC 2共 g ␭ ⫺ 1兲 2 ⫺g ␩3/ 2 exp 共⫺g ␩1/ 2␰ ␩兲共1 ⫹ g 1/␩ 2␰ ␩兲 共 g ␩ ⫺ 1兲 2 2g ␩3/ 2共 g ␩ ⫺ 3兲 1 exp 共⫺g 1/␩ 2␰ ␩兲 ⫺ 共 g ␩ ⫺ 1兲 3 共 g ␩ ⫺ 1兲 2 䡠 exp 共⫺␰ ␩兲共1 ⫹ ␰ ␩兲 ⫹ 2共3g ␩ ⫺ 1兲 exp 共⫺␰ ␩兲 共 g ␩ ⫺ 1兲 3 g ␩ ⫽ ␣␤␩ 2/D where ␭ and ␩ are the correlation scale of the exponential and the hole covariance functions, respectively. By substituting (22) in (21), the covariance function of SM distribution is derived in terms of the inverse Fourier transform of S␪(k) (see Appendix A): R ␪共 ␰ ␭兲 ⫽ ␴ 2z C 2 ⫹ 冋 (25) Using the same approach, the hole covariance-spectrum pair 册 (26) ␰ ␩ ⫽ 兩 ␰ 兩/ ␩ The variance of ␪ is obtained if ␰␩ ⫽ 0 ␴ ␪2 ⫽ 共 g 5/␩ 2 ⫺ 5g 3/␩ 2 ⫹ 5g ␩ ⫺ 1兲 2 2 ␴zC 共 g ␩ ⫺ 1兲 3 (27) The SM variances (normalized by ␴ 2z C 2 ) in (25) and (27) corresponding to the input of exponential and hole-type covariance functions of elevation, are graphed in Figure 3. Figure 4 shows the sensitivity of SM variance to the parameters ␭, ␤, ␣⫺1, and D. Moreover, the correlation function of SM (i.e., covariance function R ␪ in (24) and (26) divided by the variance in (25) and (27), respectively) corresponding to the exponential and hole-type function inputs are graphed in Figure 5 for several values of ␭. Figure 6 shows the sensitivity of the SM correlation function to ␤, ␣⫺1, and D (Figures 6a, 6b, and 6c, respectively), as well as the correlation scales of SM distribution for several typical values of ␤, ␣⫺1, and D. These results will be discussed in section 5. In order to investigate the degree of correlation between YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION 1257 topography and SM distribution, the cross-covariance function between ␪ and elevation z is derived as follows: Exponential function R ␪z共 ␰ ␭兲 ⫽ ␴ 2z C 冋 g ␭1/ 2 1 exp 共⫺g ␭1/ 2␰ ␭兲 ⫺ exp 共⫺␰ ␭兲 g␭ ⫺ 1 g␭ ⫺ 1 册 (28) Hole-type function R ␪z共 ␰ ␩兲 ⫽ ␴ 2z C 冋 2g ␩3/ 2 1 exp 共⫺g 1/␩ 2␰ ␩兲 ⫹ 共 g ␩ ⫺ 1兲 2 g␩ ⫺ 1 䡠 exp 共⫺g ␩1/ 2␰ ␩兲共1 ⫹ g 1/␩ 2␰ ␩兲 ⫺ 2共2g ␩ ⫺ 1兲 exp 共⫺␰ ␩兲 共 g ␩ ⫺ 1兲 2 册 (29) The derivations of (28) and (29) are similar to those presented in Appendix A and hence are not included in this paper. The one-point (i.e., lag zero) cross covariance between ␪ and z can be derived if ␰␭ ⫽ ␰␩ ⫽ 0 in (28) and (2): Exponential function ␴ ␪z ⫽ ␴ 2z C 1 ⫹ g ␭1/ 2 (30) Hole-type function ␴ ␪z ⫽ ␴ 2z C 2g ␩3/ 2 ⫺ 3g ␩ ⫹ 1 共 g ␩ ⫺ 1兲 2 (31) The one-point (lag zero) cross correlations derived from dividing (30) and (31) by the corresponding ␴ ␪ ␴ z , are displayed in Figure 7, and the cross-correlation functions (i.e., cross- Figure 5. (a) One-dimensional and (b) two-dimensional correlation functions of SM distribution, R ␪ / ␴ 2␪ , for several values of elevation correlation scale ␭ (␤ ⫽ 10⫺3, ␣⫺1 ⫽ 0.5 m, D ⫽ 1 m). The dotted line in the figure marks the correlation scale corresponding to the correlation of e ⫺1 . covariance functions in (28) and (29) divided by the corresponding ␴ ␪ ␴ z , respectively) resulting from the input of exponential and hole-type functions are shown in Figure 8. These figures are discussed in section 5. 4.3. Two-Dimensional Analysis The spectral approach can be extended to the twodimensional unsaturated flow case by using multidimensional spectral representation of elevation and soil moisture. A Bessel-type spectrum-covariance pair for two-dimensional isotropic spatial processes was suggested by P. Whittle [see Gelhar, 1993, p. 47] R z共 ␰ 兲 ⫽ ␴ 2z 冉冊 ␰ ␰ K ␭ 1 ␭ (32) ␴ 2z ␭ 2 S z共k兲 ⫽ ␲ 共1 ⫹ ␭ 2k 2兲 2 where k 2 ⫽ k 2x ⫹ k 2y . K 1 is the modified Bessel function of first order. The covariance function of SM can be obtained by introducing (32) into (21) (see Appendix B): R ␪共 ␰ ␭兲 ⫽ C 2␴ 2z 冋 䡠 K 0共 g ␭1/ 2␰ ␭兲 ⫹ 1 4g ␭ K 共␰ 兲 ␰⬘K 共␰ 兲 ⫺ 共 g ␭ ⫺ 1兲 2 ␭ 1 ␭ 共 g ␭ ⫺ 1兲 3 0 ␭ g ␭3/ 2 4g ␭ ␰ K 共 g 1/ 2␰ 兲 ⫹ 共 g ␭ ⫺ 1兲 2 ␭ 1 ␭ ␭ 共 g ␭ ⫺ 1兲 3 册 (33) where the dimensionless parameters g ␭ and ␰␭ are defined in (24). The corresponding SM variance is (see Appendix B) ␴ ␪2 ⫽ Figure 4. Sensitivity of (a) elevation correlation scale ␭, (b) vertical hydraulic gradient ␤, (c) capillary fringe thickness ␣⫺1, and (d) soil depth D on the variance of SM distribution for typical ranges of these parameters: ␭ ⫽ 103 m, ␣⫺1 ⫽ 0.5 m, ␤ ⫽ 10⫺3, D ⫽ 1 m, ␴ z ⫽ 30 m. The variances of SM distribution for three different values of the standard deviation of elevation are shown in Figure 4a. 冋 册 g␭ ⫹ 1 2g ␭ ln g ␭ ␴ 2z C 2 2⫺ 共 g ␭ ⫺ 1兲 共 g ␭ ⫺ 1兲 3 (34) The two-dimensional SM variance in (34) is displayed in Figure 3 as a function of dimensionless parameter g ␭ along with the SM one-dimensional variances. The correlation function of two-dimensional SM distribution in (33) is shown in Figure 5b for several values of ␭. The cross-covariance function between ␪ and z can be derived as follows: R ␪z共 ␰ ␭兲 ⫽ C ␴ 2z ⫹ 冋 ⫺2g ␭ 1 K 共 g 1/ 2␰ 兲 ⫺ ␰⬘K 共␰ 兲 共 g ␭ ⫺ 1兲 2 0 ␭ ␭ 共 g ␭ ⫺ 1兲 ␭ 1 ␭ 2g ␭ K 共␰ 兲 共 g ␭ ⫺ 1兲 2 0 ␭ 册 (35) 1258 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 6. One-dimensional correlation functions of SM for several values of (a) vertical hydraulic gradient ␤, (b) capillary fringe thickness ␣⫺1, and (c) soil depth D, and the corresponding correlation scales (␭ ⫽ 103 m). The corresponding one-point cross covariance (i.e., at lag zero) is ␴ ␪z ⫽ ␴ 2z C 冋 2g ␭ 1 ln g ␭ ⫺ 共 g ␭ ⫺ 1兲 2 g␭ ⫺ 1 册 (36) The one-point cross correlation ((36) divided by ␴ ␪ ␴ z ) is plotted in Figure 7 for comparison with the corresponding one-dimensional solutions. The cross-correlation function (i.e., (35) divided by ␴ ␪ ␴ z ) is shown in Figure 8 along with the corresponding one-dimensional solutions. 5. Discussion Three distinct mechanisms can be identified as the dominant factors in determining the steady state horizontal distribution of soil moisture assuming a bare homogeneous soil condition. First, topography forces a large-scale distribution of SM that mimics the elevation field. Specifically, the amplitude and correlation scale of SM distribution are proportional to those of elevation distribution if only topographic forcing is considered. Figure 7. One- and two-dimensional one-point (i.e., lag zero) cross-correlations between SM ␪ and elevation z as functions of dimensionless parameter g ␭ . YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 8. (a) One-dimensional and (b) two-dimensional cross-correlation functions between SM ␪ and elevation z for several values of elevation correlation scale ␭ (␤ ⫽ 10⫺3, ␣⫺1 ⫽ 0.5 m, D ⫽ 1 m). Generally, areas of topographic convergence (valleys) are favorable locations for the accumulation of SM, while areas of topographic divergence (ridges) usually experience dry SM conditions. Second, the climate forcing (evaporation) smoothes the SM distribution by differentially extracting water from the relatively wet regions. This results in a reduced variance of SM distribution. Finally, the soil properties affect the SM distribution through capillary forces, which can be quantified by the thickness of the capillary fringe, ␣⫺1. For a thick capillary fringe (␣⫺1 large, for example, clayey soils), water molecules are strongly bound to the soil and resist gravity drainage, mainly as a result of the smaller pores and hence larger surface area. Relatively, sand and gravel contain a high percentage of larger pores such that water is only weakly bound to the soil by capillary forces. Therefore SM in sandy soils with smaller capillary resistance can be removed by the topographic forcing or climatic forcing relatively easier than is the case in clayey soils. Since the vertical divergence fluxes (i.e., evaporation) have the effect of smoothing the SM distribution, this capillary resistance has an adverse effect on the vertical fluxes, and hence would tend to enhance the variance of the SM distribution. Contrarily, in the horizontal direction, the variance of SM distribution induced by topography should be reduced because of the capillary resistance. However, as will be illustrated later, the impact of capillary resistance on the vertical fluxes of water is more significant than its impact on the topography-induced horizontal fluxes, which leads to a net increase of SM variance as the capillary resistance increases. This is due to the larger hydraulic gradient in the vertical direction resulting from the disparity in scale between the vertical and horizontal directions. In the following we discuss the results presented in the Figures 3–9. In Figure 3 the normalized variance of SM is plotted as a function of dimensionless parameter g ␭ , which is 1259 defined by g ␭ ⫽ ␣␤␭ 2 /D. (Unless stated otherwise, typical values such as correlation scale of elevation ␭ ⫽ 103 m, vertical hydraulic gradient ␤ ⫽ 10⫺3, capillary fringe thickness ␣⫺1 ⫽ 0.5 m, soil depth D ⫽ 1 m, specific moisture capacity C ⫽ 0.1 m⫺1 and variance of elevation ␴ 2z ⫽ (30 m)2 are taken as the nominal values in producing these figures). Figure 3 shows that the variance of SM distribution decreases as g ␭ increases. Large g ␭ results from large ␭ and ␤, as well as small ␣⫺1 and D. The value of g ␭ ranges over several orders of magnitude mostly due to the wide range of ␭ (⬵ 10 2 –10 4 m) and ␤(⬵ 10⫺4–10⫺1). From Figure 3 it can be seen that using the exponential and the hole-type covariance function of elevation virtually leads to a similar magnitude of SM variance and a similar rate of decrease with g ␭ . Therefore the statistical properties of SM distribution are not sensitive to the functional form of elevation covariance. The two-dimensional SM variance resulting from the Bessel-type covariance is approximately one order of magnitude smaller than the onedimensional SM variance at small g ␭ , and 3– 4 orders of magnitude smaller at large g ␭ . That is, two-dimensional SM variance decreases with increasing g ␭ at a rate faster than the one-dimensional variance. This stems from the increased freedom of SM movement in two dimensions than in one dimension. In one dimension a low-conductivity zone (i.e., dry soil) can essentially prohibit SM movement along its gradient direction. However, in two dimensions, soil water movement can take place laterally, passing around the low-conductivity region and consequently causing a smoother distribution of SM. Similar conclusions have been reached in the published literature on stochastic groundwater hydrology. [e.g., Bakr et al., 1978; Gelhar, 1993]. The sensitivity of SM variance to ␭, ␣⫺1, ␤, and D is shown in Figure 4. Importantly, the correlation scale of topography has a direct impact on SM distribution. As shown in Figure 4a, a smooth elevation field with large correlation scale leads to a smooth SM distribution, while a surface with significant topographic undulation leads to a highly variable SM distribution. The variance of SM distribution is proportional to the variance of elevation distribution (see the expressions in (25) and (27)). This trend remains unaltered although not presented in Figures 4b, 4c, and 4d. In addition, a larger variance of SM distribution arises from a smaller vertical divergence of SM (small ␤), a fine-textured soil (large capillary resistance), and a deep soil layer. (The elevation variance of (30 m)2 is used in producing Figure 4). A large vertical gradient ␤ creates a substantial vertical moisture sink that eventually reduces the variance of horizontal SM distribution. Large capillary fringe thickness ␣⫺1 constrains the moisture movement in both horizontal and vertical direction. Because of the large aspect ration between the horizontal and vertical scales (i.e., ⬃103) considered, a large disparity in the hydraulic gradient exists in these two directions. The effect of ␣⫺1 on limiting the efficiency of vertical climatic fluxes and hence smoothing the SM distribution is thus more significant than the effects of the capillary resistance on the horizontal fluxes induced by topographic forcing. This explains the slowly rising trend of SM variance with the capillary fringe thickness in Figure 4b. Further, soil depth, which can be viewed as a characteristic scale of the soil reservoir, regulates the magnitude of vertical climatic (evaporative) fluxes in terms of hydraulic gradient in the vertical direction. A relatively deep soil layer alleviates the rate of vertical fluxes and as a result provides a smaller reduction of the variance of SM horizontal distribution. 1260 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION Figure 5 shows the one- and two-dimensional correlation functions of SM distribution for several values of elevation correlation scale ␭ with typical value ␤ ⫽ 10⫺3. The dotted line marks the correlation scale corresponding to the correlation of e ⫺1 . In general, smaller ␭ leads to a SM horizontal distribution with a shorter correlation scale, and the hole phenomenon occurs at small lags (i.e., 10 –30 m) in the one-dimensional correlation function. The one-dimensional SM correlation functions resulting from the exponential and hole-type covariance functions of elevation are similar. This again suggests that the statistical properties of SM distribution are not sensitive to the functional form of the input covariance. The overlap of the SM correlation for ␭ ⫽ 103 m and ␭ ⫽ 104 m indicates that the correlation structure of topography has little impact on the SM distribution for ␭ larger than 103 m (for the typical value ␤ ⫽ 10⫺3). For two-dimensional analysis the overlap of the SM correlation function for large ␭ is similar to that in one dimension. However, for the range of ␭ (103–104 m), the twodimensional SM distribution tends to correlate over longer distances (i.e., 100 –200 m) in comparison with the onedimensional distribution (i.e., 20 –30 m). Further, the hole phenomenon occurring in one dimension at small lags does not occur in two dimensions for typical values of ␭. In two dimensions, soil water has greater freedom in choosing a flow path and usually avoids the relatively low conductivity zones, which leads to a smoother SM horizontal distribution. The influences of vertical hydraulic gradient ␤, capillary resistance ␣⫺1, and soil depth D on the one-dimensional SM correlation function and its correlation scale are illustrated in Figure 6. The vertical climatic fluxes have the effect of shortening the correlation scale of SM horizontal distribution. For the homogeneous soil system investigated in this paper, if the elevation field is flat, the correlation scale of SM horizontal distribution would be infinitely large and the variance of SM distribution would be zero. Because of the absence of horizontal unsaturated flow, this situation remains unchanged regardless of the magnitude of ␤. For a rolling topography, SM tends to accumulate at the valleys, which results in a nonuniform horizontal distribution of SM. Under this situation, water is removed by vertical climatic fluxes from the high-SM regions, which leads to a smaller amplitude and a smaller correlation scale of horizontal distribution. This is demonstrated in Figure 6a such that the increase of ␤ results in the decrease of SM correlation scale. As was stated in the preceding section, ␤ is the vertical hydraulic gradient driving water away from a high-SM region to the atmosphere by evaporation. This process balances the SM differences between wet and dry regions in soils and thus has the effect of reducing the variance and correlation scale of SM distribution. On the other hand, the capillary force of soil matrix tends to resist both the horizontal topographic forcing and vertical climatic forcing; that is, it has a negative effect on the tendency of climate to reduce the correlation scale. This impact of capillary resistance is shown in the rising trend of SM correlation scale in Figure 6b. Further, the soil depth D regulates the magnitude of vertical fluxes. Large D (resulting in weaker vertical hydraulic gradient) tends to dampen the effect of ␤, whereas small D tends to enhance that effect (Figure 6c). The one- and two-dimensional one-point (i.e., lag zero) cross correlations ␴ ␪ z are plotted in Figure 7 as a function of dimensionless parameter g ␭ . As seen in the figure, the onepoint cross correlation between elevation and SM distribution decreases as g ␭ increases. This trend is similar to that of SM variance demonstrated in Figure 3. Large ␭ and ␤, as well as small ␣⫺1 and D, reduce the cross correlation of topography and SM. The one-point cross correlation approaches zero as ␭ increases to a very large value, which indicates the independence of topography and SM distribution for the situation of a uniform elevation field. For small g ␭ (small ␭ and ␤), Figure 7 implies a relatively high correlation between ␪ and z. Note that a perfect correlation relation between ␪ and z corresponds to ␤ ⫽ 0 (i.e., unit cross correlation can be derived from substituting ␤ ⫽ 0 in (25) and (30)), which describes the situation when topography is the only dominant forcing. Furthermore, for small g ␭ the one-point cross correlation in two dimensions is higher than that in one dimension. This regime is dominated by the topographically induced horizontal fluxes. Under such conditions, water has more freedom to move in two dimensions and hence tends to closely mimic topography, resulting in a higher cross correlation. As g ␭ increases (␭ or ␤ increases), the two-dimensional one-point cross correlation between ␪ and z approaches the one-dimensional cross correlation as the role of topographically induced horizontal fluxes becomes less and less important in comparison to the climate forced vertical fluxes. Figure 8 shows the one- and two-dimensional crosscorrelation function between ␪ and z for several values of ␭ with ␤ ⫽ 10⫺3. Generally, the trends of one- and twodimensional cross correlations between ␪ and z are similar except that in two dimensions, SM distribution correlates with the topography for a larger distance (i.e., hundreds of meters) than it does in one dimension (i.e., tens of meters). For a small ␭ the cross correlation between ␪ and z is relatively larger for small lags and then decays faster than the corresponding cross correlation for a larger ␭. Overall, the degree of cross correlation between ␪ and z decreases with the ␭. As was mentioned earlier, water in a hilly topography tends to accumulate at the topographic convergence areas (bottom of the valleys) which usually lie next to areas of topographic divergence (ridges). This implies the negative cross correlation between ␪ and z for the separation distance at the order of ␭. Since ␭ is representative of the average scale of topographic undulation, a small ␭ in Figure 8 describing a rolling topography suggests that the SM distribution has higher positive correlation with the elevation at lags smaller than ␭, while the same distribution has higher negative correlation with the elevation at lags larger than ␭. The distance of the transition from positive to negative cross correlation increases with ␭ at the expense of the decreasing magnitude of the cross correlation. Therefore the cross correlation between the topographic forcing and SM distribution is inversely proportional to the correlation scale of elevation. For ␭ ⬎ 103 m the topography might not have apparent effects on the horizontal redistribution of SM. Although this argument seems intuitively reasonable, it should be noted that ␭ depends on the resolution of elevation data. Finer resolutions of data (usually resulting in smaller ␭) are capable of capturing the variation of microtopography, which is shown to have significant impacts on SM distribution (see Figure 7). Finally, the combined impact of varying ␤ and ␭ on the cross correlation between ␪ and z is shown in the contour plot of Figure 9. This figure indicates the range of significant cross correlation between ␪ and z lies in ␭ ⬍ 103 m and ␤ ⬍ 10⫺3. The cross correlation between ␪ and z decreases as ␤ increases. In summary, the magnitude and scale of the cross correlation between topography and SM horizontal distribution are dependent on the correlation scale of elevation and the magni- YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION 1261 Figure 9. Contours of the one-dimensional one-point cross correlation for varying values of correlation scale of elevation ␭ and vertical hydraulic gradient ␤ (␣⫺1 ⫽ 0.5 m, D ⫽ 1 m). tude of vertical evaporative flux. Microtopography with a correlation scale smaller than 103 m has a significant impact on the horizontal SM distribution. 6. Conclusions An equation explicitly relating the variability of soil moisture distribution to the variability of topography under steady state homogeneous soil conditions is proposed and solved using the perturbation method. Both the one- and two-dimensional unsaturated flow equations are analyzed by assuming the exponential form and hole-type form of elevation covariance function in one dimension, and the Bessel-type covariance function of elevation in two dimensions. The covariance and variance of a soil moisture distribution, as well as the cross correlation between elevation and soil moisture distribution, are derived in closed form expressions for one and two dimensions. The steady state soil moisture distribution under homogeneous bare soil conditions is regulated by three distinct factors: topography, climate, and soil properties. First, topography forces a distribution of soil moisture that tends to mimic the elevation field at large scales. The other two factors are the vertical divergence of water in response to the climate forcing and the capillary resistance to water movement. The climate forcing, which is represented by a simple parameterization of vertical divergence of moisture flux, tends to smooth the spatial distribution of soil moisture. However, the capillary force exerted by the soil matrix tends to resist the displacement of water and hence exert adverse effects against the topography and climate forcings. The impact of capillary forces on the vertical fluxes of water seems to be more significant than their impact on the topographically induced horizontal fluxes. In addition to these three factors, soil depth regulates the impact of the soil moisture vertical flux such that a relatively large soil depth results in a relatively larger soil moisture variance. The variance of SM distribution decreases with elevation correlation scale ␭ and vertical hydraulic gradient ␤, and in- creases with variance of elevation ␴ 2z , capillary fringe thickness ␣⫺1, and soil depth D. Moreover, the two-dimensional analysis indicates that the SM variance is smaller than that in one dimension. Additionally, the hole phenomenon occurring in one-dimensional SM correlation function vanishes in two dimensions for typical range of elevation correlation scales. These two findings result from the greater freedom of soil moisture movement in two dimensions than in one dimension. The cross correlation between soil moisture horizontal distribution and topography decreases as the elevation correlation scale and the vertical evaporative fluxes increase. This trend corresponds to a shift of the dominant forcing on the soil moisture distribution from the topographically induced horizontal fluxes to climate-forced vertical fluxes. In contrast to the variance of soil moisture distribution, the cross correlation between soil moisture and elevation in two dimensions is higher than that in one dimension, especially for the condition of small elevation correlation scale and small vertical evaporative fluxes. This paper considered the steady state distribution of soil moisture assuming the homogeneous soil condition. Future research will focus on (1) the corresponding transient problem, (2) incorporating the variability of rainfall and (3) the impacts of soil heterogeneity and vegetation on the topographically induced soil moisture distribution. To test some of the theoretical concept proposed in this paper, a field experiment has been designed to measure soil moisture along a hillslope using neutron probe technology. The experimental site is at Harvard Forest, located in central Massachusetts. Results from this experiment will be reported in a forthcoming paper. Appendix A: Derivations of Equation (24) and (26) A1. Exponential Function From the S z (k) in (22) and the spectral relationship of soil moisture (SM) and elevation in (21), the spectrum of SM can be derived as 1262 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION S ␪共k兲 ⫽ ⫽ 2 ␴ 2z C 2␭ ␲ 2 ␴ 2z C 2 ␲␭ 冉 冋 k4 ␣␤ k2 ⫹ D 冊 Appendix B: Derivation of Equation (34) 2 By substituting the Whittle spectrum (equation (32)) into (21), we have (i.e., k 2 ⫽ k 2x ⫹ k 2y ) 共1 ⫹ ␭ 2k 2兲 a b c ⫹ ⫹ 共k ⫹ A兲 2 k 2 ⫹ A k 2 ⫹ B 2 册 (A1) S ␪共k兲 ⫽ ␴ 2z C 2␭ 2 ␲ 冉 冋 where 1 B⫽ 2 ␭ ␣␤ A⫽ D ⫺A 2 a⫽ A⫺B ⫽ k2 ⫹ 冕 A⫽ ⫹ b a e ⫺ 冑A兩␰兩共1 ⫹ 2 A 3/ 2 e ⫺ 冑A兩␰兩 ⫹ 冑A c 冑B 冑A兩 ␰ 兩兲 e ⫺ 冑B兩␰兩 册 (A2) ⫽ a b d c ⫹ ⫹ ⫹ 共k 2 ⫹ A兲 2 k 2 ⫹ A 共k 2 ⫹ B兲 2 k 2 ⫹ B where d⫽ ⫺2 AB 共 A ⫺ B兲 3 a⫽ ⫺A 3 共 A ⫺ B兲 2 ⫺B 3 c⫽ 共 A ⫺ B兲 2 b⫽ 2 ␴ 2z C 2 ␭2 e ik␰S ␪共k兲 dk 2 冑A c ␰ K 1共 冑B ␰ 兲 2 冑B ⫹ bK 0共 冑A ␰ 兲 ⫹ dK 0共 冑B ␰ 兲 册 g ␭3/ 2 4g ␭ ␰ K 共 g 1/ 2␰ 兲 ⫹ K 共 g 1/ 2␰ 兲 共 g ␭ ⫺ 1兲 2 ␭ 1 ␭ ␭ 共 g ␭ ⫺ 1兲 3 0 ␭ ␭ ⫹ 1 4g ␭ ␰ K 共␰ 兲 ⫺ K 共␰ 兲 共 g ␭ ⫺ 1兲 2 ␭ 1 ␭ 共 g ␭ ⫺ 1兲 3 0 ␭ 冕 ⬁ kS ␪共k兲 dk 0 a ⫺ 冑A兩␰兩 e 共1 ⫹ A 3/ 2 冑A兩 ␰ 兩兲 ⫹ 冑B兩 ␰ 兩兲 ⫹ d 冑B 2b 冑A ⫽ e ⫺ 冑A兩␰兩 e ⫺ 冑B兩␰兩 册 (A6) 冋 ⫺⬁ c ⫺ 冑B兩␰兩 e 共1 ⫹ B 3/ 2 a ␰ K 1共 冑A ␰ 兲 册 (A7) Since the SM variance cannot be evaluated directly from (A7), it should be derived by letting ␰␭ ⫽ 0 in the second line of (A6): B 2共3A ⫺ B兲 d⫽ 共 A ⫺ B兲 3 ␴ ␪2 ⫽ 2 冋 冋 By defining g ␭ ⫽ ␣␤␭ 2 /D and ␰ ␭ ⫽ 兩 ␰ 兩/ ␭ , (A6) can be written into the following dimensionless form as in (33): A 2共 A ⫺ 3B兲 共 A ⫺ B兲 3 ⬁ ␴ 2z C 2 ⫽ ␩ e i共kx␰x⫹ky␰y兲S ␪共k x, k y兲 dk x dk y 共 ␰ 2 ⫽ ␰ 2x ⫹ ␰ 2y 兲 R ␪共 ␰ ␭兲 ⫽ C 2␴ 2z 1 B⫽ 2 ␩ ␣␤ A⫽ D ⫹ B2 共 A ⫺ B兲 2 kJ 0共k ␰ 兲S ␪共k兲 dk ⫹ 册 (A3) 冕 c⫽ ⫺⬁ By the same procedure except taking S z (k) from (23) into (21), (A1) and (A2) becomes R ␪共k兲 ⫽ 2 AB 共 A ⫺ B兲 3 ⬁ Hole-Type Function 冋 1 ␭2 ⫺⬁3⬁ A2. 2␴ C S ␪共k兲 ⫽ ␲␩ B⫽ b⫽ 冕冕 冕 R ␪共 ␰ x, ␰ y兲 ⫽ ⫽ 2 册 The SM covariance function can be derived as By defining the non-dimensional parameter g ␭ ⫽ ␣␤␭ 2 /D and substituting A, B, a, b, c from (A1) into (A2), (24) is thus derived. 2 z ␣␤ D A2 共 A ⫺ B兲 2 e ik␰S ␪共k兲 dk 冋 共1 ⫹ ␭ 2k 2兲 2 a⫽ ⬁ ␴ 2z C 2 ␭ 2 where ⫺⬁ ⫽ 冊 b d a c ⫹ ⫹ ⫹ 共k 2 ⫹ A兲 2 k 2 ⫹ A 共k 2 ⫹ B兲 2 k 2 ⫹ B The decomposition by partial fraction in (A1) is for the convenience of integration. The SM covariance function can be obtained by taking inverse Fourier transform of (A1): R ␪共k兲 ⫽ ␣␤ D (A5) B2 c⫽ 共 A ⫺ B兲 2 A共 A ⫺ 2B兲 b⫽ 共 A ⫺ B兲 2 ␴ 2z C 2␭ 2 ␲ k4 2 ␴ 2z C 2 ␭2 冕 冋 ⬁ k b a ⫹ 共k 2 ⫹ A兲 2 共k 2 ⫹ A兲 0 ⫹ (A4) Thus (26) is derived and expressed in terms of the dimensionless parameter g ␩ ⫽ ␣␤␩ 2 /D. ⫽ c d ⫹ 共k 2 ⫹ B兲 2 共k 2 ⫹ B兲 ␴ 2z C 2 ␭2 冉 册 dk a c ⫹ b ln A ⫹ ⫹ d ln B A B 冊 YEH AND ELTAHIR: TOPOGRAPHY AND SOIL MOISTURE DISTRIBUTION ⫽ ⫽ ␴ 2z C 2 ␭2 ␴ 2z C 2 ␭2 ⫽ ␴ 2z C 2 冉 冉 冋 a c k2 ⫹ A ⫹ ⫹ b ln 2 A B k ⫹B a c B ⫹ ⫹ b ln A B A 冏冊 ⬁ 0 冊 g␭ ⫹ 1 2g ␭ ⫺ ln g ␭ 共 g ␭ ⫺ 1兲 2 共 g ␭ ⫺ 1兲 3 册 (A8) Notice that b ⫽ ⫺d from (A5). Thus (34) is obtained. Acknowledgments. The authors are grateful to the hydrology program of the National Aeronautics and Space Administration (NASA) for their support of the work under grant NAGW-4707. We also thank Guleid Artan (University of California, San Diego), Gabriel Katul (Duke University), and an anonymous reviewer. The manuscript benefited greatly from their helpful comments. References Bakr, A. A. M., L. W. Gelhar, A. L. Gutjahr, and J. R. MacMillan, Stochastic analysis of spatial variability in subsurface flows, 1, Comparison of one- and three-dimensional flows, Water Resour. Res., 14(2), 263–272, 1978. Beven, K. J., Kinematic subsurface stormflow, Water Resour. Res., 17(5), 1419 –1424, 1981. Beven, K. J., On subsurface stormflow: An analysis of response times, Hydrol. Sci. J., 27(4), 505–521, 1982. Beven, K. J., and M. J. Kirkby, A physically based, variable contributing area model of basin hydrology, Hydrol. Sci. J., 24(1), 43– 69, 1979. Bras, R. L., Hydrology; An Introduction to Hydrologic Science, 643 pp., Addison-Wesley, Reading, Mass., 1990. Entekhabi, D., and P. S. Eagleson, Land surfeace hydrology parameterization for atmospheric general circulation models, J. Clim., 2, 816 – 831, 1989. Gardner, W. R., Some steady state solutions of the unsaturated moisture flow equation with application to evaporation from a water table, Soil Sci., 85(4), 228 –232, 1958. Gelhar, L. W., Stochastic Subsurface Hydrology, 389 pp., Prentice Hall, Englewood Cliffs, N. J., 1993. Genereux, D. P., and H. F. Hemond, Naturally occurring radon 222 as a tracer for stream flow generation: Steady state methodology and field example, Water Resour. Res., 26(12), 3065–3076, 1990. 1263 Giorgini, A., M. Bergman, J. Pravia, and A. Hamidi, Lateral moisture movements in unsaturated anisotropic media, Tech. Rep. 169, 105 pp., Water Resour. Res. Cent., Purdue University, Lafayette, Indiana, 1984. Hewlett, J. D., and A. R. Hibbert, Moisture and energy conditions within a sloping soil mass, J. Geophys. Res., 68, 1081–1087, 1963. Hurley, D. G., and G. Pantelis, Unsaturated and saturated flow through a thin porous layer on a hillslope, Water Resour. Res., 21(6), 821– 824, 1985. Jackson, T. J., and T. J. Schmugge, Passive microwave remote sensing system for soil moisture: Some supporting research, IEEE Trans. Geosci. Remote Sens., 27(2), 225–235, 1989. Lowry, W. P., The falling rate phase of evaporative soil-moisture loss—A critical evaluation, Bull. Am. Meteorol. Soc., 40, 605– 608, 1959. Mantoglou, A., and L. W. Gelhar, Stochastic modeling of large-scale transient unsaturated flow systems, Water Resour. Res., 23(1), 37– 46, 1987. McCord, J. T., and D. B. Stephens, Lateral moisture flow beneath a sandy hillslope without an apparent impeding layer, Hydrol. Processes, 1, 225–238, 1987. Philip, J. R., The theory of infiltration, 1, The infiltration equation and its solution, Soil Sci., 83, 345–347, 1957. Philip, J. R., The theory of infiltration, Adv. Hydrosci., 5, 215–290, 1969. Pullan, A. J., The quasi-linear approximation for unsaturated porous media flow, Water Resour. Res., 26(6), 1219 –1234, 1990. Rodriguez-Iturbe, I., D. Entekhabi, and R. L. Bras, Nonlinear dynamics of soil moisture at climate scales, 1, Stochastic analysis, Water Resour. Res., 27(8), 1899 –1906, 1991. Solomon, A. M., and H. H. Shugart, Vegetation Dynamics and Global Change, 338 pp., Chapman and Hall, New York, 1993. Stagnitti, F., M. B. Parlange, T. S. Steenhuis, and J.-Y. Parlange, Drainage from a uniform soil layer on a hillslope, Water Resour. Res., 22(5), 631– 634, 1986. Zaslavasky, D., and G. Sinai, Surface hydrology, I, Explanation of the phenomena, J. Hydraul. Div. Am. Soc. Civ. Eng., 107, 1–16, 1981. E. A. B. Eltahir and P. J.-F. Yeh, R. M. Parsons Laboratory, Department of Civil and Environmental Engineering, Massachusetts Institute of Technology, Cambridge, MA 02139. (e-mail: eltahir@ mit.edu; patyeh@mit.edu) (Received April 18, 1997; revised January 5, 1998; accepted January 7, 1998.) WATER RESOURCES RESEARCH, VOL. 34, NO. 8, PAGE 2075, AUGUST, 1998 Correction to “Stochastic analysis of the relationship between topography and the spatial distribution of soil moisture” by Pat J.-F. Yeh and Elfatih A. B. Eltahir In the paper “Stochastic analysis of the relationship between topography and the spatial distribution of soil moisture” by Pat J.-F. Yeh and Elfatih A. B. Eltahir (Water Resources Research, 34(5), 1251–1263, 1998), the negative sign is missing in the cross-correlation equations (28), (29), (30), (31), (35), and (36). Accordingly, the sign in the ordinate of Figures 7 and 8, as well as in the contours of Figure 9, should be reversed. The cross correlation between elevation distribution and soil moisture distribution at lag zero should be negative. (Received May 20, 1998.) Copyright 1998 by the American Geophysical Union. Paper number 98WR01876. 0043-1397/98/98WR-01876$09.00 2075