Published in Vadose Zone Journal 3:1300-1308 (2004)
© 2004 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
ORIGINAL RESEARCH
The Effect of Vertically Decreasing Macropore Fractions on Simulations of Non-Equilibrium Solute Transport
Nathan W. Haws* and
P. Suresh C. Rao
Department of Agronomy, Purdue University, West Lafayette, IN 47907-2054
* Corresponding author (nhaws{at}jhu.edu)
Received 7 January 2004.
 |
ABSTRACT
|
|---|
The effects of a depth-decreasing macropore fraction and the corresponding decrease in the diffusive mass-transfer coefficient for solute exchange between soil matrix and macropore regions were simulated for a hypothetical, one-dimensional, mobileimmobile, transport domain. These effects are illustrated using an expression for a Damkhöler number, and with numerical simulations with a dual-porosity module within the HYDRUS-2D code. The depth-wise change in the macroporosity was represented with an exponentially decreasing macropore fraction, based on a pore-network model. Five simulations were run, four with various steady infiltration rates (two at water saturation and two at different levels of unsaturation) and one with transient infiltration rates. Deviations between the predicted solute breakthrough curves (BTCs) for simulations using depth-variable parameters and those using spatially averaged parameters were minor at the lowest infiltration rate, indicating little effect of non-equilibrium conditions. However, at higher infiltration rates, the BTC peak height and peak timing varied by as much as a factor of two, indicating significant influence of the variability of the macropore fractions with depth. Using depth-averaged values to represent the variable macropore fraction will closely approximate the BTC of a spatially variable case for both steady-state and transient infiltration; however, the concentration profile within the transport domain will not be accurately predicted. A sensitivity analysis conducted for the saturated hydraulic conductivity, length of the soil ped, and dispersivity parameters shows that changes in the prediction error between the spatially variable and surface-valued simulations are most sensitive to these parameters as the soil approaches water saturation.
Abbreviations: BTC, breakthrough curve
 |
INTRODUCTION
|
|---|
PORE SIZE AND GEOMETRY significantly influence water flow and solute transport in the vadose zone. This is particularly true of soils with a network of structural macropores that potentially lead to rapid solute breakthrough and a high degree of physical non-equilibrium (e.g., Brusseau and Rao, 1990; Vervoort et al., 1999; Ersahin et al., 2002). Non-equilibrium phenomena are more apparent as the soil approaches water saturation and participation of macropores in the flow and transport processes increases. When modeling transport through macroporous soils, the macropore fraction is frequently assumed to be spatially constant. However, a decreasing macropore abundance and connectivity with depth may be more realistic of some field scenarios. One likely application is agricultural soils of the U.S. Midwest, which contain macropore networks comprised of biochannels (e.g., roots and worm holes) that diminish with depth. Likewise, the abundance of inter-ped voids and cracks, which also contribute to preferential transport, also decrease with depth as the soil structure changes from granular, to blocky, to prismatic (USDA, 1998).
The spatial organization of preferential-flow networks affects the degree of transport non-equilibrium. When the fraction of preferential flow networks is spatially homogeneous across a scale of observation, solute advection and matrix diffusion will, in theory, linearly approach equilibrium with increasing transport distance. However, when the macropore fraction decreases with depth due to convergence of preferential pathways, the approach to equilibrium with increasing distance is impeded in at least two ways. First, the decreasing macropore fraction results in increasing pore-water velocities. Second, increasing spacing between macropores leads to increasing times for inter-ped solute diffusion.
Deurer et al. (2003) formally quantified the depth variability of preferential networks. Assuming that values of soil hydraulic properties could indicate connectivity of solute and water pathways, they used a geostatistical approach to reconstruct the spatial patterns of preferential networks. This approach showed that preferential drainage networks exhibit a consistent dendritic pattern with depth and that the macropore fraction decreased exponentially with depth. They then constructed a simple pore network model to demonstrate that during active macropore flow the decreasing transport volume causes water and solutes to move at increasing velocities.
The influence of a depth-decreasing macropore network is investigated here using the pore-network model of Deurer et al. (2003) to portray a hypotheticalthough realisticmacroporous transport domain. While a pore network model cannot fully represent complicated macropore sequences, such models have been shown to give simple, yet fairly accurate, depictions of the routes governing flow and transport through soils (Vogel and Roth, 2003). It is assumed here that the pore-network conceptualization of Deurer et al. (2003) satisfactorily characterizes the general behavior of non-equilibrium transport through a field soil. This allows a simple network model of a macropore network to be coupled with a flow and transport code (HYDRUS-2D; Simunek et al., 1999) for the purpose of evaluating the hypothesis that a depth-decreasing macropore network affects the degree of non-equilibrium conditions for solute transport. Because solute transport through macropore networks is highly dependent on the degree of water saturation, both (i) saturated, steady flow and (ii) unsaturated, steady and transient flow conditions are explored here. The implications of these findings can be used as a guide for future considerations of how flow and transport parameters are measured and modeled.
 |
