Published online 1 January 2007
Published in Soil Sci Soc Am J 71:26-34 (2007)
DOI: 10.2136/sssaj2005.0395
© 2007 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
SOIL PHYSICS
An Inverse Method to Estimate the Source-Sink Term in the Nitrate Transport Equation
Jianchu Shi and
Qiang Zuo*
Dep. of Soil and Water Sciences and Key Lab., of Plant-Soil Interactions, MOE, College of Resources and Environment, China Agricultural Univ., Beijing 100094, China
Renduo Zhang
School of Environ. Science and Eng., Sun Yat-Sen (Zhongshan) Univ., Guangzhou 510275, China
* Corresponding author (qiangzuo{at}cau.edu.cn).
 |
ABSTRACT
|
|---|
The source-sink term (SST) in the convection-dispersion equation (CDE) to simulate nitrate (NO3––N) transport integrates several NO3––N transformation processes in soils. The term is affected by considerably complicated micro-environmental conditions and is difficult to be measured directly. In this study, the average SST distribution was estimated by solving the CDE with an unknown SST iteratively using an inverse method. The required input information for the SST estimation was easily obtained, including soil hydraulic properties, two successively measured NO3––N concentration distributions, and the boundary and initial conditions. Numerical experiments were designed to examine the accuracy and stability of the inverse approach, considering spatial intervals of measurement data along the soil profile, time intervals between the successive measurements of soil NO3––N concentration, boundary conditions, layered soils, and measurement errors of NO3––N concentration. Comparisons with theoretical results showed that the inverse method was reliable for estimating the SST in the CDE. Data from a column experiment with winter wheat (Triticum aestivum L. ) growth were used to demonstrate applications of the inverse method. The root-nitrate-uptake (RNU) rate distributions were estimated according to the proposed inverse procedure and the soil NO3––N transport with RNU in the columns was simulated. The simulated soil NO3––N concentration distributions were comparable with the measured values. The relative errors between the simulated and measured values of the total N mass extracted by winter wheat were <10%.
|
INTRODUCTION
|
|---|
Extensive use of agricultural chemicals has caused environmental pollution problems (Schmied et al., 2000). As one of the most ubiquitous chemical contaminants, the NO3––N concentration in surface and ground waters continues to increase (Spalding and Exner, 1993). To balance the N supply for optimal plant performance and minimal losses to the environment, accurate simulations of N dynamics become critical. The CDE combined with a SST is the most commonly used process representation for predicting solute dynamics (Vasssilis and Wyseure, 1998). Enormous effort has been made on modeling N dynamics. However, the dynamics are still poorly understood because of difficulties in determining the SST accurately (Ma and Shaffer, 2001). The SST involves several complicated processes in soils, such as immobilization, nitrification, denitrification, and uptake by plant roots. Comprehensive SST models have been developed, but are still limited by inadequate descriptions of the simultaneous processes of N turnover (Hansen et al., 1995) and incomplete definitions of input parameters (De Willigen, 1991).
The SST parameters (e.g., the rates of nitrification, denitrification, immobilization, and root uptake), most of which cannot be measured directly, are influenced by many factors such as soil water content, N concentration, root distributions, organic matter content, and others (Hansen et al., 1995). The trial-error method is often used to obtain coefficients related to those rates, which may not be optimized in a strict mathematical sense (Shaffer et al., 2001). On the other hand, the inverse method optimizes the parameters related to soil water flow or solute transport in soils by utilizing data from transient transport experiments and presents an attractive alternative for parameters estimations (Schmied et al., 2000; Priesack et al., 2001). Usually, the inverse method can be successfully used to estimate some less time-independent parameters (e.g., hydraulic parameters such as the saturated hydraulic conductivity) with relatively easily measured data such as soil water contents (Zijlstra and Dane, 1996). As for dynamic or time-dependent terms, for example, the root-water-uptake (RWU) rate distributions, their average status over a time period is more meaningful. The RWU rate distributions are important for understanding the dynamic processes of soil water flow but cannot be measured directly. Recently, an inverse method to estimate the average RWU rate distribution over a time period was developed by Zuo and Zhang (2002). The method uses two successively measured soil water content distributions to estimate the SST in the soil water flow equation, namely, the RWU rate. The method can be applied to analyze the changing patterns of RWU and to optimize the parameters for setting up the RWU model and simulating soil water flow (Zuo et al., 2004). Similarly, nitrate distribution in the soil profile is relatively easy to measure, and it reflects the interrelationship among soil NO3––N transformations in soils. Therefore, it can be used to estimate the average SST distribution of the CDE with an inverse method. Nevertheless, most inverse methods have some limitations mainly related to the non-uniqueness and instability of the optimized parameter set (Hupet et al., 2003). Furthermore, in most publications on inverse methods, the estimation of unknown parameters has been limited to a homogeneous soil, with only a few studies for layered soils (Zijlstra and Dane, 1996; Priesack et al., 2001).
The objective of this study was to apply an inverse method similar to that of Zuo and Zhang (2002) to estimate the SST of NO3––N in the CDE by using input information of two successively measured NO3––N concentration distributions in the soil. Several numerical examples were set up to examine the accuracy and stability of the proposed inverse method. Experimental data were used to show the applicability of the method. Based on the estimated SST, the distribution of RNU rate was calculated, the parameter in the RNU model was optimized, and NO3––N transport in soil columns was simulated.
|
MATERIALS AND METHODS
|
|---|
Water Flow in Soils
Successful simulations of NO3––N dynamics require accurate modeling of soil water flow. One-dimensional vertical soil water flow with RWU is simulated using Richards' equation as follows (van Genuchten, 1987; Wu et al., 1999):
 | [1] |
 | [2] |
 | [3] |
 | [4] |
