VZJ sign up for etocs
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
QUICK SEARCH:   [advanced]


     

Published online 24 August 2006
Published in Vadose Zone J 5:990-1004 (2006)
DOI: 10.2136/vzj2005.0124
© 2006 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in ISI Web of Science
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via ISI Web of Science (1)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Cannavo, P.
Right arrow Articles by Renault, P.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Cannavo, P.
Right arrow Articles by Renault, P.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Cannavo, P.
Right arrow Articles by Renault, P.
Related Collections
Right arrow Biogeochemical Processes
Right arrow Coupled Flow/Transport Models

ORIGINAL RESEARCH

Modeling Seasonal Variations in Carbon Dioxide and Nitrous Oxide in the Vadose Zone

P. Cannavoa,b, F. Lafoliea,*, B. Nicolardotc and P. Renaulta

a INRA, Unité Climat Sol Environnement, Bâtiment Sol, Domaine Saint-Paul, Site Agroparc, 84914 Avignon Cedex 9, France
b UAPV, Laboratoire d'hydrogéologie, 74 rue Louis Pasteur, 84029 Avignon Cedex 1, France
c INRA, Unité d'Agronomie de Laon-Reims-Mons, Centre de Recherches en Agronomie et Environnement, 2, esplanade Roland Garros, BP 224, 51686 Reims Cedex 2, France
* Corresponding author (lafolie{at}avignon.inra.fr)

Received 18 October 2005.



   ABSTRACT
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
INTRODUCTION
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Soil CO2 and N2O concentrations were simulated with a model predicting C and N transport in the vadose zone during a 7-mo field experiment, after maize (Zea mays L.) harvesting and incorporation of maize residues into the soil. The gas transport model was based on the dusty gas theory and combined with the PASTIS model. During the experiment, soil atmosphere (CO2 and N2O), soil solution (NO3– and dissolved organic carbon [DOC]), soil water content and temperature, and potential denitrifying and aerobic respiratory activities were measured in a 2.50-m-thick soil profile. Soil gas concentrations were correctly simulated even though the model did not simulate all the biological processes that produced N2O. Nitrous oxide concentration peaks after rain were slightly overestimated, as the WFPS (water-filled pore space) was not estimated accurately enough to predict local anoxic conditions. To model CO2 concentrations, account had to be taken of DOC adsorption onto soil mineral particles and of zymogenous biomass death during the period when the ground was frozen. The model satisfactorily simulated NO3– concentrations in the top soil profile, notably during major rainfall events, and maize residue dry matter loss during the experiment. The modeling of biological processes needs to be improved to provide a better simulation of C and N transport in the vadose zone. In particular, the use of WFPS was not sufficient to predict anoxic periods; simulations should improve if soil aggregate structure is also taken into account.

Abbreviations: AUB, autochthonous microbial biomass • CEL, cellulose • DGM, dusty gas model • DOC, dissolved organic carbon • DOY, Day of Year • FOM, fresh organic matter • HCE, hemicellulose • HOM, humified organic matter • RDM, rapidly decomposable material • WFPS, water-filled pore space • ZYB, zymogenous microbial biomass


   INTRODUCTION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
THE spatial distribution of several gases closely interacts with various soil biogeochemical cycles. For instance, O2, CO2, N2, N2O, NH3, H2, H2S, and CH4 are produced and consumed during various microbial processes, including aerobic and anaerobic microbial respiration, N2 symbiotic or nonsymbiotic fixation, fermentations, and acetogenic transformations (Sierra and Renault, 1994; Pelmont, 1993; Dassonville and Renault, 2002). Their distribution affects abiotic geochemistry (Stumm and Morgan, 1996) and microorganism behavior (Richards and Webster, 1999; Glinski et al., 2000; Young and Ritz, 2000). To our knowledge, few models that combine C and N biotransformations to water, solutes, and heat transport in soils (Ersahin and Karaman, 2001; Garnier et al., 2001; Neergaard et al., 2002; Garnier et al., 2003) include a description of gas transport, despite its importance and the feasibility of using simulated vs. experimental gas concentrations as a means to check and improve models of C and N behavior. Leffelaar (1988) developed a model at the scale of the aggregate to study the dynamics of partial anaerobiosis and denitrification, combining multinary gas transport with water and solute transport and N transformations.

Gas transport in soil results from both diffusion and convection (Rolston, 1986). Convection cannot be neglected (Elberling et al., 1998), even if air pressure gradients are often close to zero as a consequence of high air permeability values (Reinecke and Sleep, 2002). Convection results from (i) a CO2 diffusion coefficient lower than that of O2 (Marrero and Mason, 1972), (ii) respiratory quotients (CO2 production/O2 consumption ratio) differing from 1, (iii) atmospheric pressure changes (Massmann and Farrier, 1992), (iv) temperature fluctuations (Sierra and Renault, 1998), and (v) water flow in air-filled pores. Two main approaches have been used for gas diffusion modeling. One is based on the use of Fick's law to describe the diffusive flux of one gas regardless of the others (Baehr and Bruell, 1990). Fick's law can be regarded as a reduced form of the Stefan–Maxwell equations, which describe gas diffusion in a system containing more than two species. It has sometimes been combined with Darcy's law to account for convection (Webb and Pruess, 2003). The other approach is the DGM (dusty gas model), which gives a description of overall multicomponent gas diffusion in porous media. This may also be combined with Darcy's law to describe convection (Thorstenson and Pollock, 1989). The DGM is based on Stefan–Maxwell equations for representing diffusion due to gas molecules colliding with each other. In addition, it assimilates solid and liquid phases of the porous medium to immobile macromolecules (i.e., the dusts; Massmann and Farrier, 1992), and thus accounts for Knudsen diffusion (i.e., diffusion of gas molecules due to collisions with the solid particles). Modeling gas diffusion using Fick's law has been widely used and often shown to satisfactorily represent gas transport in the vadose zone (Simunek and Suarez, 1993; Moldrup et al., 1998, 2000, 2003). The DGM has been widely used in chemical engineering (diffusion in porous catalysts) and the first applications in environmental research were described some years ago in Thorstenson and Pollock (1989) and Massmann and Farrier (1992). Since then, many studies have been dedicated to comparing the model based on Fick's law and the DGM. The aim in these investigations was to quantify differences in predictions between the two approaches and to decide when Fick's law models could be used in place of the more complicated DGM (Jaynes and Rogowski, 1983; Leffelaar, 1987; Thorstenson and Pollock, 1989; Baehr and Bruell, 1990; Abriola et al., 1992; Webb and Pruess, 2003; Fen and Abriola, 2004). Modifications of the Fick's law model to account for viscous flows induced by unequimolar diffusion were also proposed (Freijer and Leffelaar, 1996; Webb and Pruess, 2003) and provide an attractive approximation of the DGM. In summary, Fick's law models do not differ markedly from the DGM when the diffusing gas is at trace concentration levels, when diffusion fluxes are equimolar, and when the porous medium air permeability is sufficiently high, but significant differences can be observed in other situations. For this reason, we selected the DGM.