THEORY AND METHODS
|
|---|
Conceptual Modeling of the Macropore Network and Physical Non-Equilibrium
The Damkhöler number,
, is a nondimensional parameter commonly used to assess chemical non-equilibrium. When the rate-limiting transfer process is diffusion into immobile regions, rather than soilwater partitioning, the Damkhöler number can be redefined as a ratio of the characteristic advection time, Ta*, to the characteristic diffusion time, Td*:
 | [1] |
where z (L) is the depth below the surface, Da (L2/T) is the effective solute diffusion coefficient through a given soil matrix, v (L/T) is the advective solute velocity through the soil macropores, a (L) is the effective path length for diffusive transfer between the soil matrix and macropore regions, and the symbol E (unitless) is used here to distinguish the physical non-equilibrium formulation of the Damkhöler number from its usual definition for chemical non-equilibrium. Similar to
, smaller values of E indicate a greater relative degree of non-equilibrium transport.
Assuming that advective water flow occurs only in the macropores (i.e., mobileimmobile concept), the mean transport velocity, v, for a specified flux density across the entire mobileimmobile domain, q (L/T), at a vertical location, z (L), is computed as:
 | [2] |
where
f(z) (L3/L3) is the depth-variable water content in the macropores. Approximating the subangular-blocky structure of Midwestern soils (USDA, 1998) as a rectangular (cubic) soil ped, the mean diffusion path length, a, can also be expressed as a function of
f(z) using:
 | [3] |
where
t(z) (L3/L3) is the total volumetric water content for the entire soil matrixmacropore region, rt (L) is the total length of the soil ped, and the overbar indicates a depth-averaged quantity. When the soil is water saturated, the macropore water content equals the macroporosity and can be directly related to macropore network fraction at a given depth using the model of Deurer et al. (2003):
 | [4] |
where the subscript s denotes that the parameter is at water saturation. The parameter Nf(z) (unitless) is the macropore network fraction that exponentially decreases with depth according to:
 | [5] |
where c (unitless) is the asymptotic value of the macropore network fraction, d (unitless) is a fitting parameter, zo (L) is the starting depth for the network, and L (L) is the characteristic length of the network (Deurer et al., 2003). Taking depth-averaged values of all spatially variable quantities and substituting Eq. [25] into Eq. [1] gives the Damkhöler number as a function of depth:
 | [6] |
The pore-network model proposed by Deurer et al. (2003) is adopted here as only one approximation of a macropore network. Other representations of depth variations are possible and could be similarly investigated.
Numerical Simulations of a Spatially Variable Macropore Network
Numerical simulations of the hypothetical macropore network were performed using the HYDRUS-2D flow and transport model (Simunek et al., 1999). Although HYDRUS-2D is a two-dimensional, finite-element code, the modeled transport domain only had a width of one element, so that all transport was one-dimensional in the vertical downward (z) direction. The version of HYDRUS-2D used in these simulations had been modified from the standard HYDRUS-2D code to include a mobileimmobile module to represent the macropore and soil matrix regions (Simunek et al., 2003). The mobileimmobile model assumes that advective water flow occurs only in the mobile domain, yet allows both advective water movement and solute diffusion between the mobile (macropore) and immobile (matrix) regions. Intra-ped transfer for advection and diffusion are represented as first-order processes driven by the pressure and concentration gradients, respectively. The equations for water flow and solute transport for the mobileimmobile model are (Simunek et al., 2003):
- Water:
 | [7a] |
 | [7b] |
 | [7c] |
- Solute:
 | [8a] |
 | [8b] |
 | [8c] |
