Published online 20 November 2007
Published in Vadose Zone J 6:890-898 (2007)
DOI: 10.2136/vzj2007.0124
© 2007 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
ORIGINAL RESEARCH
Analytical Advection–Dispersion Model for Transport and Plant Uptake of Contaminants in the Root Zone
T. H. Skaggsa,*,
N. J. Jarvisb,
E. M. Pontedeiroc,
M. Th. van Genuchtena and
R. M. Cottad
a USDA-ARS, U.S. Salinity Lab., 450 W. Big Springs Rd., Riverside, CA 92507
b Swedish Univ. of Agricultural Sciences, Uppsala, Sweden
c Brazilian Nuclear Energy Commission, Rio de Janeiro, Brazil
d Federal Univ. of Rio de Janeiro, Rio de Janeiro, Brazil
* Corresponding author (tskaggs{at}ussl.ars.usda.gov).
All rights reserved. No part of this periodical may be reproduced or transmitted in any form or by any means, electronic or mechanical, including photocopying, recording, or any information storage and retrieval system, without permission in writing from the publisher.
Received 5 July 2007.
 |
ABSTRACT
|
|---|
In regulatory and risk management analyses of environmental contaminants, the vadose zone may be treated as a subcomponent within a larger environmental modeling framework. For the complexity of the larger system model to remain at manageable levels, it is desirable that subcomponent models be relatively simple and require few input parameters. In this work, we develop an advective–dispersive solute transport equation that includes plant uptake of water and solute and present an analytical solution. Assumptions underlying the transport model include linear solute sorption, first-order uptake, and a uniform soil water content. We examine the latter assumption in detail and demonstrate the effects of rooting depth, soil texture, and leaching fraction on the uniformity of the root-zone water content. The new analytical advection–dispersion model should be useful for estimating the transport and uptake of strongly sorbing and persistent contaminants, where the timescale relevant for assessing environmental impacts is long (decades) and short-term fluctuations caused by, for example, precipitation can be averaged. As an illustration, model predictions are made for the uptake of cadmium (Cd) by wheat (Triticum aestivum L.) grown in sludge-amended soil. The predictions are compared with those of a "one-compartment" model that has been proposed previously for risk analysis and regulatory studies. The comparison shows that the one-compartment model overestimates the long-term, steady-state Cd concentration in harvested wheat grain. The analytical advection–dispersion model is recommended as a tool for environmental risk assessment of strongly sorbing, persistent contaminants.
Abbreviations: ADE, advection–dispersion equation • GITT, generalized integral transform technique
|
INTRODUCTION
|
|---|
The movement of solutes in the root zone involves complex interactions of physical, chemical, and biological processes. Understanding and modeling these processes, as well as their impacts on plant growth and water and soil quality, remains a challenging area of research for scientists and engineers. Numerical computer simulation models such as HYDRUS (
im
nek et al., 2005), UNSATCHEM (imnek et al., 1996), MACRO (Larsbo et al., 2005), and RZWQM (Ahuja et al., 2000) permit comprehensive analyses of root zone processes by allowing for a wide range of soil properties, chemical reactions, and root functions. Although these models provide great flexibility, their complexity has a downside. Using the models requires specifying values for a large number of parameters, many of which are often unknown or unobtainable. Long-term simulations, such as may be required when evaluating the transport and plant uptake of strongly sorbing (i.e., slow-moving) contaminants such as trace metals and radionuclides, are difficult because of the need to specify time-varying boundary conditions over many years, and because of lengthy execution times. Also, the complexity of numerical codes is such that many of the details of the models are often not apparent to users other than the model developers, thus making it difficult for both technical experts and regulatory agencies to interpret simulation results.
Alternatively, solute transport and plant uptake in the root zone can be studied using less-complex models, which typically take the form of a partial differential equation that can be solved analytically. While analytical models have the advantage of being relatively transparent with respect to model inputs and outputs, this simplicity often comes at the expense of significant assumptions and approximations, such as requiring uniform soil conditions and/or representing nonlinear chemical reactions with linear models. In some instances, however, these simplifications may not be too restrictive. For example, in long-term simulations (i.e., over tens or hundreds of years), short-term boundary fluctuations caused by precipitation can be averaged and an analytical model with time-averaged boundary conditions and root-zone average soil and transport properties may closely approximate a more complex, transient model (e.g., Destouni, 1991). Such an approximation may be poor in the case of quickly degrading solutes (e.g., pesticides) but should work well for strongly sorbing and persistent contaminants (e.g., trace metals, radionuclides) when the relevant timescale for assessing environmental impacts is measured in decades or centuries rather than seasons and years.
Among analytically solvable solute transport models, abundant literature exists on the classical advection–dispersion equation (ADE). While analytical solutions to the ADE are widely available (e.g., van Genuchten and Alves, 1982; Skaggs and Leij, 2002), most do not explicitly incorporate water and solute uptake by roots, thus limiting their usefulness in analyzing root zone processes. Raats (1975) used the method of characteristics to obtain an analytical transport model that included root water uptake but neglected diffusion and dispersion processes. Generalizations and extensions of the Raats (1975) model have been provided by Ginn and Murphy (1997) and Schoups and Hopmans (2002). The latter study included plant uptake of both water and solute. Hoffman and van Genuchten (1983) reviewed various steady-state transport solutions that include root water uptake while neglecting diffusion and dispersion, whereas Raats (1981) considered the effect of a constant solute diffusion coefficient on the shape of steady solute profiles. To date, none of the analytical modeling efforts have included a general (i.e., possibly nonconstant) description of diffusion and dispersion. In this work, we develop a solute transport model that accounts for advection and dispersion diffusion, as well as for the uptake of water and solutes by roots, and we present an analytical solution. Additionally, we discuss general characteristics of the proposed transport model and consider an example application involving the uptake of cadmium by wheat (Triticum aestivum L.) growing in sludge-amended soil.
|
Theory
|
|---|
Steady Water Flow with Plant Uptake
Conservation of mass for one-dimensional steady water flow can be written as
 | [1] |
