Published online 21 June 2006
Published in Soil Sci Soc Am J 70:1262-1271 (2006)
DOI: 10.2136/sssaj2005.0247
© 2006 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
Soil Physics
Fitting Uniaxial Soil Compression Using Initial Bulk Density, Water Content, and Matric Potential
D. D. Fritton*
Dep. of Crop and Soil Sciences, Pennsylvania State Univ., University Park, PA 16802
* Corresponding author (ddf{at}psu.edu)
 |
ABSTRACT
|
|---|
The existing bulk density and wetness of a soil influence soil compaction, but the ability to quantitatively describe these effects is limited. The purpose of this study was to generalize an empirical nonlinear equation for uniaxial soil compression, presently limited to fitting one compression curve at a time, by adding a capability to describe a single soil material at variable initial bulk density and wetness. The original relationship between the compressed soil bulk density and the applied stress was first simplified by reducing the number of coefficients required to fit a single experimental curve from three to two. The generalized relationship combined the original equation, an equation that removed a high correlation that existed between two of its coefficients, and two equations that reflected the impact of the bulk density, water content, and matric potential on the two remaining coefficients. This generalized relationship was tested for soil materials at various initial bulk density, water content, and matric potential values where texture and organic matter content were held constant. The study used 21 datasets representing three to five horizons of four soils and one soil mixed with four different amounts of sand. For each horizon or soil material studied, the generalized relationship was fitted simultaneously using nonlinear regression to 2–14 compression curves per horizon, including disturbed and undisturbed soil samples in nine cases. The generalized relationship fit each dataset with an R2 
0.932 (P < 0.001) and was judged superior to the best existing equation for multiple curves.
|
INTRODUCTION
|
|---|
MANY attempts to quantify soil compaction rely on a semi-log linear equation of soil compression. In this equation, the mechanical properties of the soil are represented by the compression index (the slope of the straight line relationship between the measure of compaction and the log of the applied stress) and the intercept (the value of the measure of compaction at the point where the log of the applied stress is zero). Once the soil has experienced any applied stress, it is necessary to define the preconsolidation stress (the value of stress below which the soil no longer follows the initial semi-log linear equation) and the elastic rebound/recompression parameter (the slope of a semi-log linear relationship between the measure of compaction and the log of the applied stress where the stress is less than the preconsolidation stress).
In the majority of cases, approaches using the semi-log linear equation have dealt with soil that is loose before compaction, such as a sieved soil or a tilled soil being prepared as a seedbed or very wet soils of low soil strength, and have recognized the need to formulate the compression index and the intercept as functions of soil wetness. In a recent example, Défossez et al. (2003) fitted soil compaction during seedbed preparation for numerous dates with different soil water contents. They used a fourth-order polynomial function to represent the effect of water content on the compression index and the intercept using 10 coefficients to characterize the soil mechanical properties. They also found that there was a correlation between the compression index and the intercept, indicating that the semi-log linear equation becomes over-parameterized when extended to multiple water contents.
This semi-log linear equation and its related mechanical properties are also part of the critical-state approaches to soil behavior (Rashid and Abedin, 2000; Wulfsohn and Adams, 2002; Toll and Ong, 2003; Gallipoli et al., 2003). Kirby (1994) and Eco (2002) have used critical state soil mechanics to model uniaxial compression. The critical state assumption of the linear semi-log equation, which causes a sharp transition from elastic to plastic behavior at the preconsolidation stress, caused problems in both simulations. Kirby (1994) needed to modify the Modified Cam Clay model to allow a gradual transition from elastic to plastic behavior. Eco (2002) simulated compression using eight coefficients for one unsaturated (suction of 300 kPa) sample of Sainte-Rosalie clay. Matric suction, degree of saturation, and the soil water characteristic curve were used as input into the simulation model to continuously update soil mechanical properties at the preconsolidation stress to generate the curved region around the preconsolidation stress area.
To overcome the limitations of the linear semi-log equations for uniaxial soil compression, several authors have introduced nonlinear equations. These nonlinear equations (Bailey et al., 1986; Assouline et al., 1997; Fritton, 2001) introduce the initial bulk density as an important soil property influencing soil compression and all three describe data below the preconsolidation stress and data above the preconsolidation stress and provide a gradual transition between the two stress regions. Recently, Assouline (2002) modified the Assouline et al. (1997) nonlinear equation to improve its fit to uniaxial soil compression data. These three equations allow soil mechanical parameters to be derived to characterize sites that have experienced prior compactive events, such as no-till fields, pastures, turfgrass areas, and subsoils.
In the earliest of the three nonlinear equations, Bailey et al. (1986) plotted the initial bulk density and all three coefficients in their equation against water content for a clay soil, indicating that all four coefficients were sensitive to water content but did not attempt to provide equations for the relationships. McBride (1989) investigated uniaxial compression for undisturbed core samples from 34 horizons in 14 soil series in the Niagara area of Ontario, Canada. Using the Bailey et al. (1986) equation, McBride used water content at sampling and plasticity properties to partition the data into four groups. Within all but one group, plastic soils with a low initial water content and a plastic index greater than 15, the water content of individual core samples significantly affected one or more of the equation coefficients. A more useful generalization of the Bailey et al. (1986) equation was developed by McNabb and Boersma (1993, 1996). They modified the Bailey et al. (1986) equation to account for the effects of the initial bulk density and the water content at the time of sampling for three Andisols and one Ultisol. The resulting equation used six coefficients. They showed that the sampling water content significantly improved the prediction of soil compression over and above a significant effect of the initial bulk density. In contrast, Imhoff et al. (2004) found no effect of water content for a toposequence of Oxisols after the initial bulk density effect had been removed using the McNabb and Boersma (1993) equation.
In general, attempts to quantitatively describe multiple compression curves that result from soil samples at different initial conditions increases the number of coefficients needed over that used for a single curve. This results from the need to describe the effect of several soil physical properties on the coefficients of the underlying equation for the compression curve for a single sample. Single curves can be described by as few as two regression coefficients (Assouline, 2002 and Eq. [3] with n1 and n2 fixed in this article). Multiple curves have required six regression coefficients (McNabb and Boersma, 1996), ten coefficients (Défossez et al., 2003), eight coefficients (Eco, 2002), and the two to six coefficients used in this article.
The Fritton (2001) and Assouline (2002) nonlinear equations provide accurate fits to uniaxial soil compression data for individual disturbed or undisturbed soil samples at a known initial bulk density and water content. They fit uniaxial soil compression data for individual soil samples more accurately than the Bailey et al. (1986) equation. However, they have not been generalized using functional relationships between their coefficients and initial soil properties. The purpose of this article is to generalize the Fritton (2001) equation for soil compression by adding a capability to describe a single soil material at variable initial bulk density and wetness.
|
STATISTICAL CONSIDERATIONS AND HYPOTHESES
|
|---|
Fritton (2001) represented soil compression for a single soil sample at an initial bulk density of 
o (Mg m–3) and an initial level of water content by the sigmoid curve of van Genuchten (1980):
 | [1] |
where is the soil bulk density (Mg m–3) at an applied stress of 
(kPa), m is the particle density (Mg m–3), 
is an empirical coefficient (kPa–1), and n and m are unitless empirical coefficients. This relationship fit (R2 > 0.97) bulk density data for 99 disturbed and 21 undisturbed individual soil samples undergoing uniaxial compression at applied stresses ranging from 0 to 2971 kPa. Based on the points made by Assouline (2002), the ( + 1) term in the original formulation of Eq. [1] has been replaced by as the independent variable. Statistical results from Fritton (2001) are used unchanged for this article. The mass of the brass pressure distribution plate and porous stone placed on each sample at the start of the compression process and considered negligible in the previous paper added 1.45 kPa to each applied stress value. This value was rounded to 1 kPa and added to the input applied stress values and results in identical values for , n, and m for each sample and avoids the problem of using a zero value of stress in the statistical calculations.
The values for , n, and m generated by nonlinear regression for Eq. [1] are not unique when the values are highly or even moderately correlated (McNabb and Boersma, 1996; Scheinost et al., 1997). Based on the work of McNabb and Boersma (1993), an approach was used that builds on the correlation between coefficients contained in Eq. [1]. It is hypothesized that n = f(m), = f(o), and m = f(o). In addition, it is hypothesized that and m are functions of water content—
, kg kg–1 (McNabb and Boersma, 1996; McBride, 1989; Défossez et al., 2003) and matric potential—
, kPa (Eco, 2002; Gallipoli et al., 2003). Specifically, the empirical relationships to be tested are
 | [2] |
where n1 and n2 are dimensionless regression coefficients, and
 | [3a] |
where
 | [3b] |
and
 | [3c] |
where a, b, c, cc, bb, d, e, f, ff, and ee are regression coefficients with the remaining coefficients and independent variables as previously defined. Equations [2], [3b], and [3c] resulted from preliminary statistical regressions using a variety of other possibilities, including linear and log transformed values for and m evaluated as functions of the individual independent variables o, , and in various polynomial relationships. In Fritton (2001), the coefficient was shown to be a direct measure of the inverse of the preconsolidation stress for a soil. The coefficients in Eq. [3b] are direct measures of a hypothesized quadratic effect of the initial bulk density (coefficients b and bb), linear effect of the initial water content (coefficient c), and linear effect of initial matric potential (coefficient cc), with an intercept (coefficient a) on the log-transformed value of . Similarly, the coefficient m was shown (Fritton, 2001) to be directly related to the compression index of a soil. The coefficients d, e, f, ff, and ee are direct measures of the hypothesized effect of the initial bulk density, water content, and matric potential on the compression index and more generally on the shape of the compression curve. These relationships are similar to those summarized by Imhoff et al. (2004).
In this study, the sample preparation procedure involved the satiation and subsequent dewatering of loose, air-dried soil for the 99 sieved (disturbed) samples, and it is hypothesized that this caused a shrink-swell/hydroconsolidation/densification response. The empirical relationship to be tested is a sigmoidal equation (Thornley, 1976) similar to one used by Peng and Horn (2005)
 | [4] |
where o is the initial soil bulk density (Mg m–3) at a water content of (kg kg–1), wet is the soil bulk density at a high water content, dry is the soil bulk density at a low water content, c is the water content for a half-maximal response (used as an empirical coefficient), and s is an empirical coefficient. In addition, Eq. [3] is compared with the best existing equation for multiple compression curves (McNabb and Boersma, 1996).
 | [5] |
where c is the compressed bulk density (Mg m–3); o is the bulk density at zero stress, which is used as a best fit regression coefficient in the equation (Mg m–3); 
i = i a–1 normalizes o for differences in initial bulk density where i is the initial bulk density and a is the average initial bulk density; is the applied stress (kPa); c = (i – 1)o adjusts the compression curve for differences in initial bulk density of each sample; m = (1 – s) with s being the saturation ratio; and A, B, C, D, and E are regression coefficients.
For the comparison of empirical equations fitted to the same set of data, the Akaike Information Criterion (AIC) was selected.
 | [6] |
It is sensitive to the number (N) of data points being fit, the number (M) of optimized coefficients in the empirical equation, and the residual sum of squares (RSS) of deviations from the fitted equation (
im
nek and Hopmans, 2002). The best equation is the one that minimizes AIC.
The possibility of multicollinearity and over-parameterization was evaluated using the concept (GraphPad Software, 2005) of dependency (DEP)
 | [7] |
where SEF is the standard error of a given coefficient in the equation when all other coefficients have been fixed at their optimum values and nonlinear regression is run with only the given coefficient being optimized, DFF is the number of degrees of freedom when the given coefficient is the only one being optimized, DF is the number of degrees of freedom when all coefficients are being optimized simultaneously, and SE is the standard error of the given coefficient when it is being optimized with all other coefficients simultaneously. The standard error for a coefficient decreases when all other coefficients are fixed, so DEP takes on values ranging from 0 to 1. Values of DEP > 0.99 indicate excessive multicollinearity (GraphPad Software, 2005).
|
MATERIALS AND METHODS
|
|---|
The experimental procedures and data used in this article are identical to those described in Fritton (2001). A brief description is included in this article.
Soil Properties
Four deep, well drained soils (Rayne silt loam [fine-loamy, mixed, mesic Typic Hapludult] formed in gray shale residuum, Bucks silt loam [fine-loamy, mixed, mesic Typic Hapludult] formed in red shale residuum, Glenelg silt loam [fine-loamy, mixed, semiactive, mesic Typic Hapludult] formed in mica schist residuum, and Hagerstown silt loam [fine, mixed, semiactive, mesic Typic Hapludalf] formed in limestone residuum) were used for this study. For one set of measurements, a commercial white quartz sand was mixed with the sieved (<2 mm) Hagerstown (B horizon, 0.85–1.00 m) soil in dry mass (sand/soil) ratios of 1:4, 2:3, 3:2, and 4:1 to extend the range of particle size distributions studied.
For the tested materials (see Table 1 in Fritton [2001] for details), the organic C content is somewhat higher in surface horizons than the underlying horizons but is generally low. The clay (<0.002 mm) content varies from 91 to 486 g kg–1, the silt (0.002–0.05 mm) content is relatively high in most samples and varies from 87 to 594 g kg–1, and the sand (0.05–2 mm) content is rather evenly distributed across the various sand fractions in most samples and varies from 110 to 822 g kg–1. The rock fragment content varies from 0 to 347 g kg–1 but is applicable only to undisturbed (core) samples because all disturbed (sieved) samples were run on <2.0-mm material.
Soil Preparation
Sieved soil (99 of the 120 samples) was poured into 63.5-mm-diameter rings that were approximately 5 mm taller than the sample height (25.4 mm) needed for compression and leveled gently. Replicated undisturbed samples (21 of the 120 total samples) were taken on five Glenelg and four Bucks soil horizons. These undisturbed samples were trimmed with a wire saw until the ends were flush with the ring. All soil samples were placed on a wet ceramic water extraction plate and satiated from the bottom by ponding water on the plate surface. After wetting for 1–2 d, pairs of soil samples were placed in pressure chambers, and water was extracted for 2 or 3 d at constant air pressures (nominally at 10, 50, 100, 300, and 500 kPa) to obtain one to four different water contents. It was assumed that the negative of the applied air pressure was a measure of the matric potential of the water in the sample. The upper four horizons of the Hagerstown series were run at only one air pressure (500 kPa). For the 1 part Hagerstown/4 parts sand dataset, one sample of the four was removed early from the pressure apparatus to obtain a sample with higher water content. One of each pair of samples was oven dried at 105°C to estimate the water content (kg kg–1), and the other was used for compression determinations because water content decreased during compression in the wetter samples.
Compression Apparatus and Procedure
Uniaxial drained compression using a triple-beam apparatus (Model C-320; ELE International, Pelham, AL) was run on each sample. A dial gauge (Model LC-3M; ELE International) with a measurement precision of 0.025 mm was read after each sample had been loaded at each level of stress for at least 30 min and used to calculate the bulk density (estimated error < 0.02 Mg m–3; see Fritton [2001]). Applied stress levels of 0, 31, 62, 93, 186, 371, 557, 743, 1114, 1485, and 2971 kPa were used for the Bucks, Glenelg, and Rayne soils by one individual. A second individual ran all Hagerstown and Hagerstown/sand mixtures several years later using fewer applied stress levels (0, 31, 62, 186, 557, 1114, and 2971 kPa).
Statistical Calculations
All statistical calculations were made using the NonlinearFit package (Boyland et al., 1992) in Mathematica (Wolfram, 1991). The Levenberg-Marquardt method to minimize the error sum of squares was chosen within the NonlinearFit package.
Equation [2] was tested using the entire database of 120 samples. Values of n and m were obtained from fits of Eq. [1] for each sample from Fritton (2001).
The data to test Eq. [3] was partitioned by soil series, horizon, and material resulting in 21 datasets representing soils with constant texture and organic matter content. Disturbed and undisturbed data for a horizon were used together as one dataset. For the Glenelg (0.30–0.53 m) dataset, three samples that had been used in the Eq. [2] test were removed to eliminate confounding between data generated before and after a long period of storage. For the upper four horizons of the Hagerstown series, the matric potential coefficients cc and ff were not used in Eq. [3] because only one air pressure had been used to prepare the samples. For the 1 part Hagerstown/4 parts sand dataset, the matric potential coefficients cc and ff were not used in Eq. [3] because one sample of the four had been removed early from the pressure apparatus to obtain a sample with higher water content.
A stepwise approach was used to determine the coefficients of Eq. [3]. The coefficient m in Eq. [3a] was first set at 2.65 Mg m–3 for each regression. The bulk density at zero applied stress, o, was set equal to the bulk density measured under the initial stress of 1 kPa and used along with the initial water content, (kg kg–1), and the matric potential, (kPa), as the input data for each soil sample (see Table 1 for an example of the input data used for the Bucks 0.76- to 1.22-m horizon for its seven samples). Nonlinear regression was run using Eq. [3] with n1 and n2 fixed (see Results) for the set of bulk density versus applied stress data at stresses greater than the initial stress of 1 kPa because the bulk density at the initial stress was input as data for a given dataset.
The results of the nonlinear regression were examined. If any coefficients were not statistically significant at the 5% level and/or if any coefficients were correlated >0.99± (McNabb and Boersma, 1996), one coefficient at a time was removed from Eq. [3], starting with the coefficient judged to have been least effective in the regression and the nonlinear regression was rerun. This iterative procedure was continued until all remaining coefficients were statistically significant at the 5% level and until none were correlated > 0.99±.
At this point in the procedure, the resulting equation was checked for over-parameterization using Eq. [7]. If the dependency (DEP) of any coefficient was >0.99, the coefficient judged to contribute most to the dependency was removed from Eq. [3], and the nonlinear regression was rerun. The removal of coefficients continued until all remaining coefficients were statistically significant at the 5% level, correlated <0.99±, and had a dependency of <0.99.
At times during this procedure, it became evident that the three criteria for a successful equation could not be met. When this happened, the process was restarted using Eq. [3], with n1 and n2 being allowed to vary instead of being fixed, and the procedure was repeated.
There were 10 non-zero stress levels per compression curve for the Bucks, Glenelg, and Rayne soils and six non-zero stress levels per compression curve for the Hagerstown and Hagerstown/sand mixtures. For two datasets, there was a missing data point for one stress level on one of the compression curves.
Equation [4] was evaluated using nonlinear regression on the 99 disturbed samples in the total dataset and on the disturbed samples for the Bucks 0.76- to 1.22-m horizon.
The coefficients for Eq. [5] were determined by nonlinear regression for selected datasets. When a coefficient was not significant at the 5% level, it was removed from Eq. [5], and the regression was rerun. The asymptotic correlation matrix was examined, and there were no cases of coefficients correlated >0.99±.
Best fit regression coefficients reported in Tables 1, 2, and 3 were rounded to three significant digits. Test cases indicate that three significant digits were adequate to control any increase in the mean square error for the equations to less than 3% over that resulting from using six significant digits. The mean absolute difference was the average of the absolute values of the difference between the experimental bulk density for each point for all stress levels greater than 1 kPa and the value calculated from the final regression equation for that point. The coefficient of determination (R2) and mean absolute differences were also rounded to three significant digits but were based on results using six significant digits.
View this table:
[in this window]
[in a new window]
|
Table 2. Best-fit coefficients determined by nonlinear regression for Eq. [3]. Dashes indicate coefficients were not significant at the 0.05 level, correlated >0.99 ± with other coefficients, and/or had a dependency > 0.99 and were removed from Eq. [3] for the regression results shown. Under the n coefficients column, fixed means n1 was set at 0.520, and n2 was set at 0.120, and neither was allowed to vary during the regression. All equations were highly significant (P < 0.001).
|
|
View this table:
[in this window]
[in a new window]
|
Table 3. Best fit coefficients determined by nonlinear regression for Eq. [5] for selected datasets from Table 2. Dashes indicate coefficients that were not significant at the 0.05 level and were removed from Eq. [5] for regression results shown.
|
|
|
RESULTS
|
|---|
General Patterns
Figure 1
illustrates the bulk density measured for the Bucks soil at the 0.76- to 1.22-m depth at the initial stress of 1 kPa and at the 10 applied stress levels of the uniaxial compression test. The initial bulk densities used to generate the smooth curves from Eq. [1] are shown in Table 1 along with values for , n, and m for each soil sample. A total of 21 regression coefficients are required (three for each sample) for the curves shown in Fig. 1.
View larger version (28K):
[in this window]
[in a new window]
|
Fig. 1. Bulk density of seven Bucks soil samples from the 0.76- to 1.22-m depth plotted as a function of the applied stress. The legend indicates the initial water content (kg kg–1) for each soil sample. Symbols represent experimental data. The smooth lines are nonlinear regression curves based on Eq. [1] best fit separately (see Table 1 for values of the input data and coefficients) to each individual soil sample.
|
|
The three sets (open symbols) of data in Fig. 1 for undisturbed cores had an initial bulk density near 1.60 Mg m–3 (Table 1) and are typical of the 21 sets of data on undisturbed cores from each of the other Bucks and Glenelg datasets except that the initial bulk density was more variable in most other datasets.
The four disturbed (sieved) sets of data (filled symbols) in Fig. 1 are broadly representative of the disturbed data from other datasets. In this case, the soil sample with the lowest initial bulk density (1.20 Mg m–3) had the highest water content (0.235 kg kg–1) and the highest matric potential (–10 kPa), and each curve for a soil sample at a higher initial bulk density had a lower water content and a lower matric potential (Table 1). This relationship resulted from shrinkage of the samples during water extraction in the pressure chamber and was fitted (R2 = 0.986, P < 0.001) by Eq. [4]. It is evident that statistically significant coefficients for the initial bulk density terms in Eq. [3] potentially explained variation caused by water content and/or matric potential and vice versa.
The shrinkage pattern indicated previously for the Bucks soil from Fig. 1 is evident in 6 of the 21 datasets. Eight additional datasets show the same trend for about three out of four individual samples. The remaining seven datasets, which included the four Hagerstown horizons where no matric potential variation was imposed, show no evident relationship between initial bulk density and water content. When Eq. [4] was used to generate a general relationship between initial soil bulk density and the initial water content for all 99 disturbed soil samples, the resulting regression was highly significant (P < 0.001) but explained less than half the variation (R2 = 0.487). The fact that only 49% of the variation is explained by Eq. [4] for the 99 disturbed samples likely results because the soils used varied in texture and organic matter content. Variability of the bulk density and the higher preconsolidation stress of the undisturbed cores masked or prevented any equivalent variation in those samples.
In general, disturbed samples with a low initial bulk density had the lowest preconsolidation stress and the highest final bulk density, whereas those starting at a higher bulk density had a higher preconsolidation stress similar to that shown in Imhoff et al. (2004) and a lower final bulk density resulting in the crossing over of curves illustrated in Fig. 1.
Figure 2
shows n plotted against m (see Table 1 for some example values) for disturbed and undisturbed data. The points (120) represent the experimental data and the smooth curve represents the equation (R2 = 0.956, P < 0.001):
 | [8] |
The inset in Fig. 2 shows more detail for the curved region of the plot. The four circled points belong to the 1Hagerstown/4sand mixture (see Table 1 for data for one sample from that mixture) and are discussed elsewhere in this article.
View larger version (16K):
[in this window]
[in a new window]
|
Fig. 2. Plot of the coefficient n versus the coefficient m (data points) from Eq. [1]. The smooth curve is the best fit (Eq. [8]) for the 120 soil samples. The inset shows greater detail for 0.02 < m < 0.22. The four circled points belong to the soil sample mixture containing one part Hagerstown and four parts sand (1 Hagerstown/4 sand) and are discussed in the text.
|
|
Effect of Initial Soil Conditions
The results of fitting Eq. [3] are tabulated in Table 2. Values are shown for all significant coefficients, the R2 values, and the mean absolute difference, and the number of data points included in each determination for the 21 datasets. As an example, consider the Bucks soil used in Fig. 1 at the 0.76- to 1.22-m depth. Six (a, c, cc, bb, d, and ff) of the 10 coefficients in Eq. [3b] and [3c] were significant (Table 2) where n1 and n2 were fixed by Eq. [8]. Initial bulk density, water content, and matric potential terms were included in the final equation that explained 99.2% of the variance. The smooth curves shown in Fig. 1, although not generated by this equation, adequately represent the output from Eq. [3] using the coefficients shown in Table 2. Comparison of the fit of Eq. [3] to the experimental data as compared with the smooth curves shown in Fig. 1 indicate subtle differences in the way the curves fit the data, except for the sample with an initial bulk density of 1.35 Mg m–3 (filled triangles) where the smooth curve from Eq. [3] passes close to but not through 3 of the 11 points where the curvature is greatest.
Comparison with an Existing Equation
Table 3 contains the regression coefficients, R2 values, and mean absolute difference obtained when Eq. [5] was used to fit selected datasets. In two cases, coefficient D was not significantly different from zero at the 5% level, and the term containing D was removed from Eq. [5] for the regression results shown.
Figure 3
shows data for three samples selected to simplify the graph from the six samples available for the Glenelg (0.00–0.30 m) dataset. These three samples, one core (undisturbed) sample with a relatively high initial bulk density and two sieved (disturbed) samples with relatively low initial bulk densities, were selected to facilitate comparisons among equations by focusing on samples with the same density and widely different water contents and, at the same time, between samples with a similar water content but widely different initial bulk density (see Table 1). The experimental data, shown as filled squares, triangles, and stars, are the same for all four parts of Fig. 3. In Fig. 3a, 3b, and 3c, coefficients used for all smooth lines were optimized for all six samples available for the Glenelg (0.00–0.30 m) dataset but are plotted only for the three selected samples.
View larger version (35K):
[in this window]
[in a new window]
|
Fig. 3. Bulk density of three Glenelg soil samples from the 0.00- to 0.30-m depth plotted as a function of the applied stress. The data are for an undisturbed core sample (stars) and for two disturbed samples at different initial water contents (triangles and squares). See Table 1 and legends for initial conditions and additional information. The smooth lines result from nonlinear fits based on (a) Eq. [5], (b) Eq. [3] with n1 and n2 as in Eq. [8], (c) Eq. [3] with optimized n1 and n2 values, and (d) Eq. [1]. Optimized values for fits in (a), (b), and (d) are in Table 3, Table 2, and Table 1, respectively, and for (c) are b = –1.43, c = –1.26, cc = 0.00121, f = 0.753, and n2 = 0.165, with all other coefficients set to zero.
|
|
|
DISCUSSION
|
|---|
Correlation of n and m
The hypothesis that the coefficient n in Eq. [1] is related to the coefficient m in Eq. [1] by the inverse relationship shown in Eq. [2] is accepted. The specific relationship is fitted by Eq. [8], is plotted in Fig. 2, and was highly significant. The substitution of Eq. [8] into Eq. [1] changes Eq. [1] from a three-coefficient (, n, and m) equation to a two-coefficient ( and m) equation for soil compression on an individual soil sample. As Assouline (2002) points out, the definition of relationships between equation coefficients and soil physical properties is easier when fewer coefficients are used in a soil compression equation. The experience of using this revised equation indicates that there are times (9 out of 21 cases in this study; see Table 2) when better results occur when n1 and n2 are used as regression coefficients rather than given universal fixed values.
Effect of Initial Soil Conditions
The hypothesis that and m are functions of the initial bulk density, water content, and matric potential of a soil sample (Eq. [3b] and [3c]) is also accepted. This decision is based on the significance of the coefficients for these three independent variables when Eq. [3b] and [3c] were used in the nonlinear regression procedure (Table 2). At least one of the initial bulk density coefficients (b, bb, e, and ee) in Eq. [3] was significant in 19 of the 21 datasets shown in Table 2. At least one of the water content coefficients (c and f) was significant in 9 of the 17 datasets shown in Table 2 where variation in water content was imposed on the samples. At least one of the matric potential coefficients (cc and ff) was significant in 12 of the 16 datasets where matric potential was introduced as an independent variable (Table 2). Coefficients representing all three independent variables were significant in 4 of the 16 datasets shown in Table 2 where all three independent variables were used in the regression. Thus, it is concluded that the initial bulk density, water content, and matric potential as represented in Eq. [3b] and [3c] were useful for describing multiple compression curves for soil samples in this dataset where these factors vary before compression.
Equation [3] fit the top four Hagerstown horizons that had no deliberate variation in water content using only the constant term (a) and the bulk density terms (b, e, and ee) for the log and m coefficients. All random variation in water content was well represented by the variation in the initial bulk density, and there were no undisturbed samples to complicate the relationships. Once a variation in water content was imposed on disturbed samples, better fits were obtained when the effect of the initial water content or the matric potential or both were retained in the regression in addition to or instead of the effect that occurred due to the initial bulk density. The Rayne horizons, the deepest Hagerstown horizon, and the Hagerstown/sand mixtures represent this situation. Addition of undisturbed core data to the mix of disturbed samples with both sets of samples at an induced initial water content and matric potential caused no difficulty for the equation but typically required one or two additional coefficients to provide the best fit. The Bucks and Glenelg horizons fell into this category.
The water content (c or f) and the matric potential (cc or ff) coefficients remained in the final regression equations in 3 out of 16 cases for the coefficient and in 1 out of 16 cases for the m coefficient (Table 2). When the water content or the matric potential coefficients were removed from the relationship shown in Table 2 for these four cases, AIC increased, indicating that the equation shown in Table 2 gave the best statistical fit. One might expect that the unique water retention relationship between water content and matric potential for each horizon or soil material would be enough to eliminate one or the other in the regressions. However, because the relationship between the water content and the matric potential is nonlinear and the water content versus log effect is linear, the matric potential versus log effect is equivalent to a nonlinear effect of water content on log . The inclusion of these two variables is likely a convenient way to capture a mixed linear/nonlinear response of water content on log (Eq. [3b]) and m (Eq. [3c]). McBride (1989) used the water content, the water content squared, and the log of the water content in their study, whereas McNabb and Boersma (1996) used one minus the degree of saturation squared. Although these relationships were used in different equations to describe compression curves, they suggest the nonlinear relationship between compression characteristics and the water content of the soil suggested previously. Because Gallipoli et al. (2003) argue that the measure of the binding between the water and the solids in the soil represented by the matric potential is an important property for the compression process while the water content determines the proportion of the soil experiencing these bonding forces, using water content and matric potential or using an equation for the water retention curve will likely become more important when fitting compression curves across soils where the water retention relationship is no longer constant.
The ability to describe disturbed and undisturbed soils using Eq. [3] is somewhat surprising. It would seem that the rock fragments contained in the undisturbed cores would cause them to respond differently than the sieved fraction from the same horizon. In considering the fact that there seems to be no decrease in the goodness of fit of the equations to the data that can be attributed to the rock fragment content, consider the Hagerstown/sand mixtures. The average bulk density of the Hagerstown (0.85- to 1.00-m depth), 4 Hagerstown/1 sand, 3 Hagerstown/2 sand, 2 Hagerstown/3 sand, and 1 Hagerstown/4 sand at an applied stress of 2972 kPa was 1.94, 1.99, 2.01, 1.87, and 1.64 Mg m–3 (see Table 1 for an example from one sample of each material) for sand contents of 110, 288, 466, 644, and 822 g kg–1 (Fritton, 2001). This indicates that the bulk density at the highest applied stress of these mixtures increases slightly as sand content increases and declines once the sand content exceeds 466 g kg–1. This is consistent with other data in the literature (Taylor and Blake, 1981; van Wesemael et al., 1995). The rock fragment contents for the Bucks and Glenelg soils (see Fritton, 2001) were in the range of values that would be expected to raise the bulk density at high applied stresses by <0.07 Mg m–3. This effect seems to have been described adequately by using the initial bulk density as an independent variable in Eq. [3] because the increase in bulk density due to the rock fragments was contained in the initial bulk density values. Thus, once the preconsolidation stress was accurately represented (by in Eq. [3b]), the undisturbed compression curves were also fit accurately.
Comparison with an Existing Equation
The McNabb and Boersma (1996) generalization of the Bailey et al. (1986) equation represents the best equation for uniaxial soil compression presently available in the literature to account for variation in the initial bulk density and in the initial water content in a single equation. The McNabb and Boersma (1996) equation was originally developed to allow experimenters to quantitatively describe samples taken at different times in the field and to help describe the effect of the spatial variability in soil samples. The results shown in Table 3 indicate that the McNabb and Boersma (1996) equation, Eq. [5], fit all these data with a lower R2 and higher mean absolute difference using a larger number of coefficients than the equivalent fits shown in Table 2 using Eq. [3], with the exception of the 1 Hagerstown/4 sand material. In this last case, the McNabb and Boersma (1996) equation fit (AIC = –131) better than Eq. [3] (AIC = –124). Thus, the equation developed in this article (Eq. [3]) provides a fit with a smaller mean absolute difference than the McNabb and Boersma (1996) fit in most cases, with the caveat that there are times when the coefficients in Eq. [3a] cannot be assigned the fixed values shown in Eq. [8].
To further elaborate on the differences between the McNabb and Boersma (1996) equation and the equation in Eq. [3], consider the example shown in Fig. 3. The McNabb and Boersma (1996) fit shown in Fig. 3a has not only a higher mean absolute difference than the fits shown in Fig. 3b, 3c, and 3d but also shows an oscillating behavior as the applied pressure increases. This oscillating behavior is most evident in the soil sample with the low initial bulk density and high initial water content (the filled triangle data). In this case, the bulk density increases with an increase in applied stress, and the slope of the curve declines before increasing again at greater applied stresses. This oscillating behavior was present in all the examples (Table 3) attempted and is most evident in disturbed samples and least evident in undisturbed samples. In all cases, this oscillating behavior is considered an undesirable artifact of the equation rather than a response to a true oscillation in the data.
In Fig. 3b, Eq. [3] with coefficients shown in Table 2 was used to generate the smooth curves. The fit to the two samples at low initial water content (the filled stars and squares) is quite good, but the fit to the high initial water content sample (the filled triangles) leaves much to be desired. As was the case in the 1 Hagerstown/4 sand material, the problem was identified in the use of fixed values for n1 and n2 in Eq. [3a]. When these values were allowed to float, the results improved (AIC = –265 versus –237) as shown in Fig. 3c. The best fit still occurs when each sample is fit individually (Fig. 3d), but these individual fits do not allow extrapolation to other initial bulk density or initial water content situations, and that is the benefit of the generalized fits available with Eq. [3].
Implications
In most cases, the coefficients shown in Table 2 cannot be interpreted as soil properties. The relative size of the coefficients is partly determined by the partial correlation among the coefficients retained for any given soil material. The coefficients are also adjusting for any lack of fit of Eq. [1] to the compression curves and to the lack of fit of the compression curves to Eq. [8]. This was especially evident in the 1 Hagerstown/4 sand dataset where the lack of fit from Eq. [8] was so severe (see the circled points in Fig. 2) that the remaining coefficients were not capable of providing an acceptable result. A similar, but less severe, lack of fit occurred in eight other horizons where n1 and n2 had to be allowed to vary to meet the three criteria (significant coefficients, correlation < 0.99±, and dependency < 0.99) for a successful regression result. Thus, it is unlikely that the coefficients presented in Table 2 can be directly related to the properties of the various soil materials, and they are presented mainly for their ability to quantitatively represent multiple compression curves for a given horizon or soil material.
|
CONCLUSIONS
|
|---|
The empirical equation (Eq. [1]) proposed by Fritton (2001) to fit the uniaxial compression of a single soil sample at one water content has been further developed (Eq. [3]) to fit soil samples at various water contents simultaneously. In the process, the original relationship between the compressed soil bulk density and the applied stress was simplified by reducing the number of coefficients required from three to two. This revised formulation allowed the remaining two coefficients to be described as functions of the initial conditions of the soil (Eq. [3b] and [3c]). The initial bulk density, the initial water content, and the initial matric potential were all useful in describing the compressed bulk density as a function of the applied stress for multiple samples from a given soil material (Table 2). The new generalized empirical equation (Eq. [3]) fit data for multiple disturbed (sieved) and undisturbed (core) samples from 21 horizons from four soils and soil/sand mixtures with an R2 0.932 (P < 0.001). The new equation had to be optimized on each data set as a whole to establish significant relationships. The new equation was shown to be superior to the best previously available equation (McNabb and Boersma, 1996). The results are consistent with the reasoning of Gallipoli et al. (2003), indicating that the amount of water and its potential are significant in modifying the compression process. In addition, the results further confirm the importance of the initial bulk density as a controlling factor in soil compaction consistent with the results of McNabb and Boersma (1993) and Imhoff et al. (2004).
|
ACKNOWLEDGMENTS
|
|---|
The Pennsylvania Department of Agriculture, the Pennsylvania Agricultural Experiment Station Project 4079, and the Penn State Fund for Research funded this work. I thank the reviewers who helped tremendously to improve this paper.
Received for publication July 26, 2005.
|
REFERENCES
|
|---|
- Assouline, S. 2002. Modeling soil compaction under uniaxial compression. Soil Sci. Soc. Am. J. 66:1784–1787.[Abstract/Free Full Text]
- Assouline, S., J. Tavares-Filho, and D. Tessier. 1997. Effect of compaction on soil physical and hydraulic properties: Experimental results and modeling. Soil Sci. Soc. Am. J. 61:390–398.[Abstract/Free Full Text]
- Bailey, A.C., C.E. Johnson, and R.L. Schafer. 1986. A model for agricultural soil compaction. J. Agric. Eng. Res. 33:257–262.[CrossRef]
- Boyland, P., J. Keiper, E. Martin, J. Novak, M. Petkovsek, S. Skiena, I. Vardi, A. Wenzlow, T. Wickham-Jones, and D. Withoff. 1992. Guide to standard Mathematica packages. Ver. 2.1. Wolfrom Research, Inc., Champaign, IL.
- Défossez, P., G. Richard, H. Boizard, and M.F. O'Sullivan. 2003. Modeling change in soil compaction due to agricultural traffic as function of soil water content. Geoderma 116:89–105.[CrossRef][Web of Science]
- Eco, R.M. 2002. Validation of an extended critical state model using the mechanical behaviour of an agricultural soil. Soil Tillage Res. 68:1–16.
- Fritton, D.D. 2001. An improved empirical equation for uniaxial soil compression for a wide range of applied stresses. Soil Sci. Soc. Am. J. 65:678–684.[Abstract/Free Full Text]
- Gallipoli, D., A. Gens, R. Sharma, and J. Vaunat. 2003. An elasto-plastic model for unsaturated soil incorporating the effects of suction and degree of saturation on mechanical behavior. Geotechnique 53:123–135.[CrossRef]
- GraphPad Software. 2005. Quick answers database, FAQ#931. Available at http://www.graphpad.com/faq/viewfaq.cfm?faq=931 (verified 7 Apr. 2006).
- Imhoff, S., A.P. Da Silva, and D. Fallow. 2004. Susceptibility to compaction, load support capacity, and soil compressibility of Hapludox. Soil Sci. Soc. Am. J. 68:17–24.[Abstract/Free Full Text]
- Kirby, J.M. 1994. Simulating soil deformation using a critical-state model: I. Laboratory tests. Eur. J. Soil Sci. 45:239–248.
- McBride, R.A. 1989. Estimation of density-moisture-stress functions from uniaxial compression of unsaturated, structured soils. Soil Tillage Res. 13:383–397.[CrossRef]
- McNabb, D.H., and L. Boersma. 1993. Evaluation of the relationship between compressibility and shear strength of Andisols. Soil Sci. Soc. Am. J. 57:923–929.[Abstract/Free Full Text]
- McNabb, D.H., and L. Boersma. 1996. Nonlinear model for compressibility of partly saturated soils. Soil Sci. Soc. Am. J. 60:333–341.[Abstract/Free Full Text]
- Peng, X., and R. Horn. 2005. Modeling soil shrinkage curve across a wide range of soil types. Soil Sci. Soc. Am. J. 69:584–592.[Abstract/Free Full Text]
- Rashid, M.A., and M.Z. Abedin. 2000. Constitutive models for triaxial consolidation and shearing of unsaturated soil. p. 137–142. In Rahardjo et al (ed.) Unsaturated Soils for Asia. Proceedings of the Asian Conference on Unsaturated Soils, Singapore, 18–19 May 2000. Balkema, Rotterdam.
- Scheinost, A.C., W. Sinowksi, and K. Auerswald. 1997. Regionalization of soil water retention curves in a highly variable soilscape. I. Developing a new pedotransfer function. Geoderma 78:129–143.[CrossRef][Web of Science]
- imnek, J., and J.W. Hopmans. 2002. Parameter optimization and nonlinear fitting. p. 139–157. In Dane, J. H. and G. C. Topp (ed.). Methods of soil analysis. Part 4. SSSA, Madison, WI.
- Taylor, D.H., and G.R. Blake. 1981. Laboratory evaluation of soil mixtures for sports turf. Soil Sci. Soc. Am. J. 45:936–940.[Abstract/Free Full Text]
- Thornley, J.H.M. 1976. Mathematical models in plant physiology. Academic Press Inc., London.
- Toll, D.G., and B.H. Ong. 2003. Critical-state coefficients for an unsaturated residual sandy clay. Geotechnique 53:93–103.
- van Genuchten, M.Th. 1980. A closed form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892–898.[Abstract/Free Full Text]
- van Wesemael, B., J. Poesen, and T. de Figueiredo. 1995. Effects of rock fragments on physical degradation of cultivated soils by rainfall. Soil Tillage Res. 33:229–250.
- Wolfram, S. 1991. Mathematica: A system for doing mathematics by computer. 2nd ed. Addison-Wesley Publishing Company, Inc., Redwood City, CA.
- Wulfsohn, D., and B.A. Adams. 2002. Elastoplastic soil mechanics: Part II. Critical-state soil mechanics for unsaturated agricultural soils. Adv. Soil Dynamics 2:71–116.
TOP
ABSTRACT
INTRODUCTION