where
(L3/L3) is a volumetric water content (the subscript f denotes the macropore, or fracture, domain and the subscript m the matrix domain), h (L) is the matric potential, K (L/T) is the hydraulic conductivity (a function of the matrix potential in the macropore region), C (M/L3) is the aqueous solute concentration, R (unitless) is the retardation factor, D (L2/T) is the hydrodynamic dispersion coefficient, and
f (unitless) is the ratio of the macropore water content to the total water content (
f/
t). The terms
w (1/T) and
s (M/L3T) are the water and solute inter-region transfers terms, respectively, where C* = Cf for
w > 0 and Cm for
w > 0 (Simunek et al., 2003). The parameters
w (1/[LT]) and
d (1/T) modifying the water and solute transfer terms are first-order transfer coefficients for advective and diffusive exchange between the mobile and immobile regions. In HYDRUS-2D, the hydraulic conductivity, K, and the unsaturated macropore water content,
f, are computed from the saturated hydraulic conductivity, Ks, and the macroporosity,
fs, based on the
and n parameters in the van Genuchten (1980) pressure-saturation equation.
For these simulations, the diffusive transfer was assumed to be responsible for all intra-ped solute exchange, and the water transfer term,
w, was neglected. This assumption is strictly valid for the simulations that were conducted under saturated and/or unsaturated steady-state water flow conditions because the pressure heads in the matrix and macropore regions are at equilibrium. For the transient simulations, water transfer between the two domains during the wetting and draining periods is possible. However, to better highlight the effect of the diffusive transfer on transport non-equilibrium, the advective transfer term was again not considered.
The model system was a 200-cm soil profile. This length, as well as the specific soil hydraulic parameters listed in Table 1, were selected to reflect measured and calibrated properties of cropland [corn, Zea mays L., and soybean, Glycine max (L.) Merr.] soils at the Purdue University Agronomy Farm. The total porosity,
t, was assumed constant, and the depth decrease in the macropore fraction and the increase in the diffusion path length were step-wise approximated using 10 soil layers, with the thickness of each layer set to correspond to a 10% change in the macropore fraction. The macropore fraction was modeled using Eq. [4] and [5] with parameters selected so that
fs varied from 0.1 at the surface (z = 0 cm) to 0.01 at L (z = 200 cm). This surface value corresponded to the calibrated values for intact soil cores taken from the Purdue Agronomy Farm, and the value at z = L was estimated based on a two-dimensional modeling study of solute transport in cornfield plots at the farm. The value of rt (2.35 cm) was chosen so that the diffusion path length, a, at the surface, computed using Eq. [3], would match that estimated for an intact soil core and several soil aggregates obtained from the same site (Haws et al., 2004).
Because HYDRUS-2D represents intra-ped solute diffusion with a first-order diffusion term,
d, rather than a true Fickian diffusion process, the values of an apparent first-order term,
dapp, were estimated from the value of Da using the relationship:
 | [9] |