where q is the water flux (L3 L–2 T–1 ), rw(z) is a sink term specifying the rate of water extraction by roots (L3 L–3 T–1 ), and z is the vertical coordinate (L), oriented positive downward with the soil surface located at z = 0. The sink term can be expressed as (e.g., Feddes et al., 1978; Skaggs et al., 2006)
 | [2] |
where T is the transpiration rate (L3 L–2 T–1) and b(z) is the normalized uptake density (L–1). Integration of Eq. [1] gives the steady-state water flux profile
 | [3] |
where q0 
q(0) and B(z) is the cumulative uptake distribution,
 | [4] |
We are interested only in the case of q0 > T, such that water flow is vertically downward (q > 0).
Advective–Dispersive Solute Transport and Plant Uptake
Conservation of mass for one-dimensional solute transport can be expressed as
 | [5] |
where C is the solute concentration (M L–3), 
is the volumetric water content (L3 L–3), 
b is the soil bulk density (M L–3), S is the sorbed solute concentration (M M–1), js is the solute mass flux (M L–2 T–1), and rs is a solute sink term (M L–3 T–1) that accounts for the uptake of solute by roots. The solute flux is assumed to follow the standard advection-dispersion model,
 | [6] |
where D is the effective diffusion–dispersion coefficient (L2 T–1), which may vary with depth due to the dependence of dispersion on solute velocity. In this study, we assume that solute uptake is a first-order process such that
 | [7] |
where 
is the uptake coefficient ( 
0). A mechanistic interpretation of in terms of transport across root membranes is provided by Dalton et al. (1975); additional discussion of the first-order uptake model in given below in "First-Order Plant Uptake." We also assume that the soil water content is approximately uniform with depth [(z) 
] and that sorption is linear with a partitioning coefficient Kd. With these assumptions, introducing Eq. [2], [3], [6], and [7] into Eq. [5] leads to the following transport equation:
 | [8] |
where R is the retardation coefficient, v0 q0/, and vT T/. Note that because q0 > T (see "Steady Water Flow with Plant Uptake"), it follows that v0 > vT. Formally, R is defined by R = 1 + Kdb/, although in practice R is commonly treated as an effective, empirical parameter. The assumption that the water content is uniform with depth is important to obtaining Eq. [8]; we examine that approximation in detail in the next section.
The left-hand side of Eq. [8] and the first term on the right-hand side are familiar terms common to many advection–dispersion transport models (e.g., Skaggs and Leij, 2002). The second term on the right-hand side is an advection term that contains a depth-dependent transport velocity, which decreases with depth from a maximum value of v0 at the soil surface to a minimum value of v0 – vT at the bottom of the root zone. The final term in Eq. [8] is a combined sink–source term that accounts for the concentrating effect of a decreasing water flux and the countering effect of solute removal. When = 0, no solute is taken up and the concentrating effect in the root zone is maximal. When solute is taken up passively with water (0 < 
1), the concentrating effect is reduced, becoming nonexistent when = 1. Specifying > 1 corresponds to active solute uptake (Schoups and Hopmans, 2002), such that solute concentrations are reduced as a result of plant uptake.
Uniformity of Soil Water Content
The derivation of the solute transport equation in the previous section involved the assumption that the soil water content is approximately uniform with depth. Previous theoretical analyses of water flow with root uptake have found that the pressure head depth-profile is generally very uniform during steady flow (Raats, 1974) or quasi-steady flow such as may be achieved with high frequency irrigation (Rawlins, 1973). Steady-state water content profiles computed by Yuan and Lu (2005, their Fig. 4) similarly showed a high degree of uniformity in the root zone. To further justify the uniform water content approximation invoked above, we present in this section an analysis of steady flow that extends the work of Raats (1974).
View larger version (14K):
[in this window]
[in a new window]
|
FIG. 4. Illustration of the Michaelis–Menten and two-line uptake models. Also shown is a dilute solution approximation of the Michaelis–Menten model.
|
|
Consider the specific case of a uniform soil where the steady water flux is specified by the Buckingham-Darcy law, in which case Eq. [1] becomes
 | [9] |
where h is the pressure head (L) and K is the hydraulic conductivity (L T–1). The hydraulic conductivity is assumed to be of the form (Gardner, 1958)
 | [10] |
where Ks (L T–1) is the saturated hydraulic conductivity and 
(L–1) is the Gardner parameter that ranges approximately from 0.1 cm–1 for coarse-textured soils to 0.01 cm–1 for fine-textured soils (Raats, 1974). Using the root water uptake term of Raats (1974), rw(z) is given by
 | [11] |