The aim of this study was to better assess C and N behavior throughout the entire vadose zone by: (i) combining a gas transport module based on the dusty gas theory with the PASTIS (Predicting Agricultural Solute Transport in Soils) model initially proposed by Lafolie (1991) to describe water flow and N transport and later modified to improve the prediction of biological processes (Lafolie et al., 1997; Garnier et al., 2001; Garnier et al., 2003); and (ii) comparing the observed soil atmosphere CO2 and N2O and the soil NO3– content with the model's predictions, using experimental data obtained during a 7-mo period from the first 2-m depth of a fluvic hypercalcaric Cambisol.


   MATERIALS AND METHODS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Experimental Site
The selected experimental site (43°58'16'' N, 5°0'25'' E) is near Avignon (France), in a small region, labeled "vulnerable zone" with regard to major groundwater NO3– pollution due to agricultural practices (EEC 91–676 directive). The soil is a fluvic hypercalcaric Cambisol (FAO classification; see soil characteristics, Table 1). Four soil layers were defined: 0 to 30, 30 to 60, 60 to 100, and 100 to 200 cm. The water table remained below the 230-cm depth throughout the experiment. The soil was cultivated with maize in 2000 and 2001 and the experimental site has a Mediterranean climate. A 7-mo experiment was conducted from 25 Oct. 2001 (Day of Year [DOY] 297), after maize harvesting and soil tillage with maize residues incorporation (20-cm-depth plowing), to 6 May 2002 (DOY 126). The soil was bare throughout the experiment. More details about the experimental site are presented in Cannavo et al. (2004).


View this table:
[in this window]
[in a new window]
  		Table 1. Main characteristics of the soil profile.

 
In Situ Measurements
Soil Densities and Porosities
Bulk densities from each soil layer were estimated with soil cores (15-cm diameter by 7 cm high). Particle density was not measured, but a mean value of 2.7 g cm–3 was chosen. The density of the soil aggregates hereafter referenced as textural density was estimated from the clay content of the soil (Fiès and Stengel, 1981).

Soil Water Suction, Moisture, and Temperature
Soil water potential ({psi}) was monitored by first calibrating automatic tensiometers (SKT 850, SDEC, France) using a pressure calibrator (Druck DPI 602) at four temperatures (5, 15, 25, and 35°C) and five suctions (0, –2, –4, –6, and –7.5 m H2O) in a climatic chamber. Calibration parameters were estimated by fitting a multilinear regression to experimental data. The relative error was ~0.5%. The tensiometers were then installed every 20 cm from 20- to 200-cm depth to monitor soil water suction. They were connected to a datalogger (Campbell Scientific, CR10) and values were recorded every 20 min. The {psi}({theta}) relationship was determined from the volumetric soil water content ({theta}) by sampling soil at different dates and under different weather conditions (dry and rainy periods) at the same depths as the tensiometers. Additional points (water content and soil water potential) were obtained in the laboratory.

Relationships between {psi} (m H2O) and {theta} (m3 m–3) were obtained using the van Genuchten relationship (van Genuchten, 1980) fitted to experimental data obtained in the field and in the laboratory (Cannavo et al., 2004). Soil water suction was also recorded, at least three times a week, using Hg tensiometers (STM 2150, SDEC, France) installed at the same depths as the electronic tensiometers. These probes were used to control the automatic tensiometers.

Soil temperature was also measured at depths of 5, 20, 60, 100, and 170 cm (probe calibration has been described in Cannavo et al., 2004).

Gas and Soil Solution Sampling
Soil atmosphere probes were installed every 20 cm from 20- to 200-cm depth to monitor O2, N2, CO2, and N2O partial pressures. Gas samples were collected in each probe once a week, or more frequently after rain, and analyzed thereafter as described by Cannavo et al. (2004). The PVC (polyvinyl chloride) probes measured 25 mm in diameter and consisted of a gas cell and an extension piece. The gas cell was 120 mm high with a 3-mm-diameter hole and was linked to the soil surface by a 5-mm-diameter tube placed in the extension piece. The tube was sealed at the soil surface by a shutoff valve to avoid atmospheric contamination. The probes were installed by drilling the soil with a 45-mm-diameter auger and placing the probes at predetermined depths. Coarse sand was packed around the gas cell to favor gas diffusion, bentonite was poured onto the coarse sand to prevent gas atmosphere penetration along the tube, and native soil was placed above the bentonite. The tubes were first purged with a volume adapted to the length and diameter of the capillary tube connected to the gas cell, then soil atmosphere samples were collected in 3-mL evacuated Venoject tubes. The samples were analyzed in the laboratory no later than 1 h after sampling. Nitrous oxide was analyzed by ECD (electron capture detector) gas chromatography, and TCD (thermal conductivity detector) gas chromatography was used for O2 and CO2.

An area close to the place where soil water potential, temperature, and soil atmosphere composition were being monitored was specifically used to sample the soil for NO3– and DOC concentration measurements. Soil samples were taken with an auger in 30-cm-thick layers from 0- to 210-cm depth. Soil solution extraction and measurement protocols can be found in Cannavo et al. (2004).