where ß is a dimensionless shape factor that accounts for the influence of soil geometry on mass transfer (e.g., van Genuchten and Dalton, 1986). The exponential decrease in the first-order diffusion coefficient corresponding to the change in diffusion path length was from 0.02 h1 (z = 0) to 0.013 h1 (z = L).
A series of simulations was conducted once for the layered model. Then, to investigate the differences in solute flux prediction using depth-decreasing versus spatially constant parameters, the same simulations were run again with a constant macropore network fraction throughout the profile based on the surface-layer parameters (
fs = 0.1,
d = 0.02 h1). The series of simulations was also run a third time using as the constant parameters the depth-averaged values of the layered simulation (
fs = 0.043,
d = 0.015 h1).
The simulation series consisted of five model runs. The first four simulations considered the effects of the degree of water saturation and flow-rate on non-equilibrium transport. These simulations were conducted under steady water flow and constant head boundary conditions. Three of the four simulations had a unit gradient with boundary pressure heads of h = 50, 5, and 0 cm, giving flow rates, q* = q/Ks, of 0.0066, 0.24, and 1.0, respectively. The fourth simulation had a vertical gradient of 1.25 (htop = +50 cm and hbottom = 0), giving a flow rate of q* = 1.25. These four boundary conditions represent two fully saturated conditions and two unsaturated flow conditions and cover a range from essentially no macropore flow (h = 50 cm boundary) to a ponded domain with a fully active macropore fraction (htop = +50 cm). At the beginning of each simulation, a 0.1-PV solute pulse was applied at the top of the column. Breakthrough curves (BTCs) were then generated at relative distances, z/L, of 0.25, 0.5, 0.75, and 1.0.
The fifth simulation had a variable-flux boundary at the top and a free drainage boundary at the bottom. This simulation had an initial condition of hydrostatic equilibrium and was run with a sequence of four infiltration events (see Fig. 1)
. The solute pulse was applied during the first infiltration event, a period of 7 h. The cumulative solute flux at z/L = 1.0 was used to evaluate the layered and spatially constant model predictions.
Test of Parameter Sensitivity
A test of the sensitivity of the discrepancy between the layered and surface-valued simulations with respect to Ks, rt, and
(dispersivity) was also performed. In turn, one test parameter was increased, then decreased, by 50% from its reference value (value used in the original simulations) while leaving the other test parameters at their respective reference values. For each test combination, simulations were run for the four different boundary conditions for the layered and surface-value cases. Then, the root mean square error (RMSE) between the layered simulation and the corresponding surface-value simulation was evaluated as an indication of the difference between the two cases. The correlation coefficient (R2) was also computed and used as a relative indication of the degree of correlation between the layered and surface-valued cases. For consistency among the various simulations, the RMSE and R2 calculations were based on equal pore volume increments for a period of 3 pore volumes. Though arbitrarily established for convenience in modeling, the 3 pore volume period allowed for between 97 and 100% mass recovery, and extending recovery periods beyond this point had only minor effects on the RMSE and R2 values.
 |
RESULTS
|
|---|
Damkhöler Number
Figure 2
shows the Damkhöler number, E, plotted as a function of depth for the four boundary conditions of the steady flow numerical simulations using the parameters listed in Table 1. The overall value of E decreases as the flow rates increase, and for the higher flow rates (q* = 0.24, 1.0, and 1.26), the value of E remains well below 1 at all depths. For the non-equilibrium transport parameters derived from the network model, the increases in E with depth are initially rapid near the surface, but quickly level at greater depths. This leveling off is contrasted with the linear rise (with z) in E values based on constant non-equilibrium transport parameters (Fig. 2). The E values based on the surface parameters (Fig. 2) are greater than the network model's E values at all depths except the surface depth, and it is more than twice that of the E values of the network model at the bottom of the transport domain. On the other hand, the E values derived from depth-averaged parameters underpredict the network model's E values at all depths except at the end of the transport domain, where the values match. The underprediction of E values based on surface parameters compared with the E values based on the network model parameters results from not considering the increase in diffusion path length and flow velocities in the macropores with depth. This is further illustrated in the numerical simulations.

View larger version (26K):
[in this window]
[in a new window]
|
Fig. 2. Damkhöler number, E, as a function of depth for the q* = 0.0066 and 1.25 boundary conditions used in the steady flow numerical simulations. Also shown are the E values for the spatially constant representation using surface values (surface) and depth-averaged values (depth-avg).
|
|
Steady Flow Numerical Simulations
Layered vs. Surface Simulation
When the layered network model simulation (labeled "layered" in Fig. 3)
is compared with the simulation with constant values based on the surface parameters (labeled "surface" in Fig. 3), the discrepancy in the BTCs at the domain outlet (z/L = 1) increases with increasing infiltration rate (Fig. 3). For the slowest flow rate (q* = 0.0066), the layered and surface BTCs are nearly identical with neither BTC manifesting significant non-equilibrium transport. However, as the flow rates increases, differences between the layered and surface BTCs amplify, with the layered case exhibiting much greater non-equilibrium transport (i.e., earlier breakthrough and longer tailing). When infiltration intensities approach or exceed Ks (q* = 0.24, 1, 1.25) the BTC peak arrival times at z/L = 1 may be overpredicted and peak heights underpredicted by a factor of two or greater in comparison with the spatially constant network fraction (Fig. 3). These observations coincide with the computed Damkhöler number values of E > 1 for q* = 0.0066 and E << 1 for the other flow rates (Fig. 2).
Layered vs. Depth-Averaged Simulations
The success of the simulations using depth-averaged values in reproducing the layered cases is dependent on both the depth where the BTC is observed and the infiltration rate. To illustrate the ability of depth-averaged values to represent the layered simulations at different depths, Fig. 4
shows the BTCs for the q* = 0.24 case at depths of z/L = 0.25, 0.50, and 1.0. The depth-averaged BTC matches the BTC of the layered case only at z/L = 1.0. This result was predicted by the matching E values at z/L = 1.0 in Fig. 2, and the small discrepancies between the simulated BTCs at this depth are attributed to the layered approximation of the continuous exponential network model. Despite the matching BTC at the end of the transport domain, the depth-averaged simulations are unable to replicate the layered BTC at relative transport distances less than z/L = 1, with discrepancies increasing as the relative depths decrease.