Equation [11] prescribes an exponential uptake distribution, with the distribution parameter 
(L) corresponding to the depth of an equivalent uniform root system having the same uptake and transpiration rates (Raats, 1974). Although Eq. [11] accommodates an infinite rooting depth, 99.8% of the uptake occurs within a depth of 6.
For the soil surface boundary condition,
 | [12] |
and the lower boundary condition at a depth z = L,
 | [13] |
the solution of Eq. [9] for L 
may be expressed as
 | [14] |
Raats (1974, Eq. [21]) presented an equivalent result expressed in terms of the matric flux potential instead of the pressure head. Equation [14] indicates that the pressure head will decrease downward, from a high value at the soil surface down to some constant value as z (or, more practically, as z 6). Profile uniformity can be evaluated by calculating the difference between the pressure head at the soil surface and at depth, 
h h(0) – h(). From Eq. [14] it follows that this pressure head difference is given by
 | [15] |
in which
 | [16] |
is termed the leaching fraction (i.e., the fraction of q0 that passes below the root zone).
Figure 1
shows a plot of h versus the dimensionless parameter for various values of LF. The limit 0 corresponds to a very fine soil and/or a shallow root distribution, whereas implies a coarse-textured soil and/or very deep rooting (Raats, 1974). Figure 1 indicates that the difference between the pressure head at the soil surface and at the effective base of the root zone is smaller (i.e., the profile becomes more uniform) when LF is large and/or is small. This means that high leaching fractions (large LF), fine-textured soils (small ), and shallow rooting depths (small ) promote uniformity in the pressure head profile.
To make these generalities more quantitative, consider that the maximum rooting depth, RD, for common agricultural crops ranges from about 100 to 300 cm (Borg and Grimes, 1986), with the global-average value being about 210 cm (Canadell et al., 1996). Given that 99.8% of uptake occurs within a depth of 6, we may characterize shallow-, average-, and deep-rooted crops as having = RD/6 20, 35, and 50 cm, respectively. Among noncrop plant species, comparable global-average values are RD = 700 cm for trees ( 120 cm), RD = 510 cm for shrubs ( 85 cm), and RD = 260 cm for herbaceous plants ( 45 cm) (Canadell et al., 1996).
Figure 2
shows h plotted versus RD (= 6) for fine-textured ( = 0.01 cm–1) and coarse-textured ( = 0.1 cm–1) soils and various values of LF. For = 0.1 cm–1 (dashed lines), h changes very little as RD increases from 100 (shallow rooted crop) to 700 cm (tree), indicating that rooting depth is of lesser importance in determining pressure head uniformity in the coarse soil. The effect of the leaching fraction is more significant, with h increasing approximately from 7 to 22 cm when LF is decreased from 0.5 to 0.1 (across rooting depths of interest, RD > 100 cm). For the fine-textured soil ( = 0.01 cm–1, solid lines), the rooting depth is far more important. For example, when LF = 0.3, h increases from 29 to 81 cm when RD is increased from 100 to 700 cm. The leaching fraction is also significant for the fine-textured soil: when RD = 300 cm (deep-rooted crop), h increases from 29 to 139 cm when LF is decreased from 0.5 to 0.1.
The largest calculated h values were approximately 180 cm for the fine-textured soil (LF = 0.1, RD = 700 cm) and 22 cm for the coarse-textured soil (LF = 0.1, RD > 100 cm). These results correspond to low leaching and deep rooting scenarios, such as might be found in semiarid and arid climates. Even for these "worst-case" scenarios, however, the variation in the water content may not be excessive, depending on the nonlinearity of the retention function, (h), and the magnitude of d/dh. The nonlinearity in (h) means that the change in water content for a given h will depend on the pressure head at the soil surface, h(z = 0). For fine-textured soils where d/dh is typically very modestly sized for all h values, a small h will correspond to a small change in water content, whereas for coarse-textured soils, d/dh may become large such that a small h can produce a relatively large change in water content if h(z = 0) is in the vicinity of the steep part of the retention curve.
The water content variation can be assessed quantitatively using the Russo (1988) model to describe the saturation:
 | [17] |