In Situ Soil Maize Residue Decomposition
In situ maize residue decomposition was studied in 24 PVC columns (20 cm in diameter by 30 cm high). Maize residue dry matter (DM), excluding grain and roots, was measured before harvest and was 49.5 and 50.5% for the leaves and stems, respectively. The N contents for leaves and stems were 0.79 and 0.99% of DM, respectively. The amount of dry matter (leaves + stem) incorporated in soil was 835 g m–2. Two soil layers were reconstructed in the columns, the lower one 10 cm thick and the upper one 17 cm thick, in which the soil was mixed with coarsely chopped dry maize residues (4–5-cm pieces), 13 g DM of leaves and 13.2 g DM of stems. The amount of residue incorporated and the leaf/stem ratio respected the averaged field conditions. The columns were inserted in the field at 10-cm intervals on the day the maize was harvested (22 Oct. 2001, DOY 294). The crop residues were immediately incorporated by the farmer on the remainder of the field and sampling area. Columns were sampled 3, 6, 9, 12, 16, 20, 24, and 28 wk later and the amounts of dry residues present in the soil were measured. Crop residues >2 mm were separated by dry sieving. The soil was then sieved (200-µm mesh) under water. Crop residues in the 0.2- to 2-mm range were removed by flotation. The two maize residue size fractions thus obtained were dried at 75°C. They were then ground (80 µm) and the total C and N contents were measured using an elemental analyzer (Carlo Erba, NA 1500, Milan, Italy). The amount of residues in the soil, without maize residue incorporation, was also quantified by flotation and drying at 75°C.

Laboratory Measurements
Gas Diffusion Coefficient and Air Permeability Measurements
Intact soil cores were sampled in triplicate using stainless steel sampling rings (7-cm internal radius by 7-cm length) on 26 Feb. 2002 (DOY 57) from five different layers: 0- 30- and 30- to 60-cm layers, clay pockets (~1-m depth), sandy layers (below 1-m depth), and again in the 0- to 30-cm layer a few days later, after heavy rains had compacted the soil surface. Air trapping was prevented by saturating the samples for at least 48 h under vacuum. They were then placed on a suction table and water suctions of 10, 20, 30, 50, 70, and 100 cm H2O were successively applied. Equilibrium was achieved after 48 h. The diffusion coefficient for each water suction was estimated as described below. Table 2 summarizes the porosity and density of the soil samples.


View this table:
[in this window]
[in a new window]
  		Table 2. Densities and porosities of sampled soils.

 
The gas diffusion coefficient was estimated using a device based on the design by Ball et al. (1981), in which Kr-85 was used as the diffusing gas, and the method of diffusion coefficient estimation proposed by Bruckler et al. (1989), later improved by Cousin et al. (1999). Several models relating the gas diffusion coefficient to porosity were assessed (Buckingham, 1904; Penman, 1940; Millington and Quirk, 1961; Millington and Shearer, 1971; Moldrup et al., 2004).

Thereafter, the air permeability of each soil sample was measured with the same device. Several gas flows were imposed using a regulating flow meter (Brooks, 5850TR Series, 10–20 mL min–1, based on N2 equivalent, ±1.0% precision) and the air pressure gradient was then measured using a differential air pressure probe (Keller, Series 41, 1.0–100.0 kPa, ±0.2% precision), to check whether the flow was laminar.

Microbial Activities
The CO2 and N2O production and consumption rates in the soil profile were estimated by determining the potential microbial activities in the laboratory as described by Cannavo et al. (2002). Soil samples were collected with an auger at ~1-mo intervals for the 0- to 30-cm layer and every 3 mo for the other layers (30–60-, 60–100-, and 100–160-cm layers). The measured values are presented in Table 3.


View this table:
[in this window]
[in a new window]
  		Table 3. Aerobic respiratory activity (ARA, 10–2 mg O2 kg–1 dry wt. h–1), and potential and semipotential denitrifying activities (PDA and SPDA, respectively, 10–5 mg N2O–N kg–1 dry wt. h–1), in the soil profile at different days of the year (DOY).

 
Carbon and N mineralization of maize residues, initially ground to 1 mm, were studied during soil incubations under controlled conditions, as described by Trinsoutrot et al. (2000). Soil respiration and N mineralization were measured on soil samples (from the 0–30-cm layer) that were incubated with or without residues (4 g DM kg–1 dry soil) for 168 d at 28 ± 0.3°C, with the soil water suction potential maintained at 3.5 x 10–3 MPa. Analytical determinations for CO2 and mineral N are also detailed in Trinsoutrot et al. (2000). The residue biochemistry was determined using the method described by Linères and Djakovitch (1993) for extraction of the soluble fraction and that proposed by Van Soest and Wine (1967) for the hemicellulose, cellulose, and lignin fractions. The percentages of residue C were 14.5, 36, 40.3, and 9.2% for the soluble, hemicellulose, cellulose, and lignin fractions, respectively. The C/N ratios for these fractions were 5.9, 73.3, 381.1, and 58.3, respectively.

Model Description
The one-dimensional mechanistic PASTIS model (Predicting Agricultural Solute Transport in Soils) initially combined a description of the simulated flow of water (Richards equation), solutes (convection–dispersion equation), and heat (transport by water flow and diffusion) with some N transformations in the soil (Lafolie, 1991). It has been combined with a model that describes C and N biotransformations more extensively (Garnier et al., 2001), and with a module that describes denitrification (Hénault and Germon, 2000). In this study, we added a description of gas transport.

Numerical solutions for water and heat transport, and N transport and dynamics have been described in Lafolie et al. (1997). Parameter values for water, heat, and solute transport are presented in Table 4.


View this table:
[in this window]
[in a new window]
  		Table 4. Water, heat, and solute parameter values.

 
The module for C and N biotransformations in the soil uses first-order equations for describing biological processes (Garnier et al., 2001). Soil organic matter is divided into five organic pools: fresh organic matter (FOM), humified organic matter (HOM), soluble organic compounds (SOL), autochthonous microbial biomass (AUB) and zymogenous microbial biomass (ZYB). The FOM pool is composed of four biochemical fractions: rapidly decomposable material (RDM), hemicelluloses (HCE), cellulose (CEL) and lignin (LIG). The autochthonous biomass decomposes the HOM, whereas the zymogenous biomass decomposes the FOM and SOL.