View larger version (24K):
[in this window]
[in a new window]
|
Fig. 4. Comparison of breakthrough curves between layered (solid) and depth-averaged (dashed) values at depths of z/L = 0.25, 0.50, and 1.0 for the q* = 0.24 case.
|
|
Transient Flow Simulations
Similar to the steady flow simulations, the depth-averaged transient flow simulation matches well with the layered transient flow simulation for the solute flux at z/L = 1 (Fig. 5)
, such that any discrepancies are attributed to the layered representation of a continuous preferential-flow network. The surface simulation, however, appreciably underestimates the solute flux at the domain outlet (z = L).
Sensitivity of Layered vs. Surface Simulations to Changes in Parameter Values
The RSME and R2 values (Table 2) highlight the relationship between parameter sensitivity and the boundary pressure head (or q* value). As expected, RMSE values are lowest and R2 values greater at more negative boundary heads, and quickly decrease with low boundary pressures (i.e., as the macropore network becomes more active). At the 50 cm boundary head, the flow rates are low enougheven for the higher Ksthat both the layered and surface-value simulation are at, or near, transport equilibrium. Consequently, changes in the discrepancies between the layered and surface-valued cases are relatively insensitive to variations in parameter values at the 50 cm boundary pressure (q* = 0.0066) and R2 values are greater than 0.99 for all test values. Additional simulations run with a Ks of 50 cm/h (not reported in Table 2) likewise gave R2 values greater than 0.99. An insensitivity of changes in Ks to changes in the RMSE and R2 is also manifest as the soil becomes saturated and the macropore network is fully active. For the boundary cases of 0 cm (q* = 1.0) and +50 cm (q* = 1.25) the RMSE changes by less than a factor of two from the +50% to 50% cases and there is virtually no change in the R2 values. The critical region for sensitivity to the value of Ks is at the slightly negative boundary heads (i.e., 5 cm, or q* = 0.24) where the value of Ks can dictate flow velocities and non-equilibrium responses that create the largest discrepancies between the layered and surface-value simulations. In this case, the RMSE increases by over a factor of three while the R2 value shows a twofold reduction from the 2.5 to 7.5 cm/h simulations.
View this table:
[in this window]
[in a new window]
|
Table 2. Root mean square error (RSME) and correlation coefficient (R2; in parentheses) between layered and surface-value simulations at the four tested boundary pressures (or q* values) and for variable saturated hydraulic conductivity (Ks), matrix diffusion radius (rt), and dispersivity ( ).
|
|
Variations in the rt parameter values show similar RMSE sensitivity characteristics to changes in the Ks values. This is somewhat unexpected because in Eq. [6] the rt parameter is a squared term, and Ks is linear. Another nonintuitive response is that the R2 values actually increase for both the rt = 1.18 and 3.52 cm with respect to the reference case (rt = 2.35 cm). The increase in R2 values for the rt = 1.18 cm simulations can be attributed to a more rapid mass transfer rate, and thus a decrease in the degree of non-equilibrium transport. The increase in R2 values for the rt = 3.52 cm cases results because the increasing velocities and early breakthrough times are primarily controlled by the macropore fraction,
f, which, as represented in Eq. [4] and [5], is not a function of rt. As a result, early breakthrough characteristics are similar for the reference case and the other rt cases; however, the slower mass transfer rates of the rt = 3.52 cm case causes a dampening of the peak breakthrough concentration of the layered case; thus, an overall increase in the R2 values (see Fig. 6)
.