where h is the soil matric potential (cm); C(h), the soil water capacity (cm–1); K(h), the soil hydraulic conductivity (cm d–1); h0(z), the initial soil matric potential in the profile (cm); E(t), the soil surface evaporation rate (cm d–1); L, the simulating depth (cm) and L 
Lr, in which Lr is the rooting depth (cm); hL(t), the matric potential at L (the lower boundary) (cm); z, the vertical coordinate originating from the soil surface and positive downward (cm); and S(z, t) the RWU rate (cm3 cm–3 d–1), defined by (van Genuchten, 1987; Wu et al., 1999; Musters and Bouten, 2000):
 | [5] |
 | [6] |
 | [7] |
 | [8] |
where 
(h) is a dimensionless reduction function related to the effect of water stress; Smax (z, t), the maximal specific water extraction rate under the optimal soil water conditions (cm3 cm–3 d–1); Tp, the potential transpiration rate (cm d–1); Lnrd(z), the normalized root length density distribution; Y(z), the root length density (cm cm–3); Y(0), the root length density at the soil surface (cm cm–3); h1 and h2 are threshold values of matric potential (cm); 
and ß (cm2) are fitted parameters; zr (=z/Lr), the normalized root depth ranging from 0 to 1.
Nitrate Transport in Soils
One-dimensional vertical movement of NO3––N in the unsaturated zone is characterized by the CDE combined with a SST (Lafolie, 1991; Vasssilis and Wyseure 1998):
 | [9] |
 | [10] |
 | [11] |
 | [12] |
where CN is the concentration of NO3––N, expressed as mass of NO3––N per volume of soil solution (mg cm–3); 
, the soil volumetric water content (cm3 cm–3); q, the Darcy's flux (cm d–1), q = v, in which v is the pore water velocity (cm d–1); CN0(z), the initial NO3––N concentration distribution (mg cm–3); Qs(t), the flux of NO3––N at soil surface (mg cm–2 d–1); CNL(t) the NO3––N concentration at the lower boundary (mg cm–3); SN(z, t) the SST integrating the transformation processes of NO3––N in soils (mg cm–3 d–1); D(,v), the hydrodynamic dispersion coefficient (cm2 d–1), which is a combination of chemical diffusion and mechanical dispersion (Millington and Quirk, 1961):
 | [13] |
where D0 is the diffusion coefficient for NO3––N in pure water (cm2 d–1); s, the saturated water content (cm3 cm–3); 
, the dispersivity (cm).
The SST SN(z, t) in the CDE unifies the transformation processes of NO3––N in soils, expressed as (Hansen et al., 1991):
 | [14] |
where Sn(z, t), Sm(z, t), Sd(z, t), and Su(z, t), respectively, are the rates of ammonium nitrification, NO3––N immobilization, denitrification, and root uptake per unit soil volume (mg cm–3 d–1), and defined by:where k1, k2, k3, respectively, are the optimal first-order rate coefficients for nitrification of ammonium, immobilization, and denitrification of NO3––N; C0(z, t) is the concentration of ammonium in the soil solution (mg cm–3); T'(z, t), the soil temperature (°C); Tm, the optimum temperature for these processes (°C) and chosen as Tm = 35°C in this study (Cabon et al., 1991); f(z), the field water capacity (cm3 cm–3); d(z, t), the threshold water content for denitrification (cm3 cm–3); 
, the dimensionless RNU factor ( 0).
An Inverse Method to Estimate the Source-Sink Term in Convection-Dispersion Equation
The average SST during a period from t = 0 to t = T may be calculated as:
 | [20] |