In the denitrification model, the total denitrification rate is considered as the product of three functions depending on: (i) NO3– content, (i) WFPS, and (iii) soil temperature. The ratio N2O/(N2O + N2) was assumed constant and a 0.1 value was chosen throughout the experiment, as it was not measured.

Air Flow Equation and Dusty Gas Model
One-dimensional gas-phase convective flow in a porous medium is assumed to obey Darcy's law and is given by (Massmann, 1989)

Formula 1[1]
where q (kg m–2 s–1) is volumetric flux per unit area, {rho} (kg.m–3) is density of the gas phase, Bk (m2) is air permeability, µ (kg m–1 s–1) is gas-phase viscosity, z (m) is depth, and P (Pa) is pressure, with

Formula 2[2]
where C (mol m–3) is molar concentration, R is the gas constant (8.314 kg m2 s–2 mol–1 K–1), and T (K) is temperature.

The continuity equation is (Massmann, 1989)

Formula 3[3]
where t (s) is time and {varepsilon}a is air-filled pore space (m3 m–3). Combining with Eq. [1] gives

Formula 4[4]

The relationship between the density and the pressure of the gas phase is given by the Mariotte law:

Formula 5[5]
where wm (kg mol–1) is the molecular weight of the gas phase.

Combining Eq. [4] and [5] and assuming that temperature gradients can be neglected yields

Formula 6[6]

To account for gas production or consumption and also for diffusive molar fluxes arising from the DGM, the equation becomes

Formula 7[7]
where Si and Pi (mol m–3 s–1) are source and sink terms accounting for gas production or consumption, and NiD (mol m–2 s–1) are diffusive molar fluxes for gas i in a mixture of {nu} gases.

The last term in Eq. [7] accounts for the diffusive flux of the various species, which can be unequimolar and can create pressure gradients in the gas phase and consequently convection in the porous medium. The consumption and production rates (Si and Pi) can also be unequimolar and thus can create convection. Equation [7] is nonlinear due to the term PBk/µ in the derivative with respect to z. Massmann (1989) showed that a linear equation in which PBk/µ is replaced by P0Bk/µ can be used without noticeable loss of accuracy as long as the pressure at the boundary (Pl) and the initial pressure inside the porous medium (P0) are such that 1.2 > Pl/P0 > 0.8. This condition is entirely satisfied in our case. Thus, the equation we solved is

Formula 8[8]

This equation is solved with a no-flux boundary condition at the bottom of the profile due to the presence of a water table. The gas pressure at the soil surface is considered as equal to atmospheric pressure. Equation [8] is very similar to the Richards equation, except that it is linear whereas the Richards equation is nonlinear. The similarity allows employment of a finite differences scheme similar to the one used for the Richards equation. Briefly, a three-point, centered, finite difference scheme is used to approximate the spatial derivatives, the permeabilities being calculated at internodal locations with an arithmetic mean. A fully implicit scheme is used to treat the time derivative, and iterations are not required since the equation is linearized. Air-filled porosity is taken as the difference between total porosity and water content, the latter being calculated by the model for any time and location. Air permeability, Bk, depends on the air-filled porosity. The simple relation Bk = Bkmax({varepsilon}a/{varepsilon}T)2 was used (with {varepsilon}T = total porosity).

For a given mixture of {nu} gases, the DGM results in a set of {nu} equations (Thorstenson and Pollock, 1989). For a given species i (i = 1,...,{nu}), the equation for the isothermal case is as follows (Mason and Malinauskas, 1983):

Formula 9[9]
where NT (mol m–2 s–1) is the total molar flux, Xi is the molar fraction, Dij* (m2 s–1) is the binary diffusion coefficient of species i and j in the porous medium, and DiK (m2 s–1) is the Knudsen diffusion coefficient.

In the DGM, the total molar flux of a given species (NiT = NiD + XiNV) is used. The diffusive molar flux (NiD) and the viscous molar flux XiNV do not appear explicitly (Eq. [9]). To solve these equations in the unsaturated zone, one must know, at any point in the profile, (i) the gas pressure, (ii) the diffusion coefficients, (iii) the porous medium gas permeability and gas viscosity, and (iv) the temperature. Gas pressure is obtained from the solution of the air flow Eq. [8]. The gas diffusion coefficients and porous medium gas permeability depend on the water contents and temperatures calculated by the model. The Knudsen diffusion coefficient is calculated as in Thorstenson and Pollock (1989).

Most studies deal with binary or ternary gas mixtures. In these cases, expressions can be derived to express the total flux as a function of the molar fraction gradients of the different species (Alzaydi et al., 1978; Jaynes and Rogowski, 1983; Abriola et al., 1992; Massmann and Farrier, 1992; Webb and Pruess, 2003). Our purpose here is to propose a numerical solution for the DGM that can be used when more than three gas species are involved.

The set of equations for the DGM can be written in the form of a linear system of {nu} equations with {nu} unknowns NiT:

Formula 10[10]

At any point in the soil profile, assuming that the mole fractions, the gas-phase pressure, and the transport parameters are known, the solution of this system provides the total molar flux for all the species. For each species, a mass balance equation can be used to obtain the gas concentrations or mole fractions. For any given species i (i = 1,{nu}) with molar fraction Xi, the mass balance equation is as follows:

Formula 11[11]
where Cw,i (mol m–3) is the concentration of species i in the aqueous phase. It is assumed that the equilibrium between the aqueous and gaseous phases is instantaneous and that Henry's law applies. For any species i, the relationship between its partial pressure Pi and its concentration in the aqueous phase Cw,i is

Formula 18[18]
where KH,i is Henry's constant. Using this relationship and the one in Eq. [2], the concentration in the aqueous phase Cw,i can be related to the molar fraction Xi by Cw,i = XiCKH,i RT. Eq. [11] can now be written as follows:

Formula 12[12]
where Ri is a storage coefficient defined by:

Formula 19[19]