View larger version (26K):
[in this window]
[in a new window]
|
Fig. 6. Results of sensitivity analysis for the h = 5 cm (q* = 0.24) case. The top graph shows the reference simulation (used in the original simulations). The graphs underneath show the layered and surface-value simulations after varying a specific test parameter (noted on right) down by 50% of its reference value (graphs on left) or up by 50% of its reference value (graphs on right).
|
|
The simulations conducted using different values of the dispersivity,
, reveal that this parameter can dramatically affect the correlation of the surface simulations with the layered simulations. While it is not surprising that the
parameter would influence the correlation, its importance might be overlooked as it does not explicitly appear in equation for the Damkhöler number (Eq. [6]). Previous studies (e.g., Brusseau and Rao, 1990) have demonstrated that larger values of
can mask some of the effects of an underestimation in
d. The simulations here similarly show that correlations, though not necessarily RMSE, between the layered and surface-value simulations are highly sensitive to values of
(Table 2, Fig. 6).
 |
DISCUSSION AND CONCLUSIONS
|
|---|
An exponentially decreasing pore-network fraction increases the manifestation of non-equilibrium transport compared with predictions that assume the surface macropore fraction remains constant through the soil profile. Such effects are negligible at infiltration rates several times smaller than the saturated hydraulic conductivity, Ks, where equilibrium is approached (E
1), and near the upper boundary of the transport domain where the reduction in the network fraction is minimal. However, when infiltration intensities approach or exceed Ks, and transport distances are long enough so that there is a significant decrease in the network density, the BTC peak arrival times may be overpredicted and peak heights underpredicted by a factor of two or greater in comparison with the spatially constant network fraction.
Constant network fraction parameters derived from depth-averaged values of the exponential network model exactly reproduce the solute response at z/L = 1 for both the steady-state and transient simulations. However, while the flux response at the domain outlet is captured when using the depth-averaged values, the solute fluxes interior to the transport domain are overpredicted (as opposed to the underpredictions of the non-equilibrium in the constant parameter model based on surface values), and these overpredictions increase as the actual non-equilibrium transport of the transport scenario increases (i.e., smaller E values influenced by higher infiltration rates and smaller transport distances). Thus, predictions of solute transport based on depth-averaged values are useful for predicting temporal solute fluxes at the domain outlet, but should not be used to describe the depth distribution of solute concentrations.
The importance of considering dual-porosity (mobileimmobile) transport processes, particularly during periods of high infiltration rate when preferential flow effects are most pronounced, has been demonstrated in several field studies (e.g., Richard and Steenhuis, 1988; Kung et al., 2000; Fortin et al., 2002). Resourceful techniques have been developed to estimate the preferential flow and transport parameters (
f and
d) using tension infiltrometers with either a single tracer chemical (Clothier et al., 1992, 1995), multiple tracers (Jaynes et al., 1995; Casey et al., 1997), and simultaneously at multiple surface locations (Al-Jabri et al., 2002). However, while the horizontal variability of the transport processes in field soils might be quantified, most field measurement techniques and solute transport models neglect the effects of changes in
f and
d with depth. For systems of a small enough scale so that the development and depth variability of macropore networks is minimal, the spatial dependence in non-equilibrium transport may be unimportant. However, in larger systems containing complex fracture or macropore networks, first-order diffusion terms that are assumed spatially constant at the surface-estimated value may appreciably underestimate the actual degree of transport non-equilibrium. The sensitivity to properly delineated soil and transport parameters and their spatial variation becomes most pronounced as the soil approaches saturation sufficiently to activate macropores.
This simulation study of a hypothetical soil with an exponentially decreasing macropore fraction has illustrated some of the likely effects of vertical variability in the soil structure. An underlying assumption of this investigation was that the mobileimmobile concept with first-order mass transfer was an adequate representation of a soil macroporematrix transport domain. The importance of the findings presented here are contingent on the validity of this assumption. Other approximations of macroporous transport domains exist (e.g., Simunek et al., 2003) and could influence the relative importance of the results of this study. In addition, the significance of the spatial macropore structure should also be investigated for soils exhibiting different depth-wise changes in their soil transport properties than characterized by the exponential model of Deurer et al. (2003).
 |
REFERENCES
|
|---|
- Al-Jabri, S.A., R. Horton, D.B. Jaynes, and A. Gaur. 2002. Field determination of soil hydraulic and chemical transport properties. Soil Sci. 167:353368.
- Brusseau, M.L., and P.S.C. Rao. 1990. Modeling solute transport in structured soils: A review. Geoderma 46:169192.[ISI]
- Casey, F.X.M., S.D. Logsdon, R. Horton, and D.B. Jaynes. 1997. Immobile water content and mass exchange of a field soil. Soil Sci. Soc. Am. J. 61:10301036.[Abstract/Free Full Text]
- Clothier, B.E., L. Heng, G.N. Magesan, and I. Vogeler. 1995. The measured mobile-water content of a soil as a function of hydraulic regime. Aust. J. Soil Res. 33:397414.
- Clothier, B.E., M.B. Kirkam, and J.E. McLean. 1992. In situ measurements of the effective transport volume for solute moving through soil. Soil Sci. Soc. Am. J. 56:733736.[Abstract/Free Full Text]
- Deurer, M., S.R. Green, B.E. Clothier, J. Bottcher, and W.H.M. Duijnisveld. 2003. Drainage networks in soils: A concept to describe bypass-flow pathways. J. Hydrol. (Amsterdam) 272:148162.
- Ersahin, S., R.I. Papendick, J.L. Smith, C.K. Keller, and V.S. Manoranjan. 2002. Macropore transport of bromide as influenced by soil structure differences. Geoderma 108:207223.[ISI]
- Fortin, J., E. Gagnon-Bertrand, L. Vezina, and M.G. Rompre. 2002. Preferential bromide and pesticide movement to tile drains under different cropping patterns. J. Environ. Qual. 31:19401952.[Abstract/Free Full Text]
- Haws, N.W., B.S. Das, and P.S.C. Rao. 2004. Dual-domain solute transfer and transport processes: Evaluation in batch and transport experiments. J. Contam. Hydrol. (in press).
- Jaynes, D.B., S.D. Logsdon, and R. Horton. 1995. Field method for measuring mobile/immobile water content and solute transfer rate. Soil Sci. Soc. Am. J. 59:352356.[Abstract/Free Full Text]
- Kung, K.-J.S., E.J. Kladivko, T.J. Gish, T.S. Steenhuis, G. Bubenzer, and C.S. Heling. 2000. Quantifying preferential flow by breakthrough of sequentially applied tracers: Silt loam soil. Soil Sci. Soc. Am. J. 64:12961304.[Abstract/Free Full Text]
- Richard, T.L., and T.S. Steenhuis. 1988. Tile drain sampling of preferential flow on a field scale. J. Contam. Hydrol. 3:307325.
- Simunek, J., N.J. Jarvis, M.T. van Genuchten, and A. Gardenas. 2003. Review and comparison of models for describing nonequilibrium and preferential flow and transport in the vadose zone. J. Hydrol. (Amsterdam) 272:1435.
- Simunek, J., M. Sejna, and M.Th. Van Genuchten. 1999. The HYDRUS-2D software package for simulating the two-dimensional movement of water, heat, and multiple solutes in variably-saturated media: Version 2.0. U.S. Salinity Lab., Riverside, CA.
- USDA. 1998. Soil survey of Tippecanoe County, Indiana. Purdue University Agric. Exp. Stn., Indiana Dep. of Nat. Resour., State Soil Conservation Board, and the Div. of Soil Conserv., West Lafayette.
- van Genuchten, M.Th. 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892898.[Abstract/Free Full Text]
- van Genuchten, M.Th., and F.N. Dalton. 1986. Model for simulating salt movement in aggregated field soil. Geoderma 38:165183.[ISI]
- Vervoort, R.W., D.E. Radcliffe, and L.T. West. 1999. Soil structure development and preferential solute flow. Water Resour. Res. 35:913928.
- Vogel, H.-J., and K. Roth. 2003. Quantitative morphology and network representation of soil pore structure. Adv. Water Resour. 24:233242.
This article has been cited by other articles:

|
 |

|
 |
 
T. B. Ramos, M. C. Goncalves, J. C. Martins, M. Th. van Genuchten, and F. P. Pires
Estimation of Soil Hydraulic Properties from Numerical Inversion of Tension Disk Infiltrometer Data
Vadose Zone J.,
May 26, 2006;
5(2):
684 - 696.
[Abstract]
[Full Text]
[PDF]
|
 |
|