However, SN(z, t) in the CDE is unknown. We proposed an inverse method to estimate the average NO3––N transformation rate 
N(z,T) in a soil–plant system, using two successively measured NO3––N distributions, namely CN(z, 0) and CN(z, T). The iterative procedure is performed as follows:- Soil matric potential profiles h(z,t) between t = 0 and t = T, correspondingly the distributions of soil water content (z,t), are obtained by solving Eq. [1] through [4]
. Subsequently, Darcy's flux q(z,t) and pore water velocity v(z,t) are calculated.
- The initial value of SST in Eq. [9], N(0)(z,T), is set to be 0. Then, Eq. [9] through [12] are solved using a numerical method for the soil NO3––N concentration distribution at time t = T, namely, CN(1)(z,T), incorporating the information of soil and solute properties and the soil water flow obtained from Step (1). Hence, an initial approximation of soil NO3––N concentration change resulted from the SST between t = 0 and t = T is calculated by
 | [21] |
in which the superscript represents the number of iterations.
- The average SST is rectified as:
 | [22] |
- Eq. [9] through [12] are solved numerically again for CN(2)(z,T) by substituting N(1)(z,T) for SN(z,t) in Eq. [9].
- The iteration process is repeated by:
 | [23] |
 | [24] |
where k represents the kth iterations (k 1). The iteration process continues until 
CN(k)(z,T) is smaller than a specified iterative control value. If there are M measured values of soil NO3––N concentration in the soil profile at time t = T, namely CN(z1, T), CN(z2,T), ...CN(zM, T), the convergence criterion to finish the iterative process is determined by:
 | [25] |
where 
is a specified iterative control value. As the convergence criterion reaches, N(k)(z,T) is used to approximate the average source-sink term N(z,T).
In reality, it is not practical to measure soil NO3––N concentrations at any depth. Furthermore, measurement errors exist in almost all measurement processes. However, continuous, smooth and relatively accurate distributions of soil NO3––N concentration are essential for the above iterative procedure, especially at time t = T. The step linear interpolation, which is often used in practice, can generate continuous concentration distributions but usually in zigzag patterns, which may result in irregularly estimated N(k)(z,T) and even non-convergence of the iteration process. Therefore, an algebraic polynomial is used to fit the measured values to obtain a continuous and smooth distribution of soil NO3––N concentration. The fitting polynomial equation is chosen as (Xu, 1995):
 | [26] |
 | [27] |
where Pm(z) represents an (m – 1)th algebraic polynomial of soil NO3––N concentration, and a1, a2,..., am are optimized parameters. In the fitting process, the m degree of Eq. [26] may start at 3 and increases successively until the maximal absolute value of errors between the measured and fitted data is less than a specified value p, namely,
 | [28] |
where C(zi) represents the measured data; p may be related to the precision of the instrument (e.g., continuous flow analyzer for measuring NO3––N concentration).
Procedure to Test the Accuracy and Convergence of the Inverse Method
Several numerical examples were designed to test the accuracy and convergence of the inverse method considering following factors: spatial interval (SI) and time interval (T) of soil NO3––N concentration measurement; boundary conditions and layered soils; and NO3––N concentration measurement errors. The numerical experiments were designed as follows:
- Input a set of parameters referring to water flow and NO3––N transport and transformation in soils. Table 1 lists the parameters and data used in the numerical experiments. Other parameters include = 10–4 (Eq. [25]) and p = 0.9 x 10–4 mg cm–3 (Eq. [28]).
- Solve Eq. [1] through [4] using the implicit finite difference method to obtain distributions of matric potential h(z, t) based on Eq. [5] through [8] and the input data in Step (1), hence distributions of soil water content (z, t), Darcy's flux q(z, t), and pore water velocity v(z, t).
- Solve Eq. [9] through [12] using the implicit finite difference method to obtain the theoretical distributions of soil NO3––N concentration CN(z, T) at time t = T based on Eq. [14] through [19] and the input data in Step (1). The theoretical average SST distribution N(z,T) was calculated with Eq. [14] through [20].
- Choose some values CN(zi, T) from CN(z, T) according to a specified SI (e.g., 5, 10, 20, or 30 cm, similar to the measurement strategies in the field) as the "measured" data points. In practice, measured NO3––N concentrations usually oscillate around the true values because of errors caused by the measurement instrument. To consider the measurement errors, the "measured" data of soil NO3––N concentration were generated randomly by
 | [29] |
where rand is a random number generator [–p 
rand(p) p].
- Fit the "measured" C*N(zi,T) points to a continuous and smooth NO3––N concentration curve using Eq. [26] through [28].
- Estimate the average distribution of the SST in the CDE using the generated soil NO3––N distributions in Step (4) and (5) and the proposed inverse method.
- Compare the estimated average SST distribution N(k)(z,T) with the theoretical distribution N(z,T) calculated from Step (3). The errors between N(k)(z,T) and N(z,T) were evaluated by the maximal absolute error (MAE) and the overall relative error (ORE):
 | [30] |
 | [31] |
 | [32] |