where S is effective saturation (0 S 1) and m is Mualem's (1976) tortuosity parameter, which we assume has the value m = 0.5. The difference in saturation between the surface and the base of the root zone is S S[h(0)] – S[h()] = S[h(0)] – S[h(0) – h]. Figure 3
shows plots of S as a function of h(0). These plots were made using the leaching fraction and root depth combinations that produced the largest and smallest values of h observed in Fig. 2 (where RD 100 cm). Thus, the water content variations shown in Fig. 3 correspond to the greatest- and least-uniform pressure head profiles plotted in Fig. 2. As expected, Fig. 3 demonstrates that in addition to the aforementioned effects of rooting depth and leaching fraction, the uniformity of the soil water content profile is affected by the pressure head (or water content) at the soil surface. In evaluating Fig. 3, our experience suggests that when S is less than about 0.1, the water content variability can be safely treated as negligible since that level of variability is similar to the precision that is obtainable with, for example, electromagnetic methods of measuring water contents under field conditions. While it is not possible to state precisely the general conditions for which an effective, uniform water content approximation is acceptable, we conclude that the approximation is reasonable whenever the leaching fraction is sufficiently large (e.g., 0.25). For coarse-textured soils where d/dh is large near the wet end of the retention curve, there is an additional requirement that the soil surface be drier than where the steep part of the curve is located (for the coarse soil depicted in Fig. 2, h(z = 0) must be less than approximately –75 cm).
Thus far, this section has been concerned with the internal consistency of the proposed transport model. That is, we have addressed the question of whether the water content can be reasonably approximated as uniform given the model assumptions of uniform soil properties and steady flow. Once that question is resolved, a related but separate question that may be asked is whether a model that specifies a uniform water content can be used to assess transport in fields where the water content is clearly nonuniform (e.g., in layered soil profiles). Evidence exists (e.g., Wierenga, 1977; Streck and Piehler, 1998) that solute distributions modeled using an effective uniform water content are not significantly different from those modeled using an explicit representation of the water content variability. Thus, we conclude that Eq. [8] may be useful in the case of nonuniform soil profiles if an appropriate effective water content can be determined.
First-Order Plant Uptake
The proposed transport model assumes that solute uptake is a first-order process (Eq. [7]). The plant uptake of solutes is a highly complex process that depends on plant species, solute species, and solution composition (e.g., the uptake of cadmium is affected by the presence of other divalent cations) (Tinker and Nye, 2000). Detailed, mechanistic descriptions of uptake–absorption processes are possible, but these defy a simple mathematical description. Measurements of uptake kinetics commonly show that uptake increases linearly with solution concentration at low concentrations and then levels off and becomes constant at high concentrations. Consistent with these observations, uptake is often modeled assuming Michaelis-Menten kinetics (Tinker and Nye, 2000):
 | [18] |
where KM is the Michaelis–Menten constant. Alternatively, uptake may be modeled using two intersecting lines (Tinker and Nye, 2000), one representing a linear increase in uptake below some critical or threshold concentration, and one representing a constant maximum uptake rate above the critical concentration (Fig. 4
).
When the concentration remains below the threshold concentration, the first-order uptake model is consistent with the two-line model, with the uptake coefficient being proportional to the slope of the first line segment (Fig. 4). A connection with the Michaelis–Menten representation also exists at low concentrations (KM >> C), where the Michaelis–Menten model may be approximated as rs 
C/KM. In this case, the uptake coefficient is 1/KM. Note that if is related to KM in this way, then uptake calculated for a given concentration with Eq. [7] will always be greater than or equal to that specified by the corresponding Michaelis–Menten model (Fig. 4).
Analytical Solution
We define D0 D(z = 0) and consider the solution of Eq. [8] subject to the following boundary and initial conditions
 | [19] |
 | [20] |
 | [21] |
where C0 is the concentration of the infiltrating solution and CI is the initial concentration of the soil water. An analytical solution may be obtained using the generalized integral transform technique, or GITT (Ozisik and Murray, 1974; Cotta, 1993; Liu et al., 2000). Liu et al. (2000) reported the GITT solution of a generic transport equation of which Eq. [8] is a special case. For the present problem, the solution given by Liu et al. (2000) may be written as
 | [22] |
where the infinite summation required in the formal solution has been truncated at M terms and where the eigenfunction 
n and its norm Nn are, respectively,
 | [23] |
 | [24] |
The eigenvalues ßn are the first M zeros of the equation
 | [25] |
Because of the specific form of Eq. [8] and the boundary conditions imposed, it is possible to simplify the expression presented by Liu et al. (2000) for Tn(t). Letting T(t) be a M-length vector with nth element Tn(t), values for Tn(t) in Eq. [22] can be obtained by computing (cf. Liu et al., 2000, Eq. [13])
 | [26] |
where elements of the M x M matrices A and B and the M-length vectors G and T(0) are given by
 | [27] |
 | [28] |
 | [29] |
 | [30] |
A computer program implementing the presented solution is available from the senior author. For simplicity we considered in this study a constant inlet concentration, C0(t) = C0, and a constant initial concentration, CI(z) = CI. A solution for the case of a time-varying inlet and depth-varying initial concentration can also be obtained, although in that case the computation of Tn is more complicated due to the time dependence of Gn (see Liu et al., 2000).
It is also of interest to consider the steady-state (
C/t = 0) solution of the transport equation (Eq. [8]) for the special case of no dispersion (D = 0). In this case, a solution can be obtained by separation-of-variables, with the result being
 | [31] |
where the leaching fraction is defined as before, that is, LF = (v0 – vT)/v0 = (q0 – T)/q0.
|
Results and Discussion
|
|---|
General Model Characteristics
Figures 5
and 6
contain plots that illustrate general characteristics of the advection–dispersion model with uptake of water and solute. The scenarios depicted in these figures involve, beginning at time t = 0 d, the application of water with solute concentration C0 to a soil that has an initially uniform concentration of CI = 0.25 C0. The solid lines in the plots are concentration vs. depth profiles computed with the analytical solution presented in the previous section (Eq. [22]). The dashed lines are steady-state profiles computed with Eq. [31] for D = 0. The simulations assume an exponential water uptake distribution (Eq. [11]) and a dispersion coefficient that was specified as
 | [32] |
