HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS		QUICK SEARCH:   [advanced] Author: Keyword(s): Year:  Vol:  Page:

This Article Abstract Figures Only Full Text (PDF) Free Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Similar articles in ISI Web of Science Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via ISI Web of Science (13) Citing Articles via Google Scholar Google Scholar Articles by Timlin, D.J. Articles by Rawls, W. Search for Related Content PubMed Articles by Timlin, D.J. Articles by Rawls, W. GeoRef GeoRef Citation Agricola Articles by Timlin, D.J. Articles by Rawls, W.

DIVISION S-1-SOIL PHYSICS

Use of Brooks-Corey Parameters to Improve Estimates of Saturated Conductivity from Effective Porosity

^a USDA-ARS Systems Research Lab, Bldg. 007, Rm 116, 10300 Baltimore Ave, Beltsville, MD 20705 USA
^b USDA-ARS Great Plains Systems Research Unit, P.O. Box E, Ft. Collins, CO 80522 USA
^c Dept of Botany, Duke University, Durham, NC 27708 USA
^d USDA-ARS Grazing Lands Research Laboratory, P.O. Box 1199, El Reno, OK 73036-1199 USA
^e Dept. of Environmental Sciences, Rutgers, The State University of NJ, 14 College Farm Rd., New Brunswick, NJ 08901 USA
^f USDA-ARS Hydrology Laboratory, Bldg 007, Rm 112, 10300 Baltimore Ave., Beltsville, MD 20705 USA

dtimlin{at}asrr.arsusda.gov

	ABSTRACT

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
Effective porosity, defined here as the difference between satiated total porosity and water-filled porosity at a matric potential of 33 kPa, has been shown to be a good predictor for saturated hydraulic conductivity (K_s) using a modified Kozeny-Carman equation. This equation is of the form of a coefficient (B) multiplied by effective porosity raised to a power (n). The purpose of this study was to improve the predictive capability of the modified Kozeny-Carman equation by including information from moisture release curves (soil water content vs. matric potential relation). We fitted the Brooks-Corey (B-C) equation parameters (pore size distribution index and air entry potential) to moisture release data from a large database (>500 samples). Values of K_s were also available from the same source. Inclusion of the pore size distribution index into the Kozeny-Carman equation improved the K_s estimation over using only effective porosity, but only slightly. The improvement came through a better estimation of large values of K_s. We also fit a relationship for the coefficient (B) of the Kozeny-Carman equation as a function of the two B-C parameters with a constant value of n = 2.5 for the exponent. Overall the use of Brooks-Corey parameters from moisture retention data improved estimates of K_s over using effective porosity (_e) alone. There is still considerable error in predicting individual K_s values, however. The best forms of the equation was when was included in the term for the coefficient for the modified Kozeny-Carman equation. The next best form was when was included in the exponent for _e The two best models appeared to better preserve the mean, standard deviation and range of the original data.

Abbreviations: B-C, Brooks-Corey • RA, Rawls data • RMSE, root mean square error • SR, Southern region data

	INTRODUCTION

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
SATURATED SOIL HYDRAULIC CONDUCTIVITY (K_s) is an important soil parameter in models that simulate infiltration and runoff processes. This soil parameter is difficult to measure and can be highly variable, necessitating a large number of samples. For this reason indirect methods have held promise as an alternative to making direct measurements. A further advantage of indirect methods is that they allow researchers to obtain an estimate of the variability of saturated conductivities based on the variability of an easily measured predictor variable (Ahuja et al., 1989).

A number of relationships have been developed that can be used to calculate K_s with easily measured soil properties. Some are purely empirical and are often related to soil texture (Rawls et al., 1992; Puckett et al., 1985). Other relationships use physically based equations. Ahuja et al. (1984, 1989) showed that a modified Kozeny-Carman equation

(1)

was applicable to a wide range of soils from the Southern Region of the USA, Hawaii, and Arizona. Here _e is the effective porosity calculated as saturated water content (_s) minus the water content at 33 kPa matric potential, and B₁ and n are coefficients.