where TN(T) and T*N(T) indicate, respectively, the theoretical and estimated average transformation mass of NO3––N in the simulating depth per unit area per unit time (mg cm–2 d–1).
View this table:
[in this window]
[in a new window]
|
Table 1. Soil properties, RWU and NO3––N transformation parameters, initial and boundary conditions for Examples (Ex.) 1 to 4.
|
|
View larger version (30K):
[in this window]
[in a new window]
|
Fig. 1. Randomly generated "measured" and fitted NO3––N concentration distributions at different times.
|
|
Estimating the RNU Factor in SST
If the NO3––N transformation parameters in SST are known except the RNU factor ( 0), we can calculate the distribution of average RNU rate u(zi,T) based on the estimated average SST distribution N(k)(zi,T) (i = 1, 2, ..., M):
 | [33] |
In addition, the average RNU rate distribution in the soil profile can be estimated according to the RNU model Eq. [19]:
 | [34] |
where n(zi,T), d(zi,T), m(zi,T), respectively, represent the average rates of nitrification, denitrification, and immobilization per unit soil volume at the measurement depth zi (mg cm–3 d–1), which can be calculated according to the similar approach for N(zi,T) (Eq.[20]); S(zi,0) and S(zi, T) are RWU rates at the measurement depth zi at t = 0 and t = T, respectively, calculated by Eq. [5] to [8]. The least-squares procedure by combining Eq. [33] and [34] can be used to obtain an optimized value (a single value) for the whole soil profile.
Soil Column Experiment
Data collected from a column experiment were used to validate the inverse method and demonstrate its applications in practice. The experiment (Exp. 1) was performed in a greenhouse using 32 PVC columns with a height of 53 cm and a diameter of 15 cm. The columns were packed with a sandy soil up to the height of 50 cm with a bulk density of 1.65 g cm–3. The particle-size distribution of the soil included 92.33% of sand, 7.43% of silt, and 0.24% of clay. The organic matter and total N in the soil were 0.11 and 0.07 g kg–1, respectively. The soil hydraulic parameters are summarized in Table 1 (Example 4).
Winter wheat (Triticum aestivum L. cv. Jingdong 8) was planted in the columns with a seed density of seven plants per column, similar to that in the field (400–600 plants per m2). The experiment lasted for 58 d (from 16 Dec. 2004 to 12 Feb. 2005) from planting to tillering stages of winter wheat. At the soil surface in each column, 3 cm of fine quartz sand were filled to reduce soil surface evaporation on 19 Dec. 2004 (3 d after planting, 3 DAP). The conditions for winter wheat growth in the greenhouse were kept as: a photosynthetic photon flux density of 500 µmol m–2 s–1 over the plant for 12 h d–1 (from 0800 to 2000 h), day/night air temperature within 20/12 ± 2°C, and a relative humidity of 40 ± 5%. Winter wheat was irrigated once every 6 d from 10 DAP (26 Dec. 2004) to keep an average water content in the root zone not less than 60% of the field water capacity, using the half-strength Hoagland solution.
The soil and root sampling work was conducted twice every 6 d (0.5 and 5.5 d after each irrigation), with two columns opened for soil cores each time, namely, four columns in an irrigation period. The first sampling was started on 10.5 DAP (at 0800 h on 27 Dec. 2004) and 16 sampling times in total were performed during the experimental period. At each sampling time, the soil cores were cut into 4-cm soil layers. Some samples were used to measure soil water content with the gravimetric method and soil NO3––N concentration with the continuous flow analyzer (TRAACS 2000, Bran + Luebbe, Norderstedt, Germany; precision p = 0.9 x 10–4 mg cm–3). The remaining samples were put into a meshwork (with grids of 0.05 cm in diameter) and washed out soil for roots. The collected roots in each soil layer and the aboveground biomass from each column were weighed by drying to a constant weight at 70°C for their dry weights and analyzed for their nitrogen contents with an element analyzer (CHNSO EA 1108, Carlo Erba, Italy).
The evapotranspiration rate (ET) was obtained from the water loss by weighing the columns at 2000 h daily. The soil surface evaporation rate (i.e., E(t) in Eq. [3]) was measured from a parallel column experiment (Exp. 2) with three duplicate soil columns, calculating from the water loss by weighing the columns at 2000 h daily. Experiment 2 was conducted in the same way as Exp. 1 but without plant in the columns. Since winter wheat in Exp. 1 grew only during seedling stage and all soil columns in both experiments were mulched with 3 cm fine quartz sands, the soil surface evaporation rate calculated from Exp. 2 was used to approximate the value E(t) in Exp. 1.
|
RESULTS AND DISCUSSION
|
|---|
Effect of Measurement Spatial and Time Intervals on the Estimated SST (Example 1)
Generally speaking, smaller spatial intervals (SI) of NO3––N concentration measurements along the vertical soil profile contain more information and result in better estimations. However, it is tedious and difficult to measure soil NO3––N concentrations at small spatial intervals (e.g., SI < 5 cm) in the field. Because of the dramatic change of soil water and NO3––N concentration near the soil surface, a relatively small spatial interval of measurement was usually needed in this zone. Therefore, a spatial interval SI = 5 cm was adopted near the soil surface (0 z 30 cm). For z > 30 cm, three spatial intervals (SI = 10, 20, and 30 cm) were discussed. The three SI schemes were represented, respectively, by SI = 5–10 cm (i.e., SI = 5 cm for 0 z 30 cm and SI = 10 cm for z > 30 cm), SI = 5–20 cm, and SI = 5–30 cm in the following discussions. For the same reason, smaller measurement time intervals (T) contain more information for the SST changing processes. However, the estimated SST from the proposed inverse method is an average distribution over a time interval T. When the error for soil NO3––N concentration measurement p is fixed, smaller time intervals will result in larger errors between the estimated and the theoretical SST values (Eq. [24]). For example, because p might be 0.9 x 10–4 mg cm–3, the differences between estimated and theoretical SST would reach as high as 4.1 x 10–5, 2.0 x 10–5, and 0.8 x 10–5 mg cm–3 d–1 for T = 1, 2, 5 d, respectively, under the saturated condition (s = 0.45 cm3 cm–3). Hence time intervals of T < 5 d were not considered. On the other hand, if the measurement time interval were too large (e.g., T > 30 d), NO3––N transformation process included in the SST would be masked by the averaging procedure, that is, the dynamic course of the SST could not be caught in time. Therefore, time intervals T = 5, 10, 15, and 30 d, were discussed in this study.
Figure 1
shows the "measured" random points with SI = 5–10 cm and the fitted continuous distributions for various time intervals. The "measured" NO3––N concentrations for SI = 5–20 and 5–30 cm at T = 10 d were generated from those for SI = 5–10 cm and then fitted using Eq. [26] through [28]. Based on the inverse procedure, the average SST distributions for different spatial intervals and time intervals were estimated.
The estimated and theoretical distributions of the SST for different spatial intervals with T = 10 d are shown in Fig. 2a
. In general, the errors between the estimated and theoretical SST increased with increasing SI values, but the estimation results were acceptable with MAE < 3.0 x 10–5 mg cm–3 d–1 and ORE 5.0% (Table 2). To balance measurement costs and estimation accuracy, it is recommended to choose SI = 5 cm within the depth of z = 30 cm, and SI = 10–20 cm at lower depths. If the simulation depth is large (e.g., L > 200 cm), SI = 30 cm may be used for z > 30 cm.
View larger version (24K):
[in this window]
[in a new window]
|
Fig. 2. Comparison of the estimated and theoretical distributions of the average source-sink term (SST) in the convection-dispersion equation (CDE) for (a) spatial intervals SI = 5–10, 5–20, 5–30 cm (Ex. 1a); and (b) time intervals T = 5, 10, 15, 30 d (Ex. 1b).
|
|
View this table:
[in this window]
[in a new window]
|
Table 2. Combinations of different treatments in the numerical experiments, the maximal absolute errors (MAE), the overall relative errors (ORE) of the estimated results, and the iteration number for Examples (Ex.) 1 to 3.
|
|
The estimated and theoretical SST distributions for T = 5, 10, 15, and 30 d with SI = 5–10 cm were compared in Fig. 2b. As time interval increased from 5 to 30 d, the ORE decreased from 3.3 to 0.2% (Table 2). However, the MAE increased from 0.8 x 10–5 to 4.4 x 10–5 mg cm–3 d–1 (Table 2) and occurred at z = 5 cm, which maybe resulted from the large change of soil water and NO3––N concentration near the soil surface. In general, the inverse method was stable when the time interval was chosen between 5 and 30 d.
Effect of Boundary Conditions and Layered Soils on the Estimated SST (Example 2)
In the numerical simulations, soil NO3––N concentrations at the lower boundary during the time interval were linearly interpolated using the measured NO3––N concentrations at t = 0 and t = T. Results from Example 1 (Ex. 1) showed that linear interpolation had little influence on the SST estimation. In this case, the change in NO3––N concentration range at the lower boundary from t = 0 to t = 10 d was only 0.0098 mg cm–3 for a silt loam (Fig. 1). To further investigate the effect of boundary conditions, a sandy soil was chosen (Table 1) with a larger NO3––N concentration change from 0.0250 to 0.0076 mg cm–3 within 10 d at the lower boundary. Compared with the theoretical SST distribution, the estimated results showed some disparity in the lower root zone. Nevertheless, the values of MAE and ORE were still within 3.0 x 10–5 mg cm–3 d–1 and 6%, respectively (Table 2). The estimation error resulted from the linear interpolation can be alleviated by increasing measurement frequency at the lower boundary, or choosing the simulating depth at the place where the change of NO3––N concentrations from 0 to T is small. Fortunately, for most soils, NO3––N concentrations below the rooting depth usually change slowly with time.
The above discussion was only limited to a single-layered, homogeneous soil. In reality, soil heterogeneity exists in almost all the fields and in different soil profiles with various textured soil layers in the root zone. In this example, we considered a two-layered soil, with the upper soil layer (0–90 cm) a sandy soil and 90–180 cm a silt loam (Table 1). In the iterative procedure, the "measured" points at T = 10 d were fitted into two continuous and smooth NO3––N concentration distributions for 0 to 90 and 90 to 180 cm (Fig. 1). The average distributions of theoretical and estimated SST were compared in Fig. 3
with a good agreement and the values of MAE and ORE less than 3.0 x 10– mg cm–3 d–1 and 2.2%, respectively (Table 2). The results showed that the inverse method was applicable for a layered soil profile.
View larger version (16K):
[in this window]
[in a new window]
|
Fig. 3. Comparison of the estimated and theoretical distributions of the average source-sink term (SST) in the convection-dispersion equation (CDE) with a time interval T = 10 d and a spatial interval SI = 5–10 cm for a two-layered soil profile (Ex. 2b).
|
|
Effect of Measurement Errors of NO3––N Concentration on the Estimated SST (Example 3)
Accuracy of the measured NO3––N concentration is important for estimating the SST and influenced by the precision of the instrument (p), as well as by soil heterogeneity or spatial variability. The effect of the precision of measurement instrument (p) on the estimated SST has been included in the discussions of the above examples. In most cases in the field, the actual measurement error (pa) should be larger than p due to soil spatial variability. Its effect on SST was discussed in this example. Similar to that using Eq. [29], the "measured" soil NO3––N concentrations were generated randomly by
 | [35] |