Note that the {nu} balance equations (Eq. [12]) are in fact coupled owing to Eq. [10]. If we assume that the fluxes (NiT) are known, Eq. [12] is easily solved with a finite difference approximation of the space derivative and a fully implicit scheme to advance in time. Given that the molar fractions must satisfy the relation

Formula 20[20]
one has only to solve the ({nu} – 1) transport equations.

One way to solve Eq. [10] and [12] is to derive algebraic expressions for the flux as a function of the various molar fraction gradients and hence to obtain a set of coupled partial differential equations. This is possible for binary or ternary gas mixtures only. To keep some flexibility in using the code, we decided to solve the set of equations by means of an iterative procedure. In this way we could use any number of gas species without being obliged to change the set of equations. The algorithm for solving the {nu} transport equations (Eq. [12]) coupled to the equations arising from the DGM (Eq. [10]) is as follows. Let Xj,it be the mole fraction of gas j at node i and time t. At each node of the finite difference grid, the set of equations (Eq. [10]) can be solved to provide a first estimate of gas flux Nj,iT at time t + dt. Next, these estimates are used to solve the mass balance equations (Eq. [11]). This step provides new estimates Xj,ik of gas mole fractions Xj,it+dt,which in turn are used with the dusty gas set of equations (Eq. [10]) to get new estimates of the gas flux. This iterative procedure is repeated until the following condition is met:

Formula 13[13]
where k is the iteration index, j is the gas species index, i is the finite difference grid index, and N is the number of nodes in the finite difference grid.

In our calculations, the convergence criterion {varepsilon} is 10–5. When convergence cannot be achieved within a prescribed number of iterations, the time step is reduced and the described numerical procedure is repeated.

Calibration Procedures
As a preliminary step, the water transport parameters were calibrated using data from contrasting (wet or dry) periods of the experiment. Initial conditions were provided by field water potential measurements (23 Oct. 2001, DOY 295) at 20-cm depth intervals from 20- to 220-cm depth. The soil profile was divided into five soil layers: 0 to 30, 30 to 60, 60 to 110, 110 to 140, and 140 to 220 cm. The unsaturated hydraulic conductivity K({theta}) parameters (van Genuchten, 1980) were calibrated for each soil layer to satisfactorily reproduce the measured soil water suctions (electronic tensiometers and Hg pressure gauge tensiometers).

Heat flux parameters were calibrated after water transport calibration. Initial conditions were provided by field soil temperature measurements at 5-, 20-, 60-, 100-, and 170-cm depth. During a 20-d period with snow cover (13 Dec. 2001–2 Jan. 2002, DOY 347–2), a –2°C soil surface temperature value was set in the model, because the model was not adapted for simulating energy balance with a snow pack on the soil surface.

The initial conditions for the solute transport model were the NO3––N and DOC soil solution concentrations measured on 23 Oct. 2001 (DOY 295) in the seven soil layers (i.e., 0–30- to 180–210-cm layers).

The biological parameters were estimated using the incubations described above. The C and N kinetics of the control soil were used to set the values of the parameters linked to the HOM and AUB pools. The C and N kinetics of the residue-amended soil were then used to set the values of the parameters related to the FOM and ZYB pools (Table 5). Soil respiration in the 0- to 30-cm layer was simulated by the C and N transformations module. Potential aerobic respiration and denitrifying activities were used for the deepest soil layers, corrected for soil temperature and water content (Garnier et al., 2001). In addition, the potential aerobic respiration below 30 cm was too high for CO2 concentration to be correctly simulated, probably due to physical aggregate disruption by sieving increasing the C bioavailability (Ashman et al., 2003). We therefore multiplied the measured respiration in the 30- to 60-cm layer by a factor of 0.23 and that of deeper layers by 0.01. These reduction factors were adjusted to ensure that the model gave a reasonable reproduction of the CO2 concentrations in the deepest layers.


View this table:
[in this window]
[in a new window]
  		Table 5. Values of biological parameters.

 
Denitrification was considered to be effective when the WFPS exceeded 70% and we set the {alpha} power at 2 in the soil WFPS function (Hénault and Germon, 2000). The Q10 value for the soil temperature function was set at 2. In addition, continuous N2O emissions occurred at the soil surface—generally under aerobic conditions—throughout the entire experiment and were attributed to continuous aerobic denitrification, based on use of the 15N isotope (Cannavo et al., 2004). Thus, a residual denitrifying activity was set in the model.

The calibrated model was then used to simulate the entire experiment (23 Oct. 2001–6 May 2002, DOY 295–126).


   RESULTS AND DISCUSSION
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
CONCLUSIONS
 REFERENCES
 
Soil Gas Diffusion Coefficient and Air Permeability
The measured Di*/Di ratios as a function of air-filled porosity are presented in Fig. 1. Di* and Di are the diffusion coefficients in the soil and in the air, respectively. Values for Di*/Di ranged from 0.0014 to 0.1432. Comparison of the experimental data with several models (without clay pocket experimental data), showed that the Buckingham (1904) model seemed to give the best fit. The maximal relative error was 55%. We therefore tried to reduce this error by using a power function that was fitted to the observations, with the exception of the clay pocket data. We obtained the relation

Formula 14[14]
This relation took the general form Di*/Di = A{varepsilon}aB, as already proposed by Buckingham (1904), Penman (1940), and Marshall (1959). The maximal relative error was still 55%, as the A and B values were quite close to those of Buckingham (1904). The fitted curve was used in the gas transport model to estimate the binary gas diffusion coefficient as a function of air-filled porosity. The differential pressure probe did not detect any pressure gradient in the soil samples when the air permeability coefficient Bk was estimated. This was due to a high Bk, which allowed gas fluxes without any measurable pressure gradient. A 10–10 m2 value was adopted for air permeability in the model, as indicated in Nazaroff and Nero (1988).


Figure 1
View larger version (19K):
[in this window]
[in a new window]
  		Fig. 1. Experimental Di*/Di ratio (diffusion coefficients in the soil and air, respectively) as a function of the soil's air-filled porosity. Comparison with simulated values from models by Buckingham (1904), Penman (1940), Millington and Quirk (1961), Millington and Shearer (1971), and Moldrup et al. (2004).

 