Even though the coefficients of Eq. [1] fitted to the data varied with soil type within a certain range, Eq. [1] fitted to K_s data for all nine different soils had an r² as good as for individual soils (Ahuja et al., 1989). In other words, Eq. [1] exhibited a degree of universality. In fact, the coefficients, B₁ and n obtained from the above fit of Eq. [1] to data for nine soils estimated K_s for several soils from Korea (Ahuja et al., 1989) and a variety of soils from Indiana (Franzmeier, 1991) with acceptable accuracy. Messing (1989) presented data for some Norwegian soils where Eq. [1] fitted the data for individual soils well, but the coefficients varied with soil type. Some of these soils had high clay contents and probably exhibited swelling–shrinking behavior, which could possibly affect the value of fitted coefficients. In any case, further research is needed to test and improve the universality of Eq. [1] for a variety of nonswelling soils and possibly for swelling soils.

Rawls et al. (1998) used Eq. [1] and proposed the use of the fractal dimension as an exponent (n) in an equation of the form:

(2)

Here _e is effective porosity and is the pore-size distribution index from the Brooks-Corey equation. Rawls et al. (1998) showed that the exponent (3 - ) can be considered to be a measure of the pore fractal dimension. When fit to mean values of _e and from soil textural classes, this equation was shown to give a good estimate of K_s with an r² of 0.92 and the intercept B₂ equal to 0.00053 m s^-1 (190 cm h^-1).

Earlier, Rawls et al. (1993) modified the Marshall equation to obtain an equation for matrix saturated hydraulic conductivity. These workers used a Sierpinski carpet generator to represent a two-dimensional soil matrix porosity and used water retention parameters from the Brooks-Corey equation. The modified Marshall's equation used by Rawls et al. (1993) is of the form:

(3)

Here l is a parameter related to the fractal dimension, is total porosity, x is an exponent, and R (cm) is the largest continuous pore radius for the Seirpinsky carpet. R is calculated from the capillary rise equation

(4)

In this equation h_b is the absolute value of the air-entry potential (cm). The value l is estimated from the fractal dimension D

(5)

Here D is the fractal dimension of soil porosity and is estimated as , and is the Brooks-Corey pore-size distribution index. Rawls et al. (1993) used a value of 4/3 for the exponent, x.

The modified Kozeny-Carman equation (Ahuja et al., 1984; Rawls et al., 1998) and the modified Marshall equation (Rawls et al., 1993) represent related approaches to indirect calculation of K_s. There is still a potential to improve the universality of the K - _e relationship. The objective of our study was to develop a combined method that uses information from the moisture release curve to improve Ahuja's modified Kozeny-Carman equation (Eq. [1]).

	Materials and methods

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
The data set used by Ahuja et al. (1984, 1989) is denoted here for convenience as the Southern Region (SR) data. These data were fully described in Ahuja et al. (1989). The names and taxonomic classifications of these soils are given in Table 1 . The data consist of moisture retention values (water content, , and pressure head, h), bulk density, soil texture, and saturated hydraulic conductivity measured on replicated (4–10), undisturbed soil cores taken from different soil horizons at several sites for each soil. These data are published in Southern Cooperative Series Bulletins or elsewhere (Bruce et al., 1983; Dane et al., 1983; Quisenberry et al., 1987; Nofziger et al., 1983; Green et al., 1982.). For all K_s data used in the present study, the soil cores were 60 to 85 mm in diameter and 60 to 75 mm in length. Constant-head methods were used to measure K_s and hanging water column and pressure plate procedures to measure water retention curves. The differences in water retention at saturation and at -33 kPa matric potential were used to obtain the effective porosity, _e. For some soils, where water retention at saturation was not initially measured, _s was calculated using measured soil bulk density and particle density data. The value of _s was calculated as 0.90 times total porosity. There were 571 sets of moisture retention data with associated values of K_s. The pressures for the moisture release curves ranged from 0.2 to 1500 kPa. Only curves with at least five retention values were used. The data set was also averaged by texture class as was the data of Rawls et al. (1993). The texture classes and the number of samples in each class are given in Table 2 .

Moisture Release Curve Parameters
Parameters for the Brooks-Corey equation were used to define the pore-size distribution index () and the air entry potential (h_b). These parameters were empirically derived from moisture release data. The Brooks-Corey equation is