where L is the dispersivity (L). Model and algorithm parameter values are given in the figure captions. Note that the computational domain length L was much larger than the plotted depth such that the plotted results are for an effectively semi-infinite profile. Also, for simplicity we assumed a value of R = 1 for the retardation factor; using a larger R value merely shifts the plots to later times.
View larger version (19K):
[in this window]
[in a new window]
|
FIG. 6. Plots of solute concentration vs. depth computed for various values of the transport parameters and L. Except where indicated otherwise, the concentration profiles were computed using the following parameter values: t = 1000 d, R = 1, = 20 cm, v0 = 1 cm d–1, vT = 0.75 cm d–1, = 0.1, L = 10 cm, CI = 0.25 C0, L = 1000 cm, and M = 80. Dashed lines are steady-state results computed for D = 0. (t = time; R = retardation factor; = root distribution parameter; q0 = water flux at soil surface; T = transpiration rate; = volumetric water content; v0 q0/; vT T/; = first-order solute uptake coefficient; L = dispersivity; D = dispersion coefficient; CI = initial soil solution concentration; C0 = infiltrating solution concentration; L = depth of computational domain; M = number of terms summed in analytical solution.)
|
|
Figure 5 shows concentration profiles computed for successive times. At very early times (5 d, in this example), the solute profile is qualitatively similar to that predicted with standard advection–dispersion models. The similarity ends quickly, however, as the concentration behind the front becomes larger than C0 while the leading edge retains the familiar dispersive shape. This solute build-up is due to the concentrating effects of a water flux that decreases with depth. At later times the solute profile reaches a steady state, with solute concentrations increasing from a minimal value at the soil surface to a maximal value at and below the base of the root zone. The maximum concentration that will be achieved is a function of the amount of water and solute being taken up, as determined by the leaching fraction (LF) and the solute uptake coefficient (). When 
0, the maximum concentration is also affected by the value of the dispersion coefficient.
Figure 6 illustrates the effects of the uptake coefficient () and dispersivity (L) on the solute profile. The profiles shown as solid lines are for a large time value such that solute concentrations within the root zone (z 6) have obtained approximately their steady values, whereas concentrations a short distance below the root zone are still increasing. Figure 6a shows that when no solute is excluded from the uptake process ( = 1), root zone solute concentrations do not increase above C0; the steady-state profile is simply C(z) = C0. For = 0, no solute is taken up by the roots and the build-up of solute is maximized. Also shown is the result for the intermediate case = 0.5. For = 1.5, active uptake of solute occurs. In this particular case the steady solute concentration profile decreases with depth, approaching a constant value that is less than C0 but greater than CI.
Figure 6b illustrates the effects of dispersion on solute profile development. The profiles were computed for three values of the dispersivity: 2, 10, and 30 cm. As expected, the results show that increasing the value of the dispersivity produces more solute spreading at the front, which increases the time required for the profile to reach steady state: the L = 2 cm profile at t = 1000 d has very nearly reached steady state, whereas the L = 30 cm profile is still increasing in the bottom of the root zone. Also as expected, the computed solute profile approaches the steady-state, D = 0 solution (dashed line, Eq. [31]) as the dispersivity approaches zero.
The computed profiles in Fig. 6b also illustrate a possible limitation of the advection–dispersion model as applied to transport with simultaneous water and solute uptake. The problem is that when solute builds up in the root zone, a positive concentration gradient develops (i.e., directed back toward the soil surface), which, according to the advection–dispersion model, produces an accompanying diffusive–dispersive flux back toward the soil surface. The magnitude of this flux in reality should reflect only diffusive processes rather than field-scale dispersion processes. The use of a diffusion–dispersion coefficient value that reflects field-scale dispersion processes hence causes an excessively large upward component of the solute flux, leading to near-surface root zone concentrations that are too large. The higher the dispersivity, the more the surface concentration will exceed C0. The field-scale dispersivity typically has a value in the range of 5 to 20 cm (Jury et al., 1991), although this parameter often carries a high degree of uncertainty. Given this uncertainty and the possible limitations of the advection–dispersion model, it would seem prudent to err on the side of choosing a relatively small value of L (or D) when using this model. This limitation is not exclusive to the model developed here; the same considerations should apply to any advection–dispersion model that features root uptake, such as the comprehensive numerical models noted in the introduction.
Example Application
Cadmium is one of the most mobile and bioavailable trace metals and constitutes a potential human health risk through dietary intake (Buchet et al., 1990). In many locations the cadmium content of arable soils has significantly increased over the years, primarily from application of cadmium-containing phosphorous fertilizers and from atmospheric deposition. The cadmium content in Swedish arable soils, for example, has apparently increased by about 30% since the beginning of the twentieth century (Andersson, 1992). The current soil content of cadmium can result in concentrations in harvested plant parts that are close to acceptable limits set for human consumption in many countries (e.g., Eriksson et al., 2000).
Legislative standards for allowable or critical metal loadings to arable land vary by country. One emerging method for risk assessment involves simple mass balance calculations wherein the soil is treated as a single compartment subject to first-order uptake into harvested plant parts, linear sorption partitioning in the soil, and advective leaching below the root zone (Pa
es, 1998; Tiktak et al., 1998; Blombäck et al., 2000; Bergkvist et al., 2005). Such an approach implicitly assumes uniform well-mixed conditions within the soil compartment and ignores depth variations of transport and uptake of both water and trace metals in the root zone. Given these assumptions, and using the same notation as above, the change in stored amount A (M L–2) is given by
 | [33] |
where qz is the water flow at the base of the root zone and A = ( + bKd)RdC. Noting that qz = q0 – T leads to
 | [34] |
where I = (v0C0)/(RRD) and k = [v0 + vT( – 1)]/(RRD). When C(0) = CI, integrating Eq. [34] from t = 0 to t gives
 | [35] |
In the one-compartment model, the concentration of harvested parts of annual plants Cp (M m–1) is given by
 | [36] |
where Bp is the harvested biomass yield (M L–2 T–1). This expression is appropriate only for plants whose growth cycle is short relative to the timescale of solute transport in the root zone (e.g., annual crops). Setting dC/dt = 0 in the one-compartment model gives the steady-state soil solution concentration as
 | [37] |
and the steady-state plant concentration as
 | [38] |
In our example application, we compare the analytical advection–dispersion and one-compartment models for the case of cadmium uptake into wheat grain, an important staple food in many parts of the world. The advection–dispersion model in this example was implemented using a constant dispersion coefficient, D(z) = D, and a linear uptake distribution,
 | [39] |
This uptake distribution conforms to the "40:30:20:10 rule" (e.g., Raats, 1974; Hoffman and van Genuchten, 1983), such that 40% of the uptake is in the top quarter of the root zone, 30% in the second quarter, 20% in the third quarter, and 10% in the bottom quarter. Omitting details, an expression equivalent to Eq. [36] may be developed for the advection–dispersion model to compute the harvested plant part concentration:
 | [40] |
The uptake coefficient is an important parameter since it determines the extent to which cadmium in the soil solution is excluded from the plant. Measurements of cadmium uptake have indicated that uptake is linear over a range of concentrations relevant for sludge applications (Hamon et al., 1999; McGrath et al., 2000). Measurements of cadmium concentrations in soil solutions and wheat grain (Bergkvist et al., 2005), together with typical values of the water use efficiency T/Bp for wheat (see Eq. [36]), suggest a value for of approximately 0.05, which we used in the calculations presented below. The remaining parameter values are shown in the caption of Fig. 7
, which compares predictions of the two models. The assumed cadmium loading rate (1.6 mg m–2 yr–1) reflects all sources of cadmium. In the wheat-growing areas of Europe, diffuse aerial deposition, together with impurities in lime, commercial fertilizer, and feeds, would not supply more than approximately 0.3 to 0.5 mg Cd m–2 yr–1. Thus, for the scenario depicted in Fig. 7, we implicitly assume a significant additional local cadmium source such as industrial emissions or sludge applications. The assumed loading rate is not unrealistic in that the current allowable limit in the EU for cadmium supplied with sludge is set at 15 mg m–2 yr–1 (EU directive, 86/278/EEC).
View larger version (13K):
[in this window]
[in a new window]
|
FIG. 7. A comparison of predictions of the analytical advection-dispersion and one-compartment models for cadmium uptake into wheat grain using the following parameter values: q0 = 0.7 m yr–1, T = 0.5 m yr–1, = 0.05, Bp = 0.5 kg m–2 yr–1, = 0.4 m3 m–3, R = 376, D = 0.015 m2 yr–1, CI = 1 mg m–3 and C0 = 2.286 mg m–3, equivalent to a cadmium loading of 1.6 mg m–2 yr–1. The root depth, RD, in the advection–dispersion equation model was set to 1 m, and two values of RD are shown for the one-compartment model (0.25 and 1 m). (q0 = water flux at soil surface; T = transpiration rate; = first-order solute uptake coefficient; Bp = harvested biomass yield; = volumetric water content; R = retardation factor; D = dispersion coefficient; CI = initial soil solution concentration; C0 = infiltrating solution concentration; RD = maximum rooting depth.)
|
|
Figure 7 shows that for the above parameterization, the analytical advection–dispersion model predicts a final steady-state Cd concentration in wheat grain of approximately 0.2 mg kg–1, reached after about 500 yr, which equals the allowable concentration under EU directives. The root depth in the one-compartment model is usually set to the plow depth (
0.25 m) since this best approximates the assumption of uniform mixing (e.g., Tiktak et al., 1998). This assumption ignores subsoil uptake of cadmium, which can be considerable (Johnsson et al., 2002). Figure 7 shows that with RD = 0.25 m, the one-compartment model substantially overestimates the grain cadmium concentration, although the pattern of increase is similar to the analytical advection–dispersion model, with a steady-state value approached after approximately 500 yr. Increasing RD to 1 m, to match the maximum depth of roots in the advection–dispersion model, produces good agreement with the advection–dispersion model at early times (up to 300 yr). However, the final steady-state value of the one-compartment model is independent of RD (see Eq. [38]), being about 75% larger than the value predicted by the advection–dispersion model. In terms of risk assessment, such a difference in model results is highly significant. We conclude that the simplifications in the one-compartment model can lead to erroneous decisions, although the predicted deviations in the present example would become significant only after a few hundred years depending on the value of RD used in the one-compartment model. Finally, we note that the analytical advection–dispersion model requires only a few extra parameters, which are easily estimated. The advection–dispersion model therefore should be recommended as the preferred tool for environmental risk assessment of persistent contaminants.
|
Conclusions
|
|---|
Although a variety of comprehensive numerical simulation models exist for analyzing solute transport processes in the root zone, regulatory and risk management analyses often require simpler modeling tools with fewer input parameters. Simpler analytical models may be especially useful for assessing environmental impacts of strongly sorbing and persistent solutes such as trace metals. This is because the relevant timescales for such assessments are measured in decades or centuries, in which case transport and plant uptake can be modeled effectively using annual and root-zone average parameters.
In this work we developed an analytical advection–dispersion model that includes plant uptake of water and solute. Although considerably simpler than many popular comprehensive numerical models, the analytical advection–dispersion model—with its explicit accounting of advection, dispersion, and uptake processes within the root zone—is much more realistic than one-compartment, well-mixed models that are currently being used or considered for use in risk assessments. As an example, we compared predictions made with the advection–dispersion and one-compartment models for the uptake of cadmium by wheat. The comparison showed that the one-compartment model significantly overestimated the long-term, steady-state Cd concentration in harvested wheat grain. For this reason, we recommend the analytical advection–dispersion model as a preferred tool for environmental risk assessment of strongly sorbing, persistent solutes.
Possible limitations of the presented model arise due to assumptions of linear solute sorption and first-order plant uptake. We thus suggest that model applications be restricted to cases of slight to moderate soil contamination where these assumptions are likely to be reasonable. Additionally, the model was derived based on an assumption of uniform soil water content. Among other factors, water content uniformity is affected by soil texture, rooting depth, and leaching fraction, with fine texture, shallow rooting, and high leaching all promoting uniformity. In general, a leaching fraction of greater than about 0.25 corresponds to a fairly uniform profile when steady-state water flow occurs, except perhaps in very coarse-textured soils where substantial nonuniformity could arise if the surface soil were to be maintained at a high degree of saturation. Where nonuniform soil conditions prevail, the model may be used if an appropriate, effective water content can be identified.
|
ACKNOWLEDGMENTS
|
|---|
This research was supported in part by the Swedish Environmental Protection Agency funded project "Cadmium leaching from arable land and uptake into crops: calculation of critical loads and development of a simulation tool for dynamic long-term predictions"; a fellowship from CNPq (Brazil) to E. M. Pontedeiro; and the SAHRA Science and Technology Center as part of NSF grant EAR-9876800.
|
REFERENCES
|
|---|
- Ahuja, L.R., K.W. Rojas, J.D. Hanson, M.J. Shaffer, and L. Ma (ed.). 2000. Root zone water quality model: Modeling management effects on water quality and crop production. Water Resources Publications, LLC, Highlands Ranch, CO.
- Andersson, A. 1992. Trace elements in agricultural soils: Fluxes, balances, and background values. Report 4077. Swedish Environmental Protection Agency, Solna.
- Bergkvist, P., N.J. Jarvis, J. Eriksson, and L. Rapp. 2005. Critical load of cadmium on arable soils in Sweden. Emergo no. 4, Studies in the Biogeophysical Environment. Dep. of Soil Sci., Swedish Univ. of Agricultural Sciences, Uppsala.
- Blombäck, K., M. Kalvi, J. Eriksson, and I. Öborn. 2000. Cadmium accumulation in agricultural soils and uptake by plants as an effect of cadmium in phosphate fertilizers: Model estimations. p. 11–31. In Assessment of risks to health and the environment in Sweden from cadmium in fertilizers. Report 4/00. National Chemicals Inspectorate, Stockholm.
- Borg, H., and D.W. Grimes. 1986. Depth development of roots with time: An empirical description. Trans. ASAE 29:194–197.[Web of Science]
- Buchet, J.P., R. Lauwerys, H. Roels, A. Bernard, P. Bruaux, F. Claeys, G. Ducoffre, P. de Plaen, J. Staessen, A. Amery, P. Lijnen, L. Thijs, D. Rondia, F. Sartor, A. Saint Remy, and L. Nick. 1990. Renal effects of cadmium body burden of the general population. Lancet 336:699–702.[CrossRef][Web of Science][Medline]
- Canadell, J., R.B. Jackson, J.R. Ehleringer, H.A. Mooney, O.E. Sala, and E.-D. Schulze. 1996. Maximum rooting depth of vegetation types at the global scale. Oecologia 108:583–595.[CrossRef][Web of Science]
- Cotta, R.M. 1993. Integral transforms in computational heat and fluid flow. CRC Press, Boca Raton, FL.
- Dalton, F.N., P.A.C. Raats, and W.R. Gardner. 1975. Simultaneous uptake of water and solutes by plant roots. Agron. J. 67:334–339.[Abstract/Free Full Text]
- Destouni, G. 1991. Applicability of the steady state flow assumption for solute advection in field soils. Water Resour. Res. 27:2129–2140.[CrossRef]
- Eriksson, J., A. Andersson, B. Stenberg, and R. Andersson. 2000. Current status of Swedish arable soils and cereals. Report 5062. (Summary in English.) Swedish Environmental Protection Agency, Stockholm.
- Feddes, R.A., P.J. Kowalik, and H. Zaradny. 1978. Simulation of field water use and crop yield. PUDOC, Wageningen, the Netherlands.
- Gardner, W.R. 1958. Some steady-state solutions of the unsaturated moisture flow equation with application to evaporation from a water table. Soil Sci. 85:228–232.
- Ginn, T.R., and E.M. Murphy. 1997. A transient flux model for convective infiltration: Forward and inverse solutions for chloride mass balance studies. Water Resour. Res. 33:2065–2079.[CrossRef]
- Hamon, R.E., P.E. Holm, S.E. Lorenz, S.P. McGrath, and T.H. Christensen. 1999. Metal uptake by plants from sludge-amended soils: Caution is required in the plateau interpretation. Plant Soil 216:53–64.[CrossRef][Web of Science]
- Hoffman, G.J., and M.Th. van Genuchten. 1983. Soil properties and efficient water use: Water management for salinity control. p. 73–85. In H.M. Taylor, W.R. Jordan, and T.R. Sinclair (ed.) Limitations to efficient water use in crop production. ASA, CSSA, and SSSA, Madison, WI.
- Johnsson, L., I. Öborn, G. Jansson, and D. Berggren. 2002. Evidence for spring wheat (Triticum aestivum) uptake from subsurface soils. In Y. Zhu and Y. Tong (ed.) Proc. of the Int. Workshop of Soil-Plant Interactions, Beijing, China. 27–30 Oct. 2002. Kluwer Academic, Dordrecht, the Netherlands.
- Jury, W.A., W.R. Gardner, and W.H. Gardner. 1991. Soil physics. John Wiley & Sons, New York.
- Larsbo, M., S. Roulier, F. Stenemo, R. Kasteel, and N.J. Jarvis. 2005. An improved dual-permeability model of water flow and solute transport in the vadose zone. Vadose Zone J. 4:398–406.[Abstract/Free Full Text]
- Liu, C., J.E. Szecsody, J.M. Zachara, and W.P. Ball. 2000. Use of the generalized integral transform method for solving equations of solute transport in porous media. Adv. Water Resour. 23:483–492.[CrossRef]
- McGrath, S.P., F.J. Zhao, S.J. Dunham, A.R. Crosland, and K. Coleman. 2000. Long-term changes in the extractability and bioavailability of zinc and cadmium after sludge application. J. Environ. Qual. 29:875–882.[Abstract/Free Full Text]
- Mualem, Y. 1976. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12:513–522.[CrossRef]
- Ozisik, M.N., and R.L. Murray. 1974. On the solution of linear diffusion problems with variable boundary condition parameters. J. Heat Transfer 96:48–51.[Web of Science]
- Paes, T. 1998. Critical loads of trace metals in soils: A method of calculation. Water Air Soil Pollut. 105:451–458.[CrossRef]
- Raats, P.A.C. 1974. Steady flows of water and salt in uniform soil profiles with plant roots. Soil Sci. Soc. Am. Proc. 38:717–722.
- Raats, P.A.C. 1975. Distributions of salts in the root zone. J. Hydrol. 27:237–248.[CrossRef]
- Raats, P.A.C. 1981. Residence times of water and solutes within and below the root zone. Agric. Water Manage. 4:63–82.[CrossRef]
- Rawlins, S.L. 1973. Principles of managing high frequency irrigation. Soil Sci. Soc. Am. Proc. 37:626–629.
- Russo, D. 1988. Determining soil hydraulic properties by parameter estimation: On the selection of a model for the hydraulic properties. Water Resour. Res. 24:453–459.
- Schoups, G., and J.W. Hopmans. 2002. Analytical model for vadose zone solute transport with root water and solute uptake. Vadose Zone J. 1:158–171.[Abstract/Free Full Text]
- imnek, J., D.L. Suarez, and M. ejna. 1996. The UNSATCHEM software package for simulating one-dimensional variably saturated water flow, heat transport, carbon dioxide production and transport, and multicomponent solute transport with major ion equilibrium and kinetic chemistry, Version 2.0. Res. Rep. no. 141. U.S. Salinity Lab., Riverside, CA.
- imnek, J., M.Th. van Genuchten, and M. ejna. 2005. The HYDRUS-1D software package for simulating the one-dimensional movement of water, heat and multiple solutes in variably-saturated media, Version 3.0. HYDRUS Software Ser. no. 1. Dep. of Environ. Sci., Univ. of Calif., Riverside.
- Skaggs, T.H., and F.J. Leij. 2002. Solute transport: Theoretical background. p. 1353–1380. In J. Dane and C. Topp (ed.) Methods of soil analysis. Part 4. ASA and SSSA, Madison, WI.
- Skaggs, T.H., M.Th. van Genuchten, P.J. Shouse, and J.A. Poss. 2006. Macroscopic approaches to root water uptake as a function of water and salinity stress. Agric. Water Manage. 86:140–149.[CrossRef]
- Streck, T., and H. Piehler. 1998. On field-scale dispersion of strongly sorbing solutes in soils. Water Resour. Res. 34:2769–2773.[CrossRef]
- Tiktak, A., R. Alkemade, H. van Grinsven, and K. Makaske. 1998. Modeling cadmium accumulation at a regional scale in the Netherlands. Nutr. Cycling Agroecosyst. 50:209–222.[CrossRef]
- Tinker, P.B., and P.H. Nye. 2000. Solute movement in the rhizosphere. Oxford Univ. Press, New York.
- van Genuchten, M. Th., and W.J. Alves. 1982. Analytical solutions of the one-dimensional convective-dispersive solute transport equation. USDA Tech. Bull. 1661. U.S. Gov. Print. Office, Washington, DC.
- Wierenga, P.J. 1977. Solute distribution profiles computed with steady-state and transient water movement models. Soil Sci. Soc. Am. J. 41:1050–1055.[Abstract/Free Full Text]
- Yuan, F., and Z. Lu. 2005. Analytical solutions for vertical flow in unsaturated, rooted soils with variable surface fluxes. Vadose Zone J. 4:1210–1218.[Abstract/Free Full Text]
TOP
ABSTRACT
INTRODUCTION