where RE represents the relative error for measured NO3––N concentration, and –RE rand(RE) RE. The values of RE were set as 5, 10, 15, 20, and 30%. The theoretical, the generated "measured" (for RE = 10, 20, and 30%) and the fitted NO3––N concentration distributions (for RE = 20% and 30%, only) at T = 10 d are shown in Fig. 4a
. The maximal differences between the "measured" and fitted NO3––N concentration, analogous to pa, were 2.17, 2.32, 3.27, 3.51, and 8.86 x10–3 mg cm–3 for RE = 5, 10, 15, 20, and 30%, respectively, which were much higher than p (= 0.9 x 10–4 mg cm–3). The theoretical and estimated SST distributions for different RE were compared in Fig. 4b. The values of MAE and ORE between the estimated and theoretical SST distributions increased with increasing RE (Table 2). However, the values of MAE and ORE were within 8.0 x 10–5 mg cm–3 d–1 and 10%, respectively, when the relative error RE was within 20%.
View larger version (23K):
[in this window]
[in a new window]
|
Fig. 4. (a) Theoretical, randomly generated "measured," and fitted soil NO3––N concentration distributions; and (b) comparison of the estimated and theoretical distributions of the average source-sink term (SST) in the convection-dispersion equation (CDE), at a time interval T = 10 d and a spatial interval SI = 5–10 cm for different relative error RE = 10%, 20% and 30% of measured nitrate concentrations (Ex. 3).
|
|
Estimating RNU Factor and Simulating NO3––N Transport with RNU Model (Example 4)
In the soil column experiment (Exp. 1), soil organic N was almost zero and there was not any ammonium nitrogen in the irrigated Hoagland solution. Therefore, the transformation processes of NO3––N in soils for ammonium nitrification, immobilization, and denitrification were neglected or assumed to be counteracted, thus the sum of the corresponding rates of nitrification, immobilization, and denitrification was set to be zero. Consequently, the transformation processes of NO3––N in Exp. 1 were predominated by root-nitrate-uptake (RNU). Using the average RWU rate estimated with the inverse method (Zuo and Zhang, 2002) and the measured soil NO3––N concentration distributions, we calculated the average RNU rate distributions u(z,T) (namely, the SST in the CDE in this case) during different irrigation periods according to the proposed procedure (Fig. 5a
).
View larger version (30K):
[in this window]
[in a new window]
|
Fig. 5. The average root-nitrate-uptake (RNU) rate distributions during different irrigation periods (DAP: days after planting) (a) estimated using the proposed inverse method; and (b) calculated with the established RNU models and the simulated NO3––N concentrations (Ex. 4).
|
|
The RNU factor in Eq. [19] was optimized as = 1.27 with the least-squares procedure, using the estimated average distributions of RWU and RNU rates during 10 to 16 DAP. Thereafter, soil NO3––N transport in the soil column for other irrigation periods (i.e., 16–22, 22–28, 28–34, 34–40, 40–46, 46–52, and 52–58 DAP) was simulated incorporating the established RNU model (Eq. [19] with = 1.27) and the estimated RWU rate distributions. Most of the root mean squared errors (RMSE1) between measured and simulated soil NO3––N concentration were <0.05 mg cm–3 with just one exception on 21.5 DAP (Table 3), which was resulted from the large discrepancy between measured and simulated values at the soil surface.
View this table:
[in this window]
[in a new window]
|
Table 3. The measured, estimated (with the inverse method) and simulated (with the RNU model and the simulated NO3––N concentrations) total nitrogen masses extracted by winter wheat, and the corresponding relative errors during different irrigation periods (DAP: days after planting); root mean squared errors (RMSE1) (mg cm–3) between measured and simulated NO3––N concentration distributions during different irrigation periods; and root mean squared errors (RMSE2) of the average RNU rate distributions (mg cm–3 d–1) between estimated with the inverse method and simulated with RNU model during different irrigation periods.
|
|
Based on the established RNU model (Eq. [19] with = 1.27), the simulated soil NO3––N concentration distributions and the estimated RWU rates, the average RNU rate distributions during other irrigation periods (i.e., 16–22, 22–28, 28–34, 34–40, 40–46, 46–52, and 52–58 DAP) were simulated with Eq. [34] and shown in Fig. 5b. The simulated average RNU rate distributions were compared with the estimated values by the inverse method (Fig. 5a) with the root mean squared errors (RMSE2) between them smaller than 3.0 x10–4 mg cm–3 d–1 (Table 3).
The total N mass extracted by winter wheat MRNU (mg) during each irrigation period was calculated by
 | [38] |
where A is the area of the column (cm2). On the other hand, the total N mass extracted by winter wheat during each period was obtained by measuring the N content and dry weight of collected canopies and roots in the experiment. The values of the total N mass extracted by winter wheat for measured, estimated (using the estimated average RNU rates), and simulated (using the simulated average RNU rates), and the relative errors during different periods are listed in Table 3. All the relative errors were smaller than 10%. The results showed that the hypothesis for neglecting the transformation processes of ammonium nitrification, immobilization and denitrification in Exp. 1 was reasonable and the two methods to calculate RNU were reliable and stable to guarantee the mass conservation.
|
CONCLUSIONS
|
|---|
An inverse method similar to that of Zuo and Zhang (2002) was applied to estimate the SST in NO3––N CDE using two successively measured soil NO3––N concentration distributions and other related information. Several numerical examples were designed to test the accuracy, stability, and effectiveness of the inverse method under various conditions, and possible effects of major factors on the numerical approach were addressed. Comparisons between the estimated and theoretical SST distributions showed a good agreement under different conditions when the relative measurement error of NO3––N concentration was within 20%.
A soil column experiment with winter wheat growth was performed to validate the inverse method and to show its applications in practice. The average RNU rate distributions during different irrigation periods were estimated using the inverse method, and the values during the first period were employed to optimize the RNU factor in the RNU model. Thereupon, soil NO3––N transport in the columns and the average RNU rate distributions during other periods were simulated with the established RNU model. The simulated NO3––N concentration distributions were in good consistency with the measured values, and the simulated average RNU rate distributions were also in a good agreement with the estimated values. The estimated or the simulated total N mass extracted by winter wheat during each period was well compared with the measured value with relative error < 10%.
|
ACKNOWLEDGMENTS
|
|---|
This study was supported partly by the National Natural Science Foundation of China (Grant No.: 50579072), the Programs for New Century Excellent Talents in Universities (NCET-04-0137) and for Changjiang Scholars and Innovative Research Team in Universities (IRT0412), Ministry of Education, China.
|
NOTES
|
|---|
Abbreviations: CDE, convection-dispersion equation; RNU, root-nitrate-uptake; RWU, root-water-uptake; SST, source-sink term.
Received for publication December 7, 2006.
|
REFERENCES
|
|---|
- Cabon, F., G. Cirard, and E. Ledoux. 1991. Modelling of the nitrogen cycle in farm land areas. Fert. Res. 27:161–169.[CrossRef]
- Carsel, R.F., and R.S. Parrish. 1988. Developing joint probability distributions of soil water retention characteristics. Water Resour. Res. 24:755–769.[CrossRef]
- Dalton, F.N., P.A.C. Raats, and W.R. Gardner. 1975. Uptake of water and solutes by plant roots. Agron. J. 67:334–339.[Abstract/Free Full Text]
- De Willigen, P. 1991. Nitrogen turnover in the soil-crop system: Comparison of fourteen simulation models. Fert. Res. 27:141–149.[CrossRef]
- Hansen, S., H.E. Jensen, N.E. Nielsen, and H. Svendsen. 1991. Simulation of nitrogen dynamics and biomass production in winter wheat using the Danish simulation model DAISY. Fert. Res. 27:245–259.[CrossRef]
- Hansen, S., H.E. Jensen, and M.J. Shaffer. 1995. Developments in modeling nitrogen transformations in soil. p. 83–107. In P.E. Bacon (ed.) Nitrogen fertilization in the environment. Marcel Dekker, New York.
- Hupet, F., S. Lambot, R.A. Feddes, J.C. van Dam, and M. Vanclooster. 2003. Estimation of root water uptake parameters by inverse modeling with soil water content data. Water Resour. Res. 39(11):1312 doi:10.1029/2003WR002046.[CrossRef]
- Lafolie, F. 1991. Modelling water flow, nitrogen transport and root uptake including physical non-equilibrium and optimization of the root water potential. Fert. Res. 27:215–231.[CrossRef]
- MeGechan, M.B., and L. Wu. 2001. A review of carbon and nitrogen processes in European soil nitrogen dynamics models. p. 103–171. In M.J. Shaffer et al. (ed.) Modeling carbon and nitrogen dynamics for soil management. Lewis publishers, CRC Press, Boca Raton, FL.
- Ma, L., and M.J. Shaffer. 2001. A review of carbon and nitrogen processes in nine U.S. soil nitrogen dynamics models. p. 55–102. In M.J. Shaffer et al. (ed.) Modeling carbon and Nitrogen dynamics for soil management. Lewis publishers, CRC Press, Boca Raton, FL.
- Millington, R.J., and J.M. Quirk. 1961. Permeability of porous solid. Trans. Faraday Soc. 57:1200–1207.[CrossRef]
- Musters, P.A.D., and W. Bouten. 2000. A method for identifying optimum strategies of measuring soil water contents for calibrating a root water uptake model. J. Hydrol. 227:273–286.[CrossRef]
- Priesack, E., S. Achatz, and R. Stenger. 2001. Parameterization of soil nitrogen transport models by use of laboratory and field data. p. 459–481. In M.J. Shaffer et al. (ed.) Modeling carbon and nitrogen dynamics for soil management. Lewis publishers, CRC Press, Boca Raton, FL.
- Romano, N., and A. Santini. 2002. Water retention and storage: Field–field water capacity. p. 723–729. In J.H. Dane and G.C. Topp (ed.) Methods of soil analysis. Part 4.SSSA Book Ser. No. 5. SSSA, Madison, WI.
- Schmied, B., K. Abbaspour, and R. Schulin. 2000. Inverse estimation of parameters in a nitrogen model using field data. Soil Sci. Soc. Am. J. 64:533–542.[Abstract/Free Full Text]
- 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]
- Shaffer, M.J., L. Ma, and S. Hansen. 2001. Introduction to simulation of carbon and nitrogen dynamics in soils. p. 1–26. In M.J. Shaffer et al. (ed.) Modeling carbon and nitrogen dynamics for soil management. Lewis publishers, CRC Press, Boca Raton, FL.
- Spalding, R.F., and M.E. Exner. 1993. Occurrence of nitrate in groundwater: A review. J. Environ. Qual. 22:392–402.[ISI]
- Valente, A.J.M., A.C.F. Ribeiro, V.M.M. Lobo, and A. Jimenez. 2004. Diffusion coefficients of lead (II) nitrate in nitric acid aqueous solutions at 298 K. J. Molecular Liquids 111:33–38.[CrossRef]
- 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.[ISI]
- van Genuchten, M.Th. 1987. A numerical model for water and solute movement in and below the root zone. Res. Rep. 121, USDAARS, U.S. Salinity Laboratory, Riverside, CA.
- Vasssilis, Z.A., and G.C.L. Wyseure. 1998. Modeling of water and nitrogen dynamics on an undisturbed soil and a restored soil after open-cast mining. Agric. Water Manage. 37:21–40.[CrossRef]
- Wu, J., R. Zhang, and S. Gui. 1999. Modeling soil water movement with water uptake by roots. Plant Soil 215:7–17.[CrossRef]
- Xu, S.L. 1995. Fitting and approximation, 8.1 Least squares curve fitting. p. 341–346. In Popular algorithms and their FORTRAN program codes. (In Chinese.) Tsinghua Univ. Press, Beijing, China.
- Zijlstra, J., and J.H. Dane. 1996. Identification of hydraulic parameters in a layered soil based on a quasi-Newton method. J. Hydrol. 181:233–250.[CrossRef]
- Zuo, Q., and R. Zhang. 2002. Estimating root-water-uptake using an inverse method. Soil Sci. 167:561–571.[CrossRef]
- Zuo, Q., L. Meng, and R. Zhang. 2004. Simulating soil water flow with root-water-uptake applying an inverse method. Soil Sci. 169:13–24.[CrossRef]
TOP
ABSTRACT
INTRODUCTION
NO3––N (Cabon et al., 1991)