Physical Simulations with the PASTIS Model
Rainfall events during the experiment are presented in Fig. 2 . At the 20-cm depth (Fig. 3A ), the model satisfactorily simulated the different drying and wetting periods. Soil water suction was also satisfactorily simulated deeper in the soil profile. There were some significant differences between experimental and simulated data, particularly in April (DOY 102) at 40 cm depth (Fig. 3B) and in February (DOY 52) at the 100-cm depth (Fig. 3C). Nevertheless, these differences were <0.4 m H2O suction. Below the 140-cm depth, the experimental and simulated data were very close. Soil volumetric water content was also simulated by the model and compared with soil volumetric water content measured from soil regularly sampled at various places in the field. The model correctly simulated the soil volumetric water content except during the cold winter period (results not shown). This difference between measured and calculated water contents during the cold winter period was not really a problem given that microbial activity was negligible during this period.


Figure 2
View larger version (8K):
[in this window]
[in a new window]
  		Fig. 2. Rainfall events, as amount of precipitation (P), during the experiment.

 

Figure 3
View larger version (15K):
[in this window]
[in a new window]
  		Fig. 3. Experimental and simulated soil water suction at (A) 20-, (B) 40-, and (C) 100-cm depths.

 
Whatever the depth, the model simulated soil temperature satisfactorily (Fig. 4 ). The –2°C soil surface temperature value set in the model for the 20-d snow cover period allowed correct reproduction of soil temperature during this period. The differences between experimental and simulated data were generally <1.5°C. It was important to minimize the soil temperature error simulation to correctly simulate seasonal variations in soil microbial activity, which is highly sensitive to soil temperature (Grundmann et al., 1995; Mergel et al., 2001).


Figure 4
View larger version (17K):
[in this window]
[in a new window]
  		Fig. 4. Experimental and simulated soil temperature at (A) 20-, (B) 60-, and (C) 170-cm depths.

 
Finally, soil water content and water suction, in addition to soil temperature, were satisfactorily simulated by the model, and this enabled us to model gas transport and microbial activities.

Biological Simulations with the PASTIS Model
Diffusion models of gas transport in soils have to account for (i) O2 that is consumed by aerobic respiration, and discriminate between aerobic and anaerobic activities, (ii) CO2 that is produced through organic C mineralization and can affect microbial activities either directly through its assimilation by autotrophs or indirectly through changes in soil pH (Suchomel et al., 1990), and (iii) N2O that can be produced by nitrification (Stevens and Laughlin, 1998), denitrification (Nobre et al., 2001), and the dissimilatory reduction of NO3– to NH4+ (Appello and Postma, 1994), the last pathway usually being negligible (Knowles, 1982).

Simulation of Carbon Dioxide Concentrations in the Soil Profile
Carbon dioxide concentrations were simulated by applying the CANTIS model with measured biological parameters (Table 5) in the 0- to 30-cm layer, and potential aerobic respiratory activity in deeper layers (Table 3). The model did not satisfactorily simulate in situ CO2 concentrations at the 20-cm depth (Fig. 5A ). It initially underestimated CO2 concentration up to 1 Nov. 2001 (DOY 305) and overestimated concentrations from then until the beginning of the cold winter period (14 Dec. 2002, DOY 348). This was an overestimation of the burst of zymogenous biomass stimulated by the incorporation of fresh organic matter into the soil. A second overestimation was observed after the cold winter period (19 Jan. 2002, DOY 19), when the soil temperature increased significantly (Fig. 4). In this case, the zymogenous biomass, stimulated by soil warming, continued to decompose the maize residues, but here again the stimulation was overestimated. After that, the model underestimated CO2 concentration because it had been assumed that the maize residue would become less available and so decreased zymogenous biomass activity. The behavior observed at the 20-cm depth was also present at deeper layers (Fig. 5B, 5C, and 5D).


Figure 5
View larger version (24K):
[in this window]
[in a new window]
  		Fig. 5. Seasonal variations in CO2 concentrations in the soil profile at (A) 20-, (B) 40-, (C) 100-, and (D) 200-cm depths.

 
Two main explanations may be proposed to explain the model's discrepancies. One is DOC adsorption, the other is zymogenous biomass mortality in winter. Soil DOC adsorption onto soil mineral particles was not taken into account and its incorporation in the model should reduce organic substrate bioavailability for ZYB. The DOC adsorption onto mineral particles in soil is well known and is significant (Jardine et al., 1989; Michalzik et al., 1998; Kalbitz et al., 2000). A partitioning constant Kd of 0.001 m3 kg–1 was therefore introduced into the model. This approximate value was based on measurement of the Kd coefficient in several other studies (Moore et al., 1992; Kaiser et al., 2000). Taking physical adsorption into account improved the simulated CO2 concentrations in the soil profile before the winter period (Fig. 6A and 6B). A slight overestimation was still observed in early December 2001 (DOY 337), as if adsorption delayed and reduced the ZYB burst when maize residue was incorporated in the soil. It seemed to be necessary to limit organic substrate bioavailability to better simulate CO2 concentration; however, this did not improve CO2 concentration simulation after the freeze period when this overestimation was more pronounced than in the previous simulation. This can be explained by the fact that less FOM was degraded before winter, thus increasing the degradation of maize residues, and consequently CO2 production, after the winter.


Figure 6
View larger version (13K):
[in this window]
[in a new window]
  		Fig. 6. Simulation of seasonal variations in CO2 concentration in the soil profile at (A) 20- and (B) 40-cm depths, taking into account the adsorption of soluble organic C to soil mineral particles.

 
This led us to test the second assumption. The soil was frozen for almost 3 wk during the cold winter period, from 17 December to 6 January (DOY 351–6). It has already been shown that frozen soils decrease soil biomass (Schimel and Clein, 1996). Thus, a process description of enhanced death of the ZYB fraction during the cold period was added to the model. An exponential function describing ZYB death in the frozen part of the soil was applied for this 3-wk period, as follows:

Formula 15[15]
where CZYB (g C m–2) is the zymogenous biomass C content, j (date) is time, and {varphi} is a calibration parameter.

Several {varphi} values for CO2 concentration simulations were tested, taking DOC adsorption into account (results not shown). When the {varphi} value was increased, the CO2 concentration overestimation occurred later in the experiment and was reduced. A value close to 0.25 seemed to provide a good agreement between measured and calculated CO2 concentrations. Such a value implied the death of a large part of the ZYB pool during the cold period; however, after this cold period, the ZYB was able to grow rapidly and its activity fitted better with in situ CO2 concentration measurements.

Thus, both DOC adsorption and ZYB death in winter had to be taken into account in the model to better fit the experimental CO2 concentration data. A new simulation was then performed, using a {varphi} calibration parameter value of 0.25 for ZYB death and a Kd value of 1 for DOC adsorption (Fig. 7 ). The model now satisfactorily simulated the experimental data at all depths. At the 20-cm depth, the model showed CO2 concentration peaks during rain events, although it sometimes overestimated them. During the first month of the experiment, CO2 concentration at the 20-cm depth was still underestimated, whereas the model satisfactorily predicted the CO2 concentration at the 40-cm depth and deeper.


Figure 7
View larger version (23K):
[in this window]
[in a new window]
  		Fig. 7. Simulations of soil CO2 concentration at (A) 20-, (B) 40-, (C) 100-, and (D) 200-cm depths taking into account DOC adsorption and ZYB fraction death during winter.

 
The cumulative CO2 released to the atmosphere, calculated by the model, for the 7-mo experiment was 1540 kg CO2–C ha–1. This is within the range reported by Bowden et al. (2000) and Kessavalou et al. (1998), but 1.3 times lower than that of Osozawa and Hasegawa (1995), who worked on arable soils with higher soil temperatures. Carbon dioxide production in the 0- to 20-cm layer came mainly from decomposition of FOM, with less CO2 derived from AUB. The CO2 produced by the ZYB was directly correlated with changes in its size during the experiment. With increasing depth, the HOM pool became dominant in producing CO2, as no FOM was present below the 20-cm depth (data not shown). The calculated cumulative CO2 produced by FOM and HOM decomposition in the 0- to 20-cm layer was 85.8 and 23.8 g C m–2, respectively.

In the 0- to 60-cm layer, the cumulative CO2 produced by FOM, HOM, ZYB, and AUB taken together throughout the experiment was 1640 kg C ha–1. The CO2 produced at depths below 60 cm was 245 kg C ha–1. Thus, the microbial production of CO2 below the 60-cm depth was significant, as it amounted to 13% of CO2 production from the soil profile as a whole. This clearly suggests that microbial activity in the deep vadose zone must be taken into account in models predicting NO3– leaching.

The simulations presented below take into account both DOC adsorption and ZYB death during the cold winter period.

Carbon Mineralization: Maize Residue Decomposition
Maize residue remaining in the soil was obtained by adding the >2- and 0.2- to 2-mm size fractions together. Decomposition displayed three distinct stages (Fig. 8 ). From residues incorporation to 17 Dec. 2001 (DOY 293–351), 44.3% of the initial C maize residue disappeared. From 17 Dec. 2001 to 4 Feb. 2002 (DOY 351–35), the soil temperature was very low (Fig. 4), which slowed residue decomposition and no significant decomposition was observed. From 4 Feb. to 3 May 2002 (DOY 35–123), as soil temperature increased, the maize residues decomposed but at a slower rate than that observed during the first period. During the 7-mo period as a whole, 68.5% of the maize residue C was lost.


Figure 8
View larger version (11K):
[in this window]
[in a new window]
  		Fig. 8. Cumulative maize residue C decomposition during the experiment.

 
The model also simulated these three stages satisfactorily. It predicted decomposition of 77.6% of the maize residue C, which was comparable to the in situ measurement; however, it underestimated the first stage (estimating 33.3% of residue C decomposed) and overestimated the third stage. We were able to limit the underestimation of the first stage by increasing the initial quantity of the ZYB biomass to 10 g C m–2. This value was based on the maize residue fragment size before its incorporation into the soil, and its decomposition factor KMZ (Garnier et al., 2001; Garnier et al., 2003). The larger the residue fragments, the higher the KMZ and initial ZYB values.

Concerning the decomposition of FOM biochemical pools simulated by the model, the HCE and CEL pools were predominant in the composition of the maize residue and were mainly responsible for dry matter loss. The RDM also contributed to dry matter loss. The lignin pool did not significantly contribute to dry matter loss as this fraction is more resistant to biodegradation (Sanger et al., 1997): only 4% was decomposed during the experiment. At the end of the experiment, almost all the RDM pool, 85% of the HCE pool, and 94% of the CEL pool had decomposed. The HOM pool decreased slightly from 5752 to 5651 g C m–2 during the experiment.

Simulation of Nitrous Oxide Concentrations in the Soil Profile
Seasonal N2O concentrations in the soil profile have already been presented and characterized in Cannavo et al. (2004). Below a depth of 20 cm, the model gave a good representation of N2O concentration, with no overestimations (Fig. 9B , 9C, and 9D). The model reproduced N2O increases at the 20-cm depth (Fig. 9A) during major rainfall events (24 Jan. [DOY 24], 18 Feb. [DOY 49], and 6 Mar. 2002 [DOY 65]). It overestimated the highest in situ N2O peak (18 February) and predicted a N2O peak on 6 March (a rainy day), a day when gas was not sampled. Gas samples were taken 3 d later and showed an increased N2O concentration. This suggests that the model was not able to describe the phase difference between the rainfall event and the N2O emission, as anoxic conditions progressively settle after a rainfall event (slow O2 consumption; Renault and Sierra, 1994). The heaviest rain (55 mm) fell on 2 May (DOY 122), after a dry period of 3 wk, and the model simulated an N2O peak that was not detected in situ. The anoxic conditions cannot have lasted long enough for the denitrifying microorganisms to synthesize denitrifying enzymes. The model uses a denitrification function defined by Hénault and Germon (2000), which is simply based on a WFPS threshold. Thus, it does not include biological processes and this suggests that the WFPS function should depend on soil structure and its respiration (Renault and Sierra, 1994). These limits were observed on two occasions. First, the model did not simulate the increase in N2O concentration following incorporation of maize residue into the soil, which favored denitrification by enhancing aerobic respiratory activity and creating local anoxic conditions (Cannavo et al., 2004). Second, the model failed to simulate increased N2O concentrations during daily freeze–thaw cycles. Thus, the denitrifying model gave a good simulation of N2O peaks during rainfall events, even though denitrifying enzymes may not have been synthesized earlier. The residual denitrifying activity in the surface layer was set at 6.5 10–5 mg N2O–N kg–1 dry weight h–1, which was equivalent to a daily output of N2O by denitrification of 4.2 g N2O–N ha–1.