(6)

where is pore-size distribution index, h_b is air entry potential (kPa), _r (cm³ cm^-3) is residual water content, and _s is saturated water content. , _r, and h_b were estimated using a combination of linear regression and a nonlinear optimization method (van Genuchten et al., 1981). A linearized form was obtained by taking a logarithmic transform of both sides of Eq. [6]:

(7)

and using the left-hand side of Eq. [7] as the dependent variable. A robust median fit linear regression algorithm (Press et al., 1986) was used to obtain values of and h_b, given an initial estimate of _r. A median fit regression method was used because the model is linear when _r is known. The nonlinear optimization program then iterated across a range of values of _r. New values of and h_b were fit for each new value of _r. This process continued until the value of _r that gave the smallest sum of squared differences (measured - observed) with corresponding values of and h_b was found. Only h- pairs where the absolute value of h was greater than 0.02 kPa were used. The root mean square error (RMSE) from the nonlinear optimization was <0.005 cm³ cm^-3 for 90% of the samples.

	Theory

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
An examination of Eq. [1], [2], and [3] show that they all have the form of a coefficient (or coefficients), which we will refer to as "B", multiplied by a measure of porosity (_e or ) that is raised to a power (n, 3 - , or x ). The Kozeny-Carman equation (Eq. [1]) can be parameterized by fitting the coefficient, B, and the exponent, n to measured K_s - _e data. Equations [2] and [3] are more complex since the coefficient and exponent are expressed as functions of additional variables, namely h_b and . The terms for h_b and in Eq. [3], when taken together, can be considered to be a coefficient, "B". We can consider Eq. [1] in the form

(8)

In this study, we consider two approaches to determining B(, h_b). If we consider B as a function of only, a possible candidate function for K_s is

(9)

Here A is an empirical coefficient and f() is a function to be derived later. Another possible candidate function for K_s where can be obtained by inspection of the modified Marshall equation (Eq. [3]):

(10)

Here R and l are defined as in Eq. [4] and [5] and A and m are empirical coefficients.

In order to investigate relationships among , h_b, and B, the coefficient, B can be expressed as

(11)

By expressing B in this manner we can investigate the possible functional relationships among , h_b, and B for constant n.

Statistical Calculations
Regressions were carried out using SAS software (SAS Institute, 1995). Except where noted and for model comparisons, regressions were carried out on log (base 10)–transformed input data. Comparisons of the models were accomplished by comparing residuals from a regression of the predicted K_s on measured K_s values from the SA data set. Transformed and untransformed values of K_s were used in these regressions. Comparisons of regression slopes were carried out by using a method given by Snedecor and Cochran (1980)(p. 387).

	Results and discussion

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
The Coefficient B as a Function of Alone
Figure 1 shows mean K_s as a function of for the RA and SR mean texture class data, where K_s is shown in logarithmic (base 10) scale. The relationships are similar for the two data sets except for small values of K_s and . Using all the data in the SR data set we fitted the function

(12)

View larger version (20K): [in this window] [in a new window]   Fig. 1 Pore-size distribution index vs. Ks for the Southern Region (SR) and Rawls (RA) data sets. Values are means for textural classes

Equation [12] suggests that the form of f() given in Eq. [9] is . Using all the data in the SR data set we fitted the following expression using regression on log-transformed values of _e and K_s

(13)

Here , and , RMSE of the . All the coefficients were significant (P < 0.001). Figure 2 shows the predicted and measured values for the original and modified relationships. Eq. [13] is similar to Eq. [1] where the coefficient B in Eq. [1] is replaced by C₃10^C₄. The RMSEs are given in Table 3 as a function of textural class. The largest change in RMSE is in the sand texture class, where the RMSE for estimates by the new Eq. [13] are less than the RMSEs for estimates by the original Eq. [1]. An inspection of the predicted and measured values in Fig. 2 shows that the predicted K_s values are closer to the 1:1 line for Eq. [13] (with ) for the largest values of K_s.

View larger version (19K): [in this window] [in a new window]   Fig. 2 Predicted and measured Ks for the Southern Region data set determined using (a) effective porosity (e) only as a predictor and (b) effective porosity (e) and the pore size distribution index ()

The Coefficient B as a Function of and h_b

View larger version (21K): [in this window] [in a new window]   Fig. 3 The intercept of the Kozeny-Carman equation (B) as a function of (a) , and (b) hb. The value of B is calculated as Ks (e2.5). The data are means for textural classes

With reference to the Rawls et al. (1993) modified Marshall equation (Eq. [3]) we used a combination of and h_b as predictors for B in the form of f(R/l) (Eq. [10]) where R and l are given by Eq. [4] and [5] respectively. We found that (R/l)^0.5 gave the best results with textural class mean data. Figure 4 shows the relationship between (R/l)^0.5 and the coefficient B from Eq. [11] for the RA and mean SR data sets. We fit a linear function for B(h_b,) that had an intercept of zero (Fig. 4). The slopes for the two relationships were not significantly different (F = 0.34). The relationship for the RA data set is:

(14)

View larger version (18K): [in this window] [in a new window]   Fig. 4 Relationship between the intercept B and R/l for Southern Region (SR) and Rawls (RA) data sets. The value of B is calculated as Ks (2.5e). The data are means for textural classes

(15)

Predicted K_s (m s^-1) values from Eq. [15] using Brooks-Corey parameters and _e values from the SR data set are plotted against measured K_s values in Fig. 5 . The relationship fits the data well with an for the log-transformed data that is similar to the r² for the original and modified Ahuja's relationships (Tables 3 and 4) , although the RMSE is slightly higher than for the original method. The parameters for this relationship were fit from the RA data and are completely independent of the SR data. The results could be improved using parameters fit to the SR data set. However, the differences will not be large since the slopes of the relationships for the two data sets shown in Fig. 4 are similar. The applicability of Eq. [15] to the SR data set where the coefficients were derived from the RA data set does demonstrate the generality of this relationship. It is also encouraging to note that the error is not that much larger than the error in estimates from the other models fit to the SR data set.

View larger version (22K): [in this window] [in a new window]   Fig. 5 Calculated and measured values of Ks, where the intercept (B) for Eq. [1] has been calculated using values of and hb, and the relationship between B and R/l for the Rawls data set given in Fig. 4

View larger version (23K): [in this window] [in a new window]   Fig. 6 Relationship between measured Ks and values predicted using Eq. [2] when fit to the Southern Region (SR) data set

The results given in Tables 3 and 4 suggest that a particular equation may be more applicable to a specific range of data. Regression with log-transformed variables acts like a weighted regression, where small values of K_s are given higher weights than would be given if the fitting method evaluated deviations of untransformed values. As a result, Eq. [1] and [13] may be more appropriate for estimating small values of K_s in finer-textured soils for example, and Eq. [2] for larger values of K_s in coarser-textured soils.

The probability plots in Fig. 7 indicate how well each model describes the original distribution of data. The distribution of the predicted values of the four models are not greatly different in the mid ranges of the data. Equations [2] and [13] both come closest to the distribution at high values of K_s. The distribution of K_s predicted using a calculated value of B (Eq. [15]) is quite close to the SR data distribution for lower values of K_s. This is encouraging considering this equation was parameterized with an independent data set. Of the four models, the distribution of K_s predicted by Eq. [13] appears closest to the distribution of measured K_s values. However, these differences are not large but are important to consider when an estimation is used to generate a distribution of K_s values across a field as a function of spatial variability.

View larger version (25K): [in this window] [in a new window]   Fig. 7 Probability plot of measured Ks and values estimated by the four models. Ks values are from the full Southern Region (SR) data set

There is also the question of the use of these methods to estimate saturated conductivities for use at the field scale. Ahuja et al. (1993) have shown that a harmonic mean K_s of layered soil can be estimated from a 2-d drainage of surface soil. This may extend the usefulness of the methods developed in our study. Ahuja has observed (Ahuja, 1993, unpublished data) that final infiltration rates taken in 25-cm-diameter infiltration rings are related to average effective porosity, _e, of the 1-m profile, as well as soil water content of the profile measured 2 to 3 d after the soil was fully wetted. These relationships are similar to Ahuja et al.'s (1989) K_s(_e) relationships. Soil water content measured 2 to 3 d after rainfall is probably a better estimator of drainable porosity than 33 kPa water content. It is likely that this property can be easily scaled up to larger areas. Research into this area would be a promising extension to this work.

	Summary and conclusions

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 
The modified Kozeny-Carman equation was used to calculate saturated conductivity (K_s) from effective porosity (_e). We were able to obtain better predictions of K_s when the pore size distribution index () from the Brooks-Corey equation was used along with _e. The coefficient of determination (r²) for log(K_s) increased from 0.70 to 0.73 and the RMSE of log(K_s) decreased from 0.60 to 0.57. The use of improved the fit for larger values of K_s (>2.5 x 10^-5 m s^-1).

We determined a functional relationship for the intercept B (K_s^-2.5) in the modified Kozeny-Carman equation (with an exponent of 2.5) as a function of h_b and from the Brooks-Corey equation. The equation for B vs. f(, h_b) was linear with an intercept of 0 when fit to a data set of textural mean values that was available in the literature. Independent predictions of K_s using parameters from a data set different from the one for which the B vs. f(, h_b) relationship was fit had an r² of 0.75 and RMSE = 0.59 for log(K_s).

Overall, the use of Brooks-Corey parameters from moisture retention data improved estimates of K_s, compared with using effective porosity (_e) alone. However, in spite of the improvement, there is still substantial prediction error. There was not a large difference in prediction error among the four models. The best form of the equation was when the Brooks-Corey pore-distribution factor, , was included in the term for the coefficient of the modified Kozeny-Carman equation. The next best form was when was included in the exponent for _e. The use of the air entry potential (h_b) did not measurably improve the estimates of K_s. The RMSEs for the two best models were not greatly different. The two best models appeared to retain the distribution of the original data better. The use of improved the estimate of K_s in coarse-textured soils, and so models that incorporate would be more appropriate for estimation of K_s in this texture class.van Genuchten 1981

Received for publication November 24, 1997.

	REFERENCES

TOP ABSTRACT INTRODUCTION Materials and methods Theory Results and discussion Summary and conclusions REFERENCES

 

• Ahuja L.R., Cassel D.K., Bruce R.R., Barnes B.B. Evaluation of spatial distribution of hydraulic conductivity using effective porosity data. Soil Sci. 1989;148:404-411.
• Ahuja L.R., Naney J.W., Green R.E., Nielsen D.R. Macroporosity to characterize spatial variability of hydraulic conductivity and effects of land management. Soil Sci. Soc. Am. J. 1984;48:699-702.[Abstract/Free Full Text]
• Ahuja L.R., Wendroth O., Nielsen D.R. Relationship between drainage of surface soil and average profile saturated conductivity. Soil Sci. Soc. Am. J. 1993;57:19-25.
• Bruce, R.R., J.H. Dane, V.L. Quisenberry, N.L. Powell, and A.W. Thomas. 1983. Physical characteristics of soils in the Southern Region: Cecil Southern Coop. Ser. Bull. 267. Georgia Agric Exp. Stn., Univ. of Georgia, Athens.
• Dane, J.H., D.K. Cassel, J.M. Davidson, W.L Pollans, and V.L. Quisenberry. 1983. Physical characteristics of soil in the Southern Region Troup and Lakeland series. Southern Coop. Ser. Bull. 262. Alabama Agric. Exp. Stn., Auburn.
• Franzmeier D.P. Estimation of hydraulic conductivity from effective porosity data for some Indiana soils. Soil Sci. Soc. Am. J. 1991;55:1801-1803.[Abstract/Free Full Text]
• Green R.E., Ahuja L.R., Chong S.K., Lau L.S. Water conduction in Hawai'i oxic soils. Water Resources Research Center. Technical Report 143. Honolulu: Univ. of Hawaii at Manoa, 1982.
• Messing I. Estimation of saturated hydraulic conductivity in clay soils from soil moisture retention data. Soil Sci. Soc. Am. J. 1989;53:665-668.[Abstract/Free Full Text]
• Nofziger, D.L., J.R. Williams, A.G. Hornsby, and A.L. Wood. 1983. Physical characteristics of soils of the Southern Region. Southern Cooperative Series Bull. 265. Central Mailing Services, Oklahoma State Univ., Stillwater.
• Press W.H., Flannery B.P., Teukolsky S.A., Vetterling W.T. Numerical Recipes. New York: Cambridge Univ. Press, 1986.
• Puckett W.E., Dane J.H., Hajek B.F. Physical and mineralogical data to determine soil hydraulic properties. Soil Sci. Soc. Am. J. 1985;49:831-836.[Abstract/Free Full Text]
• Quisenberry, V.L., D.K. Cassel, J.H. Dane, and J.C. Parker. 1987. Physical characteristics of soils in the Southern Region: Norfolk, Dothan, Goldsboro, Wagram. Southern Coop. Ser. Bull. 263. South Carolina Agric. Exp. Stn., Clemson Univ., Clemson, SC.
• Rawls, W.J., L.R. Ahuja, and D.L. Brakensiek. 1992. Estimating soil hydraulic properties from soils data. p. 329–340. In M. Th. van Genuchten et al. (ed.) Proc. Int. Workshop on Indirect Methods for Estimating the Hydraulic Properties of Unsaturated Soils. Riverside, CA. 11–13 Oct. 1989. USDA and Univ. of California, Riverside.
• Rawls W.J., Brakensiek D.L., Logsdon S.D. Predicting saturated hydraulic conductivity utilizing fractal principles. Soil Sci. Soc. Am. J. 1993;57:1193-1197.[Abstract/Free Full Text]
• Rawls W.J., Brakensiek D.L., Saxton K.E. Estimation of soil water properties. Trans. ASAE 1982;25:1316-1320.
• Rawls W.J., Gimenez D., Grossman R. Use of soil texture, bulk density, and slope of the water retention curve to predict saturated hydraulic conductivity. Trans ASAE 1998;41:983-988.
• SAS Institute. 1995. SAS/STAT guide for personal computers. Version 6.07 ed. SAS Inst., Cary, NC.
• Snedecor G.W., Cochran W. Statistical methods, 7th ed Ames: Iowa State Univ. Press, 1980.
• van Genuchten, M. Th. 1981. Non-equilibrium transport parameters from miscible displacement experiments. Research Rep. 119. U.S. Salinity Lab., USDA-SEA-ARS, Riverside, CA.

This article has been cited by other articles:

				L. R. Ahuja, L. Ma, and D. J. Timlin Trans-Disciplinary Soil Physics Research Critical to Synthesis and Modeling of Agricultural Systems Soil Sci. Soc. Am. J., February 2, 2006; 70(2): 311 - 326. [Abstract] [Full Text] [PDF]

				A. Nemes, W. J. Rawls, and Y. A. Pachepsky Influence of Organic Matter on the Estimation of Saturated Hydraulic Conductivity Soil Sci. Soc. Am. J., June 28, 2005; 69(4): 1330 - 1337. [Abstract] [Full Text] [PDF]

				B. V. Iversen, B. V. Iversen, P. Moldrup, P. Schjonning, and O. H. Jacobsen Field Application of a Portable Air Permeameter to Characterize Spatial Variability in Air and Water Permeability Vadose Zone J., November 1, 2003; 2(4): 618 - 626. [Abstract] [Full Text] [PDF]

				B. Minasny and A. B. McBratney The Neuro-m Method for Fitting Neural Network Parametric Pedotransfer Functions Soil Sci. Soc. Am. J., March 1, 2002; 66(2): 352 - 361. [Abstract] [Full Text] [PDF]

This Article Abstract Figures Only Full Text (PDF) Free Alert me when this article is cited Alert me if a correction is posted Services Similar articles in this journal Similar articles in ISI Web of Science Alert me to new issues of the journal Download to citation manager Citing Articles Citing Articles via HighWire Citing Articles via ISI Web of Science (13) Citing Articles via Google Scholar Google Scholar Articles by Timlin, D.J. Articles by Rawls, W. Search for Related Content PubMed Articles by Timlin, D.J. Articles by Rawls, W. GeoRef GeoRef Citation Agricola Articles by Timlin, D.J. Articles by Rawls, W.