Figure 9
View larger version (19K):
[in this window]
[in a new window]
  		Fig. 9. Seasonal variations in N2O concentration in the soil profile at (A) 20-, (B) 40-, (C) 100-, and (D) 200-cm depths.

 
According to the model, the cumulative amount of N2O gas released to the atmosphere during the 7-mo experiment was 0.9 kg N2O–N ha–1. This is seven times higher than the value reported by Bowden et al. (2000), who worked with a temperate forest soil, and 2.3 times higher than that of Kessavalou et al. (1998), who worked with a tilled arable soil.

Simulation of Nitrogen Mineralization and Nitrate Concentration in the Soil Profile
The net N mineralization during the 7-mo experiment was 215 kg N ha–1. This was 3.6 times higher than values published by Mary et al. (1999), who studied straw decomposition after disk plowing. Daily net N mineralization in the 0- to 20-cm layer was always positive and varied between 0.007 and 3.439 kg N ha–1.

The amount of NO3– in situ did not present significant temporal variations in the soil layers below 30 to 60 cm; however, considerable variations were observed in the 0- to 30-cm layer, and these were correctly reproduced by the model (Fig. 10 ). Before 16 Jan. 2002 (DOY 16), the model generally overestimated the in situ data. There was no major rainfall event during this period and the model simulated a slightly increasing NO3– amount in this layer. In situ NO3– concentration measurements were measured in triplicate, with a mean 65% variation coefficient. Thus, a large standard error measurement does not explain all the NO3– concentration variations observed before16 Jan. 2002 (DOY 16). Major rainfall events occurred after this date, however, and the model correctly simulated the decrease in NO3– concentration until 14 Feb. 2002 (DOY 45) and the increase in NO3– concentration thereafter. The model simulated NO3– leaching from the top of the 30- to 60-cm layer during the January to March period (DOY 22–82), and was relatively consistent with the measurements (Fig. 10B). The model simulated some NO3– leaching to lower layers after the heavy February rains (Fig. 10C and 10D); however, the NO3– concentration measured after this event did not increase below the 60-cm depth. This discrepancy could not be due to a fault in the model, as the PASTIS model had already correctly simulated NO3– concentrations in the deep vadose zone in earlier studies (Garnier et al., 2001, 2003). In addition, there was no water leaching to the groundwater. It would therefore appear that some biological processes had not been taken into account. Denitrification in the 0- to 60-cm layer might prevent NO3– leaching; however, this would lead to higher N2O concentrations than the measured values, since almost 300 kg N ha–1 was leached from the 0- to 30-cm layer during the February rains. Denitrification with N2 as the terminal product was a possibility, since no N2O increase was observed in the profile. Another assumption could be NO3– storage within the microbial biomass that could later be partially denitrified. Ellis et al. (1996) observed the disappearance of 330 kg NO3––N ha–1 after a 2-h incubation under anaerobic conditions, which could be stored by both denitrifiers and microorganisms capable of autotrophic denitrification, and led to an underestimation of potential gaseous N loss.


Figure 10
View larger version (23K):
[in this window]
[in a new window]
  		Fig. 10. Nitrate amounts in the (A) 0- to 30-, (B) 30- to 60-, (C) 90- to 120-, and (D) 120- to 150-cm layers during the experiment. Bars represent standard deviations.

 

   CONCLUSIONS
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
REFERENCES
 
We incorporated a gas transport module that combined the DGM describing gas diffusion with Darcy's law. Model simulations for water, heat, NO3–, CO2, and N2O were compared with measurements obtained from a 7-mo field experiment performed from autumn to spring after crop residue incorporation into the soil. The simulations of water flow and heat transport were in good agreement with field measurements. Carbon dioxide production by the biomass decomposing crop residues seemed to be overestimated unless the accessibility of the substrate (DOC adsorption on solids) was reduced and ZYB death during the frost period was taken into account. The model correctly reproduced the decomposition of crop residues. Our simple approach to N2O emissions did not allow us to accurately describe highly localized anoxic conditions favoring denitrification. Nitrate amounts in the 0- to 30- and 30- to 60-cm soil layers were correctly reproduced. The amounts of NO3– measured deeper in the profile were much lower than in the top layers and did not vary, whereas the model predicted an increase in NO3– in deeper layers due to leaching from the top layers. The fact that large amounts of NO3– seemed to have disappeared from the 90- to 120-cm layer and beneath could not be explained satisfactorily. Denitrification with N2 as the terminal product is a possibility since no N2O increase was observed in the profile, which suggests that N2O production was limited to the first 30 cm. Simulations showed that gas production in the deeper layers had to be considered to correctly reproduce gas concentration in the soil atmosphere and fluxes to the atmosphere.


   ACKNOWLEDGMENTS
 
This work was carried out in the Climate Soil and Environment Unit, INRA, Avignon (France), and the Laon-Reims-Mons Agronomy Unit, INRA, Reims (France). We wish to thank D. Warwick for reviewing the English version of the manuscript, M. Rabiet for gas diffusion measurements, T. Flacher, farmer, for lending a part of his field, M. Daniel and R. Simler for design and construction of gas probes, M. Bourlet for soil analyses, and Dr. J-P. Legros for helping with the FAO classification.


   REFERENCES
TOP
 EXECUTIVE SUMMARY
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES