JOURNAL OF GEOPHYSICAL RESEARCH, VOL. 116, G03040, doi:10.1029/2011JG001644, 2011

Impact of transient soil water simulation to estimated nitrogenleaching and emission at high‐ and low‐deposition forest sites
in Southern California
Fengming Yuan,1 Thomas Meixner,2 Mark E. Fenn,3 and Jirka Šimůnek4
Received 6 January 2011; revised 13 July 2011; accepted 16 July 2011; published 29 September 2011.

[1] Soil water dynamics and drainage are key abiotic factors controlling losses ofatmospherically deposited N in Southern California. In this paper soil N leaching and trace
gaseous emissions simulated by the DAYCENT biogeochemical model using its originalsemi‐dynamic water flow module were compared to that coupled with a finite elementtransient water flow module (HYDRUS), for two mixed conifer forests with annual
deposition rates of about 70 and 9 kg N ha−1, in the San Bernardino National Forest.
Numerical solution of the Richards equation implemented in HYDRUS water module
could improve response of surface soil water dynamics to precipitation pattern, comparedto the original, and consequently it resulted in annual N gaseous emission loss about
1.5 ∼ 2 times higher. While the two flow modules predicted similar amounts of annual
water drainage, the HYDRUS water module simulated more frequent, but smaller
drainage fluxes, which favors soil mineralization and downward transport. In normal
precipitation years, annual leaching losses predicted by the HYDRUS coupled DAYCENT
model was about 5–18 kg N ha−1 higher due to different temporal patterns of dailywater drainage. In dry and wet years, leaching losses were similar. Our analysis suggeststhat it is necessary to fully capture dynamics of transient water flow (e.g., by numericallysolving the transient Richards equation) in order to adequately estimate soil N gaseousemissions and N transport and thus leaching, although it requires more computational
resources while the uncertainty in model improvement is still large due to lackof measurements.
Citation: Yuan, F., T. Meixner, M. E. Fenn, and J. Šimůnek (2011), Impact of transient soil water simulation to estimatednitrogen leaching and emission at high‐ and low‐deposition forest sites in Southern California, J. Geophys. Res., 116, G03040,
doi:10.1029/2011JG001644.

1. Introduction[2] In southern California, natural ecosystems in regions
downwind of Los Angeles and nearby agricultural areas sufferfrom heavily polluted air, which produces the “N saturation”
phenomenon due to the highest atmospheric N deposition
rates in the United States [Fenn et al., 1996, 2003]. Besides
field investigations on the fate of anthropogenic N additions inthese ecosystems, the modeling approach has been used [e.g.,
Arbaugh et al., 1999; Fenn et al., 2008]. The modeling adopted
for these forest ecosystems is similar to applications in agricultural systems [Chen et al., 2008; Del Grosso et al., 2005,
1
Institute of Arctic Biology, University of Alaska Fairbanks, Fairbanks,
Alaska, USA.
2
Department of Hydrology and Water Resources, University of
Arizona, Tucson, Arizona, USA.
3
Pacific Southwest Research Station, USDA Forest Service, Riverside,
California, USA.
4
Department of Environmental Sciences, University of California,
Riverside, California, USA.

Copyright 2011 by the American Geophysical Union.
0148‐0227/11/2011JG001644

2006, 2009], and grasslands [Parton et al., 2001]. Although
this sort of modeling approach has been widely used toincrease our understanding of the fate of N‐fertilizer applications in these managed systems, it is less known how it performs when applied to natural ecosystems, such as forests
under continuous atmospheric N deposition.
[3] It is well known that the movement and losses of excess
soil N depend upon soil water dynamics, because water flow is
the primary soil solute transport pathway [Feyen et al., 1998].
Additionally, soil water‐filled pore space (WFPS) is a primaryfactor controlling N trace gas production [Dobbie and Smith,2003] and soil water contents near saturation favor primary
N trace gas emissions [e.g., del Prado et al., 2006], althoughaerobic conditions may produce non‐ignorable trace gaseous
N from soil nitrification in these systems. Therefore a reliable
description of soil water content and fluxes in both wetting anddrying processes is necessary in order to adequately represent
the N cycle in soil‐plant systems. What, however, are theappropriate spatial and temporal scales that need to be considered to sufficiently represent important water processes inbiogeochemical models for various purposes [Chen et al.,2008]? This problem is especially critical for the estimation

G030401 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

of gaseous emissions, which are highly variable in time; andeven pulsed for sudden wetting or drying events during andafter rainfall or snowmelt [e.g., Frolking et al., 1998; Smithet al., 1998; Stolk et al., 2011].
[4] It has been shown that dynamic soil water models basedon the numerical solution of the Richards equation can adequately simulate transient water flow and thus soil water
dynamics [de Jong and Bootsma, 1996; Ranatunga et al.,2008; Scanlon et al., 2002], and these approaches have been
extensively applied in land surface models [Shao andHenderson‐Sellers, 1996; Lee and Abriola, 1999] and agricultural and environmental modeling systems [Ranatunga
et al., 2008; Gu and Riley, 2010; Gu et al., 2009]. The disadvantage of the Richards equation based water flow models isthat they require knowledge of unsaturated soil hydraulicproperties, which are rarely available, and daily or even shortertime steps and fine soil horizon discretization for numericalsolution of the Richards equation [Lee and Abriola, 1999].
These approaches require significant computational resources
and more data inputs as well, and thus are not preferred incomplex ecosystem models for regional or global applicationor for effective and exhaustive analysis of soil, water, plant,
and climatic interactions on multiple biogeochemical processes [de Jong and Bootsma, 1996; Rodriguez‐Iturbe et al.,2007]. On the other hand, the commonly used tipping bucket
models are too simple to represent soil water dynamics [deJong and Bootsma, 1996; Shao and Henderson‐Sellers,
1996]. Thus semi‐dynamic water models were extensively
developed to allow for a robust representation of soil water
dynamics and fast computation in such circumstances, e.g.,
DSSAT [Jones et al., 2003; Ritchie, 1998], DAYCENT[Parton et al., 1998], and EU‐Rotate_N [Doltra and Muñoz,
2010]. These models adopted the field capacity concept froma tipping‐bucket model, but allow multiple soil layers withmore realistic and mechanical modifications regarding infiltration, soil internal water flow, bottom drainage, and evapotranspiration in order to catch the dynamic features of
infiltration (during rainfall or snowmelt) and redistribution ofsoil water. Pachepsky et al. [2006] and Doltra and Muñoz[2010] demonstrated that these water flow models could
achieve similar accuracy as Richards equation based models
for simulation of soil water fluxes. However, the concept offield capacity is not only ambiguous in its definition for the
purpose of mathematical modeling [Twarakavi et al., 2009]; italso essentially states that soil water above this point would
drain from the soil, which is true when soil water flow reaches
steady state in a period of a few days upon soil internal flow
[Ritchie 1998]. It implies that soils are either fully saturated
such as after an intense rainfall, or at or below field capacity[e.g., Doltra and Muñoz, 2010, Figure 1]. This assumption
could be inappropriate for estimating N trace gas emissions,since simulated soil water content would rarely be in ranges
close to saturation, a condition favoring anaerobic N gaseousemissions during and after rainfall or snowmelt [del Prado
et al., 2006]. This approach may also produce a very different drainage pattern than those based on the solution of the
Richards flow equation, especially at smaller time scales[Scanlon et al., 2002] or different ecosystems than those in thework by Pachepsky et al. [2006].
[5] There exist many modeling studies that coupled transientwater dynamics and complicated N biogeochemical processes
in order to investigate water effects on N biogeochemical cycles

G03040especially for gaseous N emission relevant processes due to its
high temporal variability [e.g., Gu and Riley, 2010; Gu et al.,2009; Hendriks et al., 2007; Stolk et al., 2011]. It appears
that model improvements of N emission process responding towater conditions and soil conditions are still in need. On theother hand, Doltra and Muñoz [2010] investigated effects oftransient water dynamic model (in HYDRUS‐2D) and semi‐
dynamic water flow module in EU‐Rotate_N model in N‐
fertilized vegetable fields, showing that N uptake difference
resulted from evapotranspiration simulation could cause
observable difference on soil nitrate and leaching losses.
However, it is still not clear to what degree the importantbiogeochemical processes (e.g., other soil N processes besides
N uptake) would be affected at different temporal scales if
biogeochemical models are coupled with either a Richards
equation based water model or a simple tipping‐bucket orsemi‐dynamic water model alone. Analysis on such an issue
would be useful to decide what sort of soil water or hydrological algorithms should be included in biogeochemical or
ecological models.
[6] In this study soil excess N losses (both leaching andemission) simulated by a well‐established biogeochemical
model (DAYCENT [Parton et al., 2001]) in mixed‐coniferous
forest ecosystem under substantial N deposition will be
investigated, by implementing a numerical solution of the one‐dimensional (1D) transient soil water flow equation (based on
the Richards equation, from HYDRUS‐1D [Šimůnek et al.,2008]) into it as an alternative to its original water flow module [Parton et al., 1998], which can be categorized as a semi‐
dynamic tipping bucket model [de Jong and Bootsma, 1996].
Comparisons of results produced by these two model versions,
hereafter referred to as DAYCENT models with the “original”
and “HYDRUS” water flow module, were conducted withregards to effects on simulated soil N accumulation andleaching and emission losses. The objective is to illustrate theimportance of adequately simulating transient soil water
dynamics and flows to estimate gaseous N emissions andN leaching and other N processes.

2. Materials and Methods
2.1. Site Description[7] The San Bernardino National Forest (SBNF) is locatedeast of Los Angeles, California. It is a mixed conifer forest,with dominant tree species ponderosa (Pinus ponderosa
Laws) and Jeffrey pine (P. jeffreyi Grev. and Balf.). Logging
has not occurred for many years and forest fires have been
suppressed since the early 1900s. Nitrogen deposition monitoring has demonstrated that in the San Bernardino Mountains there is an N deposition gradient from west to east andsouth with decreased deposition, with increasing distancefrom the urban region. Precipitation declines along this gradient as well [Fenn et al., 2003]. In this study, a high N
deposition site, Camp Paivika (CP), on the western end of thisgradient, and a relatively low deposition site, Barton Flats(BF), 46 km east of CP, were chosen. The current N
deposition levels at CP and BF were estimated at 70 and
8.8 kg N ha−1 yr−1 [Fenn et al., 2008], of which about
5.4 (CP) and 2.4 (BF) kg N ha−1 yr−1 occur in the form of
wet deposition (M. E. Fenn, unpublished data).
[8] The soils are generally sandy loam in texture, with 70–
72% sand, 8–9% silt, and 19–21% clay for CP, and 66–69%

2 of 15

G03040

G03040YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES
Table 1. Soil Physical Properties, Soil Organic Matter, and Root Fraction for Soil Layers Required by DAYCENT,
Based on the Measurement of Grulke et al. [1998] at Experimental Sites Camp Paivika, Barton Flats, and James
Mountain Reserves [Vargas and Allen, 2008] in or Near the San Bernardino National Forest, California
Silt

Sand

Clay

Soil Organic
Matter (%)

Root
Fraction

0–1, 1–4
4–15
15–30, 30–45
45–60, 60–75, 75–90, 90–105 105–120

Camp Paivika
1.32
0.09
1.43
0.09
1.53
0.09
1.55
0.09

0.70
0.71
0.72
0.72

0.21
0.20
0.19
0.19

0.21
0.20
0.19
0.10

0.6
0.2
0.1
0.1

0–1, 1–4
4–15
15–30
30–45, 45–60, 60–75, 75–90, 90–105, 105–120

Barton Flats
1.47
0.12
1.53
0.09
1.58
0.08
1.59
0.08

0.66
0.67
0.69
0.69

0.22
0.24
0.23
0.23

0.22
0.25
0.24
0.10

0.6
0.2
0.1
0.1

James Mountains Reserve
1.20
0.10
1.30
0.15

0.83
0.82

0.07
0.03

0.08
0.04

0.9
0.1

Layers (cm)

0–2, 2–8, 8–16
16–26, 26–36, 36–46

Bulk Density(g cm−3)

sand, 8–12% silt, and 22–25% clay for BF [Grulke et al.,1998]. For applying DAYCENT, it was suggested that thesurface 15 cm soil layers are divided into 3 thin layers (e.g.,
0–1, 1–4, and 4–15 cm in this study) in order to simulate therapid changes of near surface water and temperature conditions, which are critical to gaseous N emissions [Parton
et al., 1998]. The soil physical properties, SOM, and root
fractions of the soil profiles needed by the models are shownin Tables 1 and 2.
[9] In order to assess how well the newly incorporated soil
water module can simulate soil moisture dynamics, a model
simulation was conducted and compared with soil waterobservations at the University of California James San Jacinto
Mountains Reserve (http://www.jamesreserve.edu) located in
the San Jacinto Mountains in Riverside County, California[Kitajima et al., 2010]. The site is an old‐growth stand of
California semiarid mixed conifer forest, surrounded by the
San Bernardino Mountains. Soils are shallow Entisols(Mixed, mesic, shallow Typic Xeropsamment) with a loamy‐
sand texture and underlay by weathered granitic bedrock. The

soil water content was measured using soil moisture ECHO
probes (ONESET®, S‐SMA‐M005) every 5–15 min at 2, 8,
and 16 cm depth since the fall 2004 [Allen et al., 2007;
Kitajima et al., 2010]. Weather data such as air temperatureand precipitation are also monitored in situ or nearby station.
In this study soil water content during period of 2006–2009
was compared, with measured or derived soil properties at thesite (Tables 1 and 2).
[10] The two DAYCENT model versions are comparedwith regards to soil solution nitrate concentration at CP and
BF sites. In 2000 at each site 6 soil pits were dug deep enough,and 3 lysimeter soil solution collectors were installed at a
depth of 1 m. The lysimeter collectors were horizontally
put into tunnels of the pit walls, following standard procedure
to be sure the lysimeter being contacted with soils well.
From 2001–2004, the soil solutions were extracted 13 times
at CP and 11 times at BF with at least one effective sample(i.e., enough water volume for record and analysis), for whichthen the NH+4 ‐N and NO−3 ‐N concentrations were analyzed.
The NH+4 ‐N in these samples was minor and can be neglected,

Table 2. Soil Bulk Density and Texture (see Table 1) Derived Hydraulic Properties and Parameters for the 10 Layers
Used by DAYCENT (r, wp, fc, s, Ks) and HYDRUS (r, s, Ks, a, n, l) at Experimental Sites Camp Paivika,
Barton Flats, and James Mountain Reserves in or Near the San Bernardino National Forest, California
wpa

fca

sa

Ksb

a, n, lc

0–1, 1–4
4–15
15–30, 30–45
45–60, 60–75, 75–90, 90–105 105–120

Camp Paivika
0.061
0.170
0.060
0.151
0.054
0.130
0.052
0.125

0.280
0.240
0.212
0.202

0.450
0.416
0.389
0.371

57.71
46.87
40.48
34.24

0.0296, 1.292,0.5
0.0381, 1.320, 0.5
0.0417, 1.348, 0.5
0.0423, 1.357, 0.5

0–1, 1–4
4–15
15–30
30–45, 45–60, 60–75, 75–90, 90–105, 105–120

Barton Flats
0.058
0.155
0.058
0.164
0.055
0.151
0.055
0.144

0.250
0.254
0.233
0.223

0.404
0.386
0.369
0.366

31.80
24.80
24.15
23.90

0.0326, 1.302, 0.5
0.0332, 1.290, 0.5
0.0370, 1.308, 0.5
0.0386, 1.319, 0.5

James Mountains Reserve
0.019
0.099
0.153
0.020
0.060
0.123

0.475
0.461

274.4
303.9

0.0399, 1.678, 0.5
0.0456, 1.661, 0.5

Layers (cm)

0–2, 2–8, 8–16
16–26, 26–36, 36–46

rar, wp, fc and s are soil residual water content, wilting point, field capacity, and saturated water content, respectively.
Ks is soil saturated hydraulic conductivity (cm day−1).
ca (cm−1), n, l, are parameters of the van Genuchten type analytical functions.
ab

3 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

so only the NO−3 ‐N concentrations are compared to our model
simulations of deep soil nitrate concentration averaged oversoils from 60–120 cm depths. So the model may be evaluatedwith regard to its ability to capture deep soil nitrate accumulation with time, rather than N transport, which needs estimation of water flow.
2.2. Overview of DAYCENT’s
N Biogeochemical Module[11] The DAYCENT biogeochemical model [Del Grosso
et al., 2000, 2001; Parton et al., 1998, 2001] had been
applied to simulate plant growth, litter‐fall, soil organic
matter dynamics and the tightly coupled soil‐plant N cycle.
The model had also been successfully applied to estimate
N trace gas emissions for agricultural systems across theUnited States [Del Grosso et al., 2006] and globally [DelGrosso et al., 2009]. An important feature of DAYCENT
is its N cycling algorithm, which considers N2O, NOx, andN2 emissions resulting from nitrification and denitrification[Parton et al., 1988a, 2001], and detailed soil N mineralization‐immobilization processes associated with variousorganic matter pools (including microbial), which are distinguished by their decomposition rates [Parton et al.,1988b]. In this study, the model was calibrated mainly bymodifying parameters relevant to tree growth, N adsorption
by litter, and litter and SOM decomposition rates in order toachieve reasonably well‐fitted outputs with observations oftree biomass and allocation, litter and SOM accumulation,
and C/N ratios as previously reported [Fenn et al., 2008].
Otherwise, for other N processes it was not calibrated andparameters used those already included in the model asreported by Parton et al. [1998] and others [e.g., Del Grosso
et al., 2002, 2009].
[12] DAYCENT simulates N mineralization from released
N bonded with 3 SOM components (active, slow, and passive pools) and metabolic residue of litter materials when
respired, assuming that C/N ratios remain fixed, while
N immobilization from soil inorganic N pools when litter
structural residue decomposes [Parton et al., 1988b]. The net
mineralized N from SOM decomposition and/or respiration
adds into soil inorganic NH+4 ‐N pools, which is further
nitrified into NO−3 ‐N. Together with N inputs by bio‐fixation, deposition, and/or fertilization (N sources), soil NH+4 ‐N
and NO−3 ‐N can be taken up by roots. The net N addition (i.e.,mineralization + sources − uptake) stores in soil excess N
pools, which are subjected to downward transport with waterflow (even leaching out of root zones with drainage) andgaseous emissions. When water flows between layers, thenitrate portion of flowing water over the upper layer’s total
water volume transports into the next layers. Although anempirical factor in the DAYCENT model accounting fornitrate transport not associated with water flow may be
applied by setting a constant fraction of nitrate moving
downward to the lower horizons, we arbitrarily turned thisparameter off in order to consider the direct effects ofdrainage (and internal water flow) by the two contrastingwater flow modules (see below) on nitrogen leaching.
[13] DAYCENT’s N gas submodel is based on the leaky
pipe concept, in which total N gas emissions are proportional to soil N transformation rates (nitrification and denitrification) and soil gas diffusivity determines the ratios ofrespective N gas species (mainly N2O, N2, and NOx)

G03040[Parton et al., 2001]. Soil NO−3 ‐N and N2 will be final
products of nitrification and denitrification, respectively, inboth of which N2O and NOx will be produced (leaking) as
by‐products. The model first calculates N2O emission rates
from those two processes, and then the ratio of NOx:N2O is
calculated according to an empirical relationship with gas
diffusivity. The ratio of NOx:N2O in nitrification is also
controlled by a pulse multiplier counting for the effect ofrain into previous dry soil. In all those N processes(including gas diffusivity), soil WFPS (i.e., soil moisture),
and/or soil texture are the major environmental control
factors, while soil temperature only affects soil N mineralization and nitrification in the model. The details andparameters refer to the work by Parton et al. [2001].
2.3. Incorporation of HYDRUS‐1D Water Module
into DAYCENT and Parameterization[14] The DAYCENT original soil water module [Parton
et al., 1998] has a two‐stage process to describe infiltration and redistribution of water. During periods of rainfall orsnowmelt, water infiltrates into a surface soil horizon until it
is saturated and then flows into the next layer at saturated
conductivity (Ks in Table 2). After infiltration, soil water inany horizon above its field capacity will be drained into thenext layer and finally drained out of the soil bottom at saturated conductivity as well. Any extra rainfall or snowmelt
water is removed as run‐off. Following this infiltration/
percolation stage, bidirectional unsaturated flux proceeds
according to the Darcy equation, which also multiplied by aso‐called “damping multiplier” (dmpflux) (default 0.000001)[Parton et al., 1998], as following,
dmpflux * hpoti1  hpoti * avcondi
fluxi ¼
disti

ð1Þ

where, hpoti−1 and hpoti are the hydraulic potential (i.e., sumof soil matric potential and gravitational head) of two
adjacent layers, i − 1 and i, avcondi is the weighted average ofunsaturated conductivity of layer (i) and its upper (i − 1) andlower ones (i + 1), and disti is the distance of the middle
points of two layers. The success of this approach depends
upon the carefully designed vertical soil layer discretizationand calculation time steps for the two processes [de Jongand Bootsma, 1996; Ranatunga et al., 2008]. Practically,in the DAYCENT water module, drainage rates are adjustable using bottom layer’s base flow or stormflow factors(ranging from 0 to 1), and internal flow rates by changing“damping multiplier.”[15] In this study, in order to better understand the impact
of soil water contents and overall ecosystem water balanceon these N processes and its losses, the soil water module of
HYDRUS‐1D [Šimůnek et al., 2008] was successfully
incorporated as an alternative submodule to DAYCENT’s
two modules for infiltration and redistribution (unsaturatedflow). The HYDRUS water module uses the method of
linear mass‐lumped finite elements to numerically solve thehighly nonlinear Richards equation describing both saturated and unsaturated water flow,

4 of 15

 

@
@
@h
¼
K
þ1 S
@t @x
@x

ð2Þ

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

G03040Figure 1. Comparison of surface (0–4 or 2 cm) soil hydraulic property functions derived from the soil
texture: (a) soil water retention curves, and (b) unsaturated conductivity, used by DAYCENT original and
HYDRUS water flow module at Camp Paivika (CP), Barton Flats (BF), and James Mountain Reserves(JMR), in or near the San Bernardino National Forest, California.
where h is the water pressure head,  is the volumetric watercontent, t is the time, x is the distance, S is the water sink/
source term, and K is the unsaturated conductivity. Both
 and K are dependent on pressure head h, which can be
expressed in the van Genuchten functions as,ðhÞr þs  rm½1 þ jhjn 

ð3aÞ

and
h

k ðhÞ ¼ Ks Sel 1  1  Se1=m

m i2

ð3bÞ

where r and s are residue and saturated water content, Se iseffective saturation (i.e., = ( − r)/(s − r)), m = 1 − 1/n, anda, n, l are the empirical parameters. All of those independent
parameters can be derived from pedotransfer functions [e.g.,
Schaap et al., 2001; Šimůnek et al., 2008; Wösten et al., 2001].
In this water module, runoff and drainage can be estimatedcorrespondingly, without any empirical adjustment as theoriginal DAYCENT water module.
[16] Both water modules were driven by the same water
input and potential evapotranspiration, which were calculated
in the original DAYCENT water submodels; then, actual
evapotranspiration and soil water content were used to drive
plant growth and soil C and N biogeochemical processes. Due
to differences in soil water content predicted by the twomodels, actual evapotranspiration is expected to be different
as well, although the potential evapotranspiration and root
distribution are identical.
[17] The soil hydraulic property functions, characterizing
soil’s retention and conductive properties, are different in two
water modules. While the DAYCENT water module uses the
Saxon‐Rawls method, HYDRUS uses the van Genuchten
type analytical functions (see equation (3) above) or others(e.g., Brooks and Corey, or the modified van Genuchten type)[see Šimůnek et al., 2008, and references therein]. The bulk
density (Table 1), the wilting point, and field capacity were
derived based on soil texture and soil organic matter (SOM)
using the Saxon‐Rawls method [Saxton et al., 1986], except

for saturated conductivity (Table 2). The Rosetta module,
included in the HYDRUS‐1D package, was applied todetermine the van Genuchten parameters required by the
HYDRUS water module [Schaap et al., 2001; Šimůnek et al.,2008], using soil texture, bulk density, and soil water content
at −33 kPa (i.e., field capacity as determined above) (Table 2).
In this procedure of estimating soil hydraulic properties, thecommon parameters for the two modules were forced to be
the same in order to allow the hydraulic functions to be as
similar as possible so that soil water flow calculation affectedby those functions could be minimized, as reported by Doltra
and Muñoz [2010].
[18] However, the hydraulic functions are still different.
For example, while soil water retention curves and hydraulicconductivity functions for the depth of 0–4 cm at CP and BF
site intersect at around field capacity and the pressure headof −10 cm (Figure 1), they show relatively large differences
near saturation. While this may cause differences in predicted soil water dynamics and movement, it is not possible
to predict differences in modeling output by comparing
these functions directly. While HYDRUS‐1D applies thesefunctions in the entire range of soil water contents in its
water modules, DAYCENT uses them only between theresidual soil water content and field capacity. Additionally,
DAYCENT’s water flow rates in the unsaturated flow stage
are calculated using Darcy’s equation and then multiplied bya so‐called “damping multiplier” (default 0.000001) as
described above.
2.4. Model Simulations[19] Inputs required by DAYCENT include daily maximum/minimum temperatures and precipitation, current andhistorical land use, disturbances and/or management, andsite‐specific soil properties and hydraulic characteristics. For
the driving weather data, we used a 65‐year daily data set
from the weather stations in the Lake Arrowhead area (1942–
2006) downloaded from the NOAA National Climate Data
Center (http://www.ncdc.noaa.gov/oa/climate/stationlocator.html).
Because our two SBNF sites are away from the weather stations in the Lake Arrowhead area, precipitation for CP and BF
needs to be extrapolated. We extracted monthly precipitation

5 of 15

G03040

G03040YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

Table 3. Annual Precipitation and Simulated Drainage, Runoff, Actual Evapotranspiration, and Soil Water Change by HYDRUS or
DAYCENT Original Soil Water Module for Years 2001–2005 at Experimental Sites Camp Paivika and Barton Flats in the San Bernardino National Forest, Californiaa
Drainage

Actual
Evapotranspiration

Soil Water Change

D

H

D

H

D

RunoffYear

Precipitation

H

D

H

2001
2002
2003
2004
2005

67.2
29.7
76.6
29.8
141.0

42.8
12.7
52.5
9.3
116.3

44.4
11.9
49.4
11.0
123.0

Camp Paivika
0
0
0
0
0

0
0
4.2
0
2.3

18.1
19.0
22.1
18.6
8.6

15.7
19.3
20.6
18.6
8.1

6.2
−2.1
1.9
1.8
15.9

7.0
−1.5
2.3
0.1
7.5

2001
2002
2003
2004
2005

53.8
15.6
69.3
32.1
151.5

27.4
0.34
41.7
10.7
111.9

24.8
0.003
35.0
10.4
99.3

Barton Flats
0
0
0
0
16.4

1.3
1.5
6.6
0
30.1

26.3
20.2
20.7
24.9
15.3

27.4
20.1
17.1
25.0
15.7

0.1
−5.0
6.9
−3.5
7.5

0.3
−6.1
10.5
−3.3
6.3

a

HYDRUS, H; DAYCENT Original, D. Annual precipitation measured from previous June to July. All measurements are in centimeters.

for the two sites from the PRISM database [Di Luzio et al.,2008] and then downscaled to daily values assuming thesame time variation pattern as at the Lake Arrowhead station.
Annual precipitation averaged 88.2 cm at CP and 59.1 cm
at BF, with more than 70% of precipitation falling from
December to March. Since both sites exhibited a Mediterranean climate, the yearly cycle considered in this paper
started in July and ended the following June. For the period
from 2001–2005, the annual precipitation (from previous
July to June) was 67.2, 29.7, 76.6, and 29.8, and 141.0 cm at
CP, and 53.8, 15.6, 69.3, 32.1, and 151.5 cm at BF (seeTable 3 and Figures 4a and 4d for the daily patterns). Years
2002 and 2004 were dry, 2001 and 2003 were nearly average, and 2005 was wet. For the James San Jacinto Mountains

Reserve site (JMR), the air temperature data from 2006–2009
was directly obtained from the in situ observations andprecipitation data was from the nearby Keenwild station(Figure 2a). The weather driving data beyond this period was
obtained from the PRISM database and adjusted accordingly
using the same procedure as for CP and BF.
[20] Model runs were initialized and stabilized by repeatedly using site‐specific weather data of all available years for
900 years (i.e., initially started from 1000 AD). A backgroundN deposition rate of 0.1 kg N ha−1 yr−1 and fire occurrence ofevery 100 years were used. For the final model runs (2001–
2005), the forests were assumed to start growing after a fire in
1895 without any additional fire. The model was calibrated toobtain historical tree growth, litter C, SOM and C/N ratios in

Figure 2. Precipitation (cm day−1) and soil volume water contents (VWC) at 2, 8, and 16 cm depths
observed and simulated by DAYCENT original and HYDRUS water flow module at the James San
Jacinto Mountains Reserves (JMR) near the San Bernardino National Forest, California, from 2006 to2009.
6 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

G03040Figure 3. Soil solution NO−3 ‐N concentrations (mg L−1) measured in the lysimeter (mean ± standard
deviation; no error bar if only one observation) and simulated using the DAYCENT model with the original or HYDRUS water flow modules at (a) Camp Paivika (CP) and (b) Barton Flats (BF), San Bernardino
National Forest, California, from 2001 to 2004.
the study area, and some parameters, mainly tree growth andmortality, were modified in order to achieve reasonably well‐fitted outputs with these observations [Fenn et al., 2008]. The
model inputs assumed that gradually increasing N deposition
began in 1930 and in 1980 reached the current levels of
70.0 and 8.8 kg N ha−1 yr−1 for CP and BF, respectively, andremained constant thereafter [Fenn et al., 2008].

3. Results
3.1. Model Evaluation in Soil Water Dynamics
and Deep Nitrate Simulation[21] As described above, the original soil water module in
DAYCENT treats water movement as a two‐stage process.
This algorithm produces surface soil saturation during periodsof continuous and large rainfall (or snowmelt) and thereafter a
rapid drop to field capacity (Figure 2). At the James San Jacinto
Mountains Reserve (JMR), surface saturation occurred very
rarely from 2006–2009, probably due to its sandy texture withvery high saturated conductivity (Table 2). During most of therainy season, the soil surface remained at field capacity or
below. Apparently the HYDRUS water flow module produced
soil moisture patterns during the rainy season that were more
realistic and corresponded better to the precipitation pattern(Figure 2). Compared to the measurements, the HYDRUS
water module appeared to over‐estimate soil moisture during
large rain events. And both water modules showed time‐
lagged high soil water content during the late precipitationseasons, probably because DAYCENT assumed some snow
accumulation in winter if temperature dropped to zero, which
was melting later on. The observation demonstrated this
could be incorrect. Since the two water flow modules produced
different water content patterns (especially during rainy periods), it can be expected that this will result in different estimates of soil nitrogen gaseous emissions, because generally
near‐saturation conditions promote gaseous emissions in the
DAYCENT model.

[22] At CP and BF sites, simulated soil solution nitrateconcentrations, calculated as the soil NO−3 ‐N content divided
by the amount of soil water in between 60 and 120 cm depths,
are compared in Figure 3 with measured NO−3 ‐N concentrations, averaged over the sampled lysimeters at a depth of 1 m(samples varied from 1–6, because very often soil solution
samples could not be extracted from all installed lysimeters).
Both measurements and simulations showed that late fall
rainfall events transported nitrate downward so that deep soil
nitrate concentrations reached an annually high value at the
beginning of the winter rainy season. Further precipitation
then leached soil nitrate out of the profile, resulting in very
low nitrate concentrations below 60 cm. It appears that, ingeneral, both modules adequately represented deep soil
nitrate dynamics (time series), although the absolute nitrateconcentration comparison might be inappropriate because(1) there was large variance in the lysimeter collected measurements, and (2) the soil solution collected by lysimeters
installed at 1 m depth may not identically represent the simulated soil portions. However, it still provided useful information for model evaluation. For an example, the measurements
were almost carried out in the 2003 rainy season, during
which the HYDRUS water flow module simulated a higher
deep soil NO−3 ‐N accumulation prior to drainage but lower by
the end of early spring drainage period compared to theoriginal module, especially at the high‐N deposition CP site
(Figure 3a). This result might be related to more N downward
movement and final leaching simulated by the HYDRUS
module (see below). Apparently, neither model simulated the
deep soil NO−3 ‐N peaks in the dry 2004 year at both sites(Figure 3), nor high concentrations observed in the wet 2003
year (see Table 3) at the CP site, indicating that further
improvement in DAYCENT’s solute transport module may
be needed. Both simulations and lysimeter measurements
showed that although the annual N deposition rate at the CPsite (70 kg N ha−1) was nearly 8 times larger than at the BF site
(8.8 kg N ha−1), deep soil nitrate accumulation was not very

7 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

G03040Figure 4. (a–f) Daily precipitation (cm day−1) and simulated differences of actual evapotranspiration(AET) and soil water change by the DAYCENT original and HYDRUS water flow modules at Camp
Paivika (CP) and Barton Flats (BF), San Bernardino National Forest, California, from 2001 to 2005.

different from both modeling and observations. This result
implies that nitrate transported downward would eventually
be leached out of the profile, irrespective of the rates ofexternal N sources and/or mineralization at these two sites.
3.2. Comparison of Soil Water Drainage Simulation
by Two Water Modules[23] At our two study sites (CP and BF), the HYDRUS
water module did not produce runoff except for very wet 2005
year at BF, and the original water module only simulated avery small amount of runoff (Table 3). Different soil watercontent simulated by the two modules did not significantly
affect estimates of actual evapotranspiration (AET) either;less than 0.15 cm day−1 and 2.4 cm per year at most (Figure 4b
and Table 3). But it did affect the timing and pattern of soil
water change (Figure 4c) and thus drainage (Figure 5) remarkably, although the differences at yearly were small (Table 3).
After a long dry period in summer and early fall, initial precipitation events did not produce much water flow out of thesystem until early spring, when soil water storage had risen
and snowmelt or rain continued (Figure 4a versus Figure 5afor CP site, Figure 4d versus Figure 5b for BF site). Drainage predicted by the DAYCENT original water module
stopped earlier and was stronger than by the HYDRUS water
module, in which drainage lasted until May, except for the
wet year 2005 (Figures 5a and 5b).
[24] HYDRUS water module simulated small drainage
fluxes on more days. It, however, simulated fewer large
drainage events, compared to the original water module(Figures 5c and 5d). For the time period (total 1827 days) of
2001–2005 at the CP site, daily drainage of less than 0.25,
from 0.25 to 1.0, and more than 1.0 cm per day, simulated by
the HYDRUS module, occurred on 278, 238, and 19 days,

respectively. The original water module simulated drainage
for the same intervals on 67, 198, and 26 days, respectively(Figure 5c). At the BF site, the same ranges of drainage were
simulated by the HYDRUS module on 308, 142, and 20 days,
and by the original water module on 37, 108, and 27 days,
respectively (Figure 5d). The original water flow module
simulated drainage fluxes of higher intensity and smallerfrequency, i.e., more pulsed, than the HYDRUS water module. This result was especially true during drier 2002 and 2004
years, but less so during the wetter 2005 (Figures 5a and 5b).
[25] By integrating the daily drainage over monthly andyearly periods for 2001–2005 at the CP site, the coefficient ofvariation decreased from 4.752 to 2.127 (monthly) and 0.924(yearly) for the HYDRUS water module simulations, andfrom 6.633 to 2.345 (monthly) and 0.986 (yearly) for theoriginal water module simulations. At the BF site, the coefficient of variation for daily, monthly, and yearly time‐scale
was 5.477, 2.277, and 1.102 for the HYDRUS module, and7.587, 2.630, and 1.123 for the original module, respectively.
These results show that the simulated difference in drainage
volumes between the two modules decreased with increased
time averaging period, until there was no significant difference at the annual time scale (Table 3).
3.3. Difference of DAYCENT Simulated N Leaching
Losses with Two Water Modules[26] In the simulations, most N leaching in 2001 through2005 at CP and BF sites occurred in the late rainfall season(from late winter to early spring), corresponding to drainage
events (Figure 6b versus Figure 5a for CP site, and Figure 6e
versus Figure 5b for BF site) rather than net addition to soil
excess N pool (i.e., mineralization + source − uptake) (Figure 6aand Figure 6d). The rates of N loss generally followed the

8 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

G03040Figure 5. (a–d) Comparison of monthly soil water drainage rates (cm month−1) and number of days
when daily drainage rates were in specified ranges simulated by the DAYCENT original and HYDRUS
water flow modules at Camp Paivika (CP) and Barton Flats (BF), San Bernardino National Forest,
California, from 2001 to 2005.
time trends in drainage rates. Consequently, N leaching losses
in the dry years of 2002 and 2004 were low compared to highN deposition amounts (see Figure 6 and Tables 3 and 4). In the
wet year of 2005, very high N leaching in December wasconsistent with the intensive drainage during the same period,
but it was not proportional to drainage. While drainage was

2–3 times greater during the wet year compared to the normal
year, it did not produce correspondingly larger nitrate leaching losses.
[27] The 2002 and 2004 years were dry, so the leachinglosses predicted by the two water modules were very low
and only slightly different. However, during the average

Figure 6. (a–f) Monthly N losses (leaching at 120 cm soil bottom and gaseous emissions) (kg ha−1month−1) simulated by the DAYCENT model with either the original or HYDRUS water flow module
at Camp Paivika (CP) and Barton Flats (BF), San Bernardino National Forest, California, from 2001
to 2005.
9 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

Table 4. Annual N Net Addition to Soil Excess N Pools and Their
Losses Due to Leaching and Gaseous Emission, Simulated by
HYDRUS or DAYCENT Original Soil Water Module, at theHigh‐N Deposition Site Camp Paivika and the Low‐N Deposition
Site Barton Flats in the San Bernardino National Forest, California,
in 2001–2005a
Net Addition to
Soil Excess
N PoolEmission(kg N ha−1)

Leaching

Year

H

D

H

D

H

D

2001
2002
2003
2004
2005

123.3
116.4
135.0
107.6
115.7

88.3
80.3
99.2
69.9
89.7

Camp Paivika
48.6
6.24
38.7
3.35
47.5

33.3
5.26
20.2
3.31
33.0

2.42
2.94
6.85
2.51
5.39

0.84
1.42
4.23
1.30
3.42

2001
2002
2003
2004
2005

36.9
33.8
42.8
33.6
37.8

25.8
20.4
29.3
23.3
26.6

Barton Flats
9.35
0.00
21.8
2.29
25.4

8.65
0.00
10.7
1.46
20.4

1.72
1.68
3.77
2.17
1.82

0.87
0.92
2.59
1.08
0.94

a
Net addition is Mineralization plus Sources minus Uptake. Leaching is
at 120 cm soil bottom. Emissions (kg N ha−1) are measured from previous
July to June. HYDRUS, H; DAYCENT Original, D. Camp Paivika is 70 kg
N ha−1 yr−1, and Barton Flats is 8.8 kg N ha−1 yr−1.

years of 2001 and 2003 and the wet year of 2005, the rates ofsoil nitrogen leaching losses predicted by the HYDRUS flow
module were very different compared to those calculated by
the original water module (Figure 6b and 6e). In these 3 years,
model runs that used the HYDRUS water flow module simulated higher N leaching losses than those that used theoriginal water flow module. Cumulative annual N leachinglosses (from previous July to June) (Table 4) at the high‐N CP
site simulated using the HYDRUS flow module were about
39 ∼ 49 kg N ha−1, while those simulated using the original
DAYCENT water flow module were about 14 ∼ 18 kg N ha−1lower. At the low‐N deposition BF site, leaching N losses
were highly variable, ranging from 9 kg N ha−1 in 2001 to
more than 20 kg N ha−1 in both 2003 and 2005. In simulations
with the HYDRUS water flow module the leaching was about
0.7, 11, and 5 kg N ha−1 higher in 2001, 2003, and 2005,
respectively, than in simulations with the DAYCENT waterflow module (Table 4).
[28] The difference of net N additions to soil excess pool by
two water modules could affect the absolute N leaching losses
at CP and BF site (Table 4 and Figure 6), but did not primarily
contribute to the difference of leaching estimation by twowater modules. By contrast, Figure 7 showed that the simulated N leaching differences by the HYDRUS and original
water modules was mostly correlated to the water drainage
differences of two modules at daily time step (Figure 7a)rather than at monthly (Figure 7b) or yearly (Figure 7c). At
the wetter CP site, regression R2 of N leaching difference anddrainage difference between two modules decreased from
0.292 at daily to 0.01 at monthly and 0.00005 at yearly, andat the relatively drier BF site, it changed from 0.382 at daily
to 0.169 at monthly and 0.259 at yearly time step. Figure 7
also showed that soil drainage difference between two
water flow modules mostly varied within 10 cm H2O
regardless of time step, while the resulted N leaching loss

G03040estimation differences varied greatly at various time steps. It
clearly demonstrated that short‐term water flow dynamics
would remarkably affect the N downward transport and its
leaching losses.
3.4. Difference of DAYCENT Simulated N GaseousEmissions with Two Water Modules[29] Although N gaseous losses were relatively small,they occurred throughout almost the entire precipitationseason, during which two major peaks appeared; one at thestart and another at the end of the rainy season (Figure 6cand 6f). Usually the late spring N emissions were muchhigher than those in the previous fall, except for the very wet
year 2005. The low N emission between these two peaks
occurred during the period of very strong N leaching losses.
The timing and amounts of leaching and emission compensated for each other (Figure 6b versus Figure 6c for CP
site, Figure 6e versus Figure 6f for BF site). The 2005 rainyseason started earlier than in other years, and stimulated
larger production of trace N gases. The strong precipitation
events that started in the middle of the rainy season thenalmost exhausted the soil nitrate near the soil surface,leaving not enough N for an observable gaseous emission
peak in the late spring of 2005, as shown in other years.
[30] With the two‐stage water flow algorithm used in the
DAYCENT original water module, very few precipitation

Figure 7. Comparison of the DAYCENT simulated
N leaching difference and the soil drainage difference,
between the original and HYDRUS water flow module,at (a) daily, (b) monthly, and (c) yearly time step, at Camp
Paivika (CP) and Barton Flats (BF), San Bernardino
National Forest, California, from 2001 to 2005.

10 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

G03040Figure 8. Comparison of the monthly DAYCENT simulated (a) N2O ratios of via denitrification to via
nitrification, and (b) ratios of NOx to N2O, between the original and HYDRUS water flow module, at Camp
Paivika (CP) and Barton Flats (BF), San Bernardino National Forest, California, from 2001 to 2005. Notethat high NOx:N2O ratios generally imply low water‐filled pore space (WFPS) in Figure 8b [see Parton
et al., 2001].
events caused soil saturation, and even when such conditions
occurred, soil moisture would very quickly drop back to field
capacity, which would then persist or be further reduced (seeFigure 2). As expected, calculations that used the HYDRUS
water module simulated higher N gaseous emissions overthe entire study period than those that used the original
DAYCENT water module (Figures 6c and 6f and Table 4).
In this study, DAYCENT with the HYDRUS water module
predicted the annual emission of 2.4 ∼ 6.8 kg N ha−1 at CP
and 1.7 ∼ 3.8 kg N ha−1 at BF, while with the original moduleof 0.8 ∼ 4.2 kg N ha−1 at CP and 0.9 ∼ 2.6 kg N ha−1 at BF(Table 4). In other words, the difference was about 1.2 ∼ 2.6
and 0.7 ∼ 1.2 kg N ha−1 year−1 for the high‐N CP site and the
low‐N BF site, respectively, or about 1.5 ∼ 2 times higher
with the HYDRUS water module. Thus differences in trace N
emissions estimated by the two water modules cannot beignored in this heavily N deposition ecosystem.
[31] In this study, of total N emissions, NOx gas was the
major species with 80–85% at CP and 75% at BF, N2O wasabout 15–18% at CP and 22% at BF, and the rest 1–3% was inform of N2. As reasonably expected from soil water simulations, especially during late rainy season (Figures 4c and 4f),
N gaseous emission difference between two water modules
described above was mainly related to soil water controlled
N biochemical processes. As briefly described previously,
DAYCENT calculated N2O emissions as leaking from both
aerobic nitrification and anaerobic denitrification processes.
Figure 8a compared monthly N2O ratios of via denitrification:
via nitrification, indicating that the model with its original
water module only simulated large denitrification produced
N2O occasionally, while with incorporated HYDRUS water
module it produced much more N2O from denitrification. On
the other hand, ratios of NOx:N2O by DAYCENT with theoriginal water module were about 56% larger than those with
HYDRUS water module (Figure 8b). Generally the ratios of
NOx:N2O decrease with increasing WFPS, because gas diffusivity decreases with WFPS [Parton et al., 2001]. Our
results indicated that the original DAYCENT model produced higher NOx:N2O ratios because of lower soil WFPS

simulated, compared to the HYDRUS water module incorporated model.

4. Discussion
4.1. Soil Water Content Dynamic and Drainages[32] The original and new HYDRUS soil water flow
modules produced different soil water and drainage temporal
patterns because of different assumptions involved in eachof these modules, and thus different governing flow equations
being solved, and different parameters (see Table 2 andFigure 1) used to characterize soil hydraulic properties. For
example, in the tipping bucket model, field capacity is a key
parameter, but it is not used in models like HYDRUS solving
Richards equation for water flow. In many cases, thisparameter is even not rigorously defined [Twarakavi et al.,2009]. These differences resulted in soil moisture patterns
produced by these two modules that were quite different,
especially after rainfall events (Figure 2). The field observations showed that water contents higher than field capacity
would last for a longer time period as simulated by the
HYDRUS module, and that quick infiltration of waterfrom top to bottom as predicted by the original flow module
didn’t capture observed water content variations with time
(Figure 2). This difference in transient soil water status nearsaturation would likely have large effects on soil N statesand dynamics of soil biogeochemical processes, such as
mineralization and nitrification occurring mainly under aerobic conditions and denitrification under anaerobic conditions, both of which would affect N downward transport andgaseous emissions. In DAYCENT this issue could be overcome by: (i) dividing the topsoil layer into very thin layers(e.g., 0–1 and 1–4 cm, or more detailed horizon discretion),and (ii) adjusting the precipitation event period (e.g., from4 to 2 h), so that the surface soils could be saturated even withsmall rainfall events [Parton et al., 1998]. Obviously, thesemodel tuning techniques may not be appropriate or workablefor various situations. For example, this study demonstrated
that the very thin top layer division did not improve much in

11 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

the two study sites with coarse soils. For fine soils, thisapproach may be viable because of very low infiltration rates,
but it may remain a challenge even for an experienced
modeler. A more realistic approach with this sort of watermodels may be to modify the hydraulic conductivity derived
from the hydraulic functions rather than assuming saturated
flow.
[33] As for soil water redistribution or the drying process, theoriginal water module appeared to dry soil too quickly compared to the HYDRUS water module (Figure 2). It is generally
accepted that physically based flow models that numericallysolve the Richards equation [Scanlon et al., 2002] or otheranalytical water flow equation [e.g., Bastiaanssen et al., 2007]
are superior to models not based on numerical solutions,especially when applied to arid regions. A specific difficultywith models that do not use the Richards equation may lie in
how to continuously move water between layers underunsaturated conditions. DAYCENT water module applies
Darcy’s Law, but the fluxes are adjusted by a “damping multiplier” probably in order to reasonably fit well with observation, which adds challenges to the model calibration process.
[34] The yearly water drainage volumes simulated by thetwo modules were similar (Table 3), implying that the two‐stage (infiltration‐redistribution) water flow module (in theoriginal DAYCENT) is appropriate for long‐term waterbalance calculations. Doltra and Muñoz [2010] also had
proved that both semi‐dynamic and transient dynamic water
modules could achieve very similar drainage simulation aslong as using hydraulic functions derived from same pedotransfer functions. However, the semi‐dynamic water modules cannot adequately describe the short‐term dynamics(Figure 2). In our study, the two‐stage water flow module of
DAYCENT produced more pulsed water drainage fluxes,
while the HYDRUS flow model simulated more continuous,
but smaller fluxes. Further comparison of internal soil water
fluxes also showed that, similarly as for bottom fluxes (i.e.,drainage), the HYDRUS module simulated more frequentand less intensive water fluxes than the DAYCENT waterflow module.
4.2. Soil Water Flow (and Drainage) and N Transport(and Leaching)[35] In this study, actual evapotranspiration (AET) from twowater modules was very similar (see Figure 4b and Table 3)
because the same potential ET and root distributions for water
uptake were used. The runoff from two modules was a small
portion of water balance (Table 3). So the difference of soil
water change was mainly caused by water drainage, whichremarkably occurred during late precipitation period (seeFigure 4). Therefore the characteristics of continuous, frequent,and less intensive water flow produced by the HYDRUS water
module would primarily be responsible to the increase of thetransport of soil surface nitrate downward and then its
leaching out of the soil profile in average and wet years. As
shown above, differences in the annual N leaching loss
simulated by the two water modules were largest during the
average year of 2003 (over 30%). For another average year
2001, the difference was still very large (18%) at the CP
site, while at the BF site it was much smaller, probably
because the dry conditions reduced drainage events anddrainage amounts (Table 4) to a critical threshold at whichfurther dry conditions would not significantly affect theG03040annual N leaching loss. These simulated differences decreasedin the very wet 2005 (Tables 3 and 4).
[36] Our analysis revealed that N leaching difference
between two water modules were actually related to the difference of water drainage at the daily time step rather thanlonger simulation intervals (Figure 7). In the DAYCENT
model, each water flow event will deplete the topsoil nitrateand stimulate soil mineralization. This process occurs primarily in the topsoil, and produces more soluble N, which inturn provides more nitrates for leaching. Because the plant
N uptake for these old‐growth trees was a small portion of
total N budget (less than 4 kg N ha−1 yr−1) and same amountof N deposition was received at each site, the difference of
simulated net N addition to soil was mainly the consequenceof large difference in N mineralization by two modules(Table 4). Figures 5a and 5d showed that with HYDRUS
water module, the N net addition to excess N pools was higher
than that with the original water module, partially caused bythis effect of water depletions. So relatively more frequentwater flow, predicted by HYDRUS compared to the original
water module, can transport more nitrate downward, so thatthe nitrate concentrations are higher in the bottom layer andlower in the top layer, except in dry years (see Figure 3 forcomparison of deep nitrate concentrations) when little
movement occurs due to lack of water.
[37] However, under intense precipitation conditions (e.g.,
the very wet 2005), initial very large leaching events produced by the models could exhaust soil inorganic nitrogenand leave less N for continuous loss later on (see Figures 6band 6e). Additionally, soil N mineralization, which canrestore soil mobile N pools, could be reduced due to moreanaerobic soil conditions in very wet years, as indicated byactually lower net N addition to soils in 2005. Also, in thewet 2005, the relative difference in the number of effective
drainage days produced by the two water flow modules was
significantly reduced compared to 2001 and 2003. Therefore, simulated N leaching in the very wet year 2005 wascomparable to those in the average years of 2001 (except for
the BF site) and 2003, although the drainage amount almost
doubled or tripled (Tables 3 and 4). The medium rainfall in
the average years caused soil water conditions and water
infiltration and drainage, which is favorable for N mineralization and downward movement and leaching. This isbecause generally, while N leaching increases with drainage(i.e., rainfall), it does not do so linearly [Simmelsgaard,
1998], mostly due to the limitation of available N.
[38] Lysimeter observations of high nitrate concentrations
at the beginning of the rainy season prior to large leaching
events (see Figure 3) imply that the actual nitrate leachingfrom the soil surface to deeper soil layers may be higher
than simulated. In other words, both modules mightunderestimate soil nitrate leaching losses, probably due tosoil nitrate diffusion and dispersion not explicitly includedin the current N module of DAYCENT. But generally it issmall compared to advective solute transport.
4.2.1. Transient (Top) Soil Water Dynamics
and N Gaseous Emissions[39] Nitrogen gaseous emissions are mainly produced by aseries of anaerobic processes as well as through the leaky pipeof nitrification in this study because ammonia volatilization
was not favored in these forests (pH 4–5 at CP and 5.8–6 at
BF [Grulke et al., 1998]). The NOx emission in the study sites

12 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

might be associated with nitrification and dominate the Nemission during dry season [e.g., Fenn et al., 1996], as shownin our simulations that N gaseous emission occurred
throughout the entire year (see Figures 6c, 6f, and 8b). But
apparently emissions during the rainy season (especially inthe spring) were more important, implying that simulating
soil water near‐saturation is critical in this ecosystem. The
late‐fall gaseous emissions were due to the early fall rainfallthat induced near‐saturation conditions and not much inorganic N leaching loss or plant uptake. With initial leaching in
the middle of the precipitation season, nitrate accumulation in
the topsoil layers was quickly reduced due to solute transport
to deep soil layers. Also the soil temperature, which was very
low at that time (i.e., winter period) at both study sites, constrained the production of trace N gases due to temperature
regulated soil NO−3 ‐N production through nitrification. Whenwater drainage and thus N leaching started to decline in thespring, another N emission peak (usually stronger) appeared,
because the rainfall season continued (so that soil watercontent was still high) and the weather was warmer, all ofwhich favored soil nitrification, denitrification, and othertrace N gaseous emission processes within the model.
[40] With the two‐stage water flow process considered in
the DAYCENT original water module, very few precipitation
events caused soil saturation, and N emissions were lower
than those predicted using the HYDRUS water module. Thisimplies that the approach used in the original water module
may be inappropriate for properly driving highly spontaneousand pulsed soil N gaseous processes in these two coarse‐soil
study sites. Even when the modeled soil profiles were divided
into layers as suggested (i.e., two very thin top layers 1 and3 cm thick), this modeling comparison study showed thatthe original approach was not as dynamical as the finite element algorithm used in HYDRUS for simulating water flow.
The HYDRUS approach produced moisture patterns thatallowed us to properly represent the pulses of soil N emissions at the daily time step, as theoretically and practically
expected.
4.2.2. Transient Soil Water or Semi‐dynamical Water
Modules for N Biogeochemical Cycle?
[41] As investigated in this study, many studies had alreadydemonstrated that Richards equation or other based numerical soil water flow modules are needed to capture dynamic
features and complexity of N biogeochemical processes
under highly variable water conditions [e.g., Gu and Riley,2010]. There are also many modeling studies showing that
more comprehensive representation of soil‐plant N processes[e.g., Doltra and Muñoz, 2010; Stolk et al., 2011] would be
required, especially for those N gaseous emission related. The
challenge is still remaining largely due to lack of continuous
and reliable field measurements on N losses, especially N
gaseous emissions, for model development and evaluation.
[42] Major concerns or obstacles regarding the application of
Richards equation‐based numerical soil water flow modules inlarge‐scale (both spatially and temporally) biogeochemical
models, such as DAYCENT, lie in the computational demandsof these modules, which usually require much longer computational time than tipping bucked modules [e.g., Crevoisieret al., 2009; de Jong and Bootsma, 1996; Scanlon et al.,2002]. This was also the case in our study. For examples,when running DAYCENT for 900 spin‐up years and then 109transient years at JMR sites with 6 soil layers, as shown in

G03040

Table 1 on a laptop with Intel® Core2 Duo CPU of 2.26 GHz
and RAM of 2 GB, with the original water module it took an
average of 150 s (10 duplicates), while with HYDRUS water
module it needed about 420 s. Then when running at CP
and BF sites with 10 soil layers, the average computation time
for the original semi‐dynamical water module was about 210–
230 s, but for the transient HYDRUS water module the run‐time needed rose up to 800–920 s. This could present a serious
obstacle when applying these sorts of models regionally[Del Grosso et al., 2006] and globally [Del Grosso et al.,2009]. However, considering the dependence of highly variable N biogeochemical processes on transient soil water contents and fluxes, it was necessary and worthwhile, as
demonstrated in this modeling comparison study, to incorporate a more dynamic soil water module to properly simulate
soil water flow and status. Recently, Crevoisier et al. [2009]
demonstrated a fast non‐iterative approach to numericallysolve the Richards equation. Such an approach may prove to bea very good alternative to tipping bucket water models for use
in biogeochemical model applications if experimentally
proven workable.

5. Conclusions[43] In the San Bernardino National Forest (SBNF) insouthern California, the mixed conifer forests have experienced the highest atmospheric N deposition in North America. To be able to track the fate of N, it is crucial that soil waterand its movement (and drainage) are appropriately modeled
for this Mediterranean type climate. This importance stems
from the fact that (1) soil water dynamics determine redox
conditions, an important factor controlling N gaseous emissions, and (2) downward water fluxes (drainage) transport
nitrate below the root zone and toward groundwater andultimately surface runoff, resulting in N leaching loss.
[44] After replacing the original DAYCENT’s semi‐
dynamic (two‐stage) water flow module with a finite element
based transient water flow module from HYDRUS‐1D(based on the solution of the Richards equation), daily soil
water contents were more dynamically simulated, better
representing the time series of soil water contents measured ata shallow soil site of an old‐growth mixed coniferous forests
surrounded by SBNF, especially during winter precipitation
periods. The two‐stage water flow algorithm used in theoriginal DAYCENT model did not allow the topsoil horizon
to remain above field capacity for a long period of time.
When comparing simulated deep soil nitrate concentrations
with the lysimeter measurements conducted at the two study
sites (CP and BF) from 2001 to 2004, simulations that usedthe HYDRUS water flow module produced a slightly higher
deep nitrate concentration prior to, but lower concentration after, the drainage period, compared to the original
DAYCENT module. However, both models captured thetime patterns of deep nitrate accumulation indicated by
the measurements.
[45] DAYCENT model with the original water module
estimated trace N gaseous losses 18–48% lower than when
the HYDRUS water flow module was incorporated, which
estimated soil N emissions of 2.4 ∼ 6.8 kg N ha−1 year−1 at thehigh‐N CP site and 1.7 ∼ 3.8 kg N ha−1 year−1 at the low‐N
BF site. This simulated N gaseous emission difference
between two water modules largely due to the transient water

13 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

flow module simulated more favorable water conditions (nearsaturation) for denitrification relevant N gaseous emissions.
[46] The semi‐dynamic water module of DAYCENT
generally simulated less frequent but higher volume (i.e.,
more spontaneous) water fluxes and drainage compared to
the HYDRUS water module. Although the annual waterdrainage simulated by the two modules was not significantly
different, the HYDRUS module predicted more N mineralization and nitrate downward movement and thus higher
leaching losses. This is because frequent but small waterflows promoted continuous nitrate transport and likely soil
N mineralization at topsoils. With contrasted drainage time
patterns, the difference of N leaching between the twomodels was about 14 ∼ 18 kg N ha−1 yr−1 at the high‐N CP
site, while about 5 ∼ 12 kg N ha−1 yr−1 at the low‐N BF site,
during normal years (with 60 ∼ 80 cm annual precipitation).
During either dry (annual precipitation of less than 30 cm)
or wet (about 120 cm) years, the differences in model predictions in N leaching were relatively small.
[47] This modeling study suggests that estimation of waterbalance over time periods longer than one month will not
likely be impacted whether semi‐dynamical or fully transient water modules are used. However, it appears thatnumerical solutions based on the transient Richards equation
are preferred when soil water flow and solute transport arecoupled and predictions of soil‐plant N cycling are needed
at time scales of days or less.
[48] Acknowledgments. We thank the Natural Resource Ecology Laboratory at Colorado State University for permitting use of their DAYCENT
model, training, and providing code for this work. Michael F. Allen andKuni Kitajima in the Center for Conservation Biology (CCB), University ofCalifornia at Riverside (UCR), shared their long‐term monitoring data at
James San Jacinto Mountains Reserves with us for the model comparison.
CCB at UCR also provided logistical support for this work. This project
was funded by National Science Foundation (DEB 04‐21530).

References
Allen, M. F., et al. (2007), Soil sensor technology: Life within a pixel,
BioScience, 57, 859–867, doi:10.1641/B571008.
Arbaugh, M. J., et al. (1999), Simulated effects of N deposition, ozone injury
and climate change on a forest stand in the San Bernardino Mountains,in Oxidant Air Pollution Impacts in the Montane Forests of Southern
California: A Case Study of the San Bernardino Mountains, edited by
P. R. Miller and J. R. McBride, pp. 353–372, Springer, New York,
doi:10.1007/978-1-4612-1436-6_17.
Bastiaanssen, W. G. M., et al. (2007), Twenty‐five years modeling irrigatedand drained soils: State of the art, Agric. Water Manage., 92, 111–125,
doi:10.1016/j.agwat.2007.05.013.
Chen, D. L., et al. (2008), N2O emissions from agricultural lands: A synthesisof simulation approaches, Plant Soil, 309, 169–189, doi:10.1007/s11104008-9634-0.
Crevoisier, D., et al. (2009), Evaluation of the Ross fast solution of
Richards’ equation in unfavourable conditions for standard finite element methods, Adv. Water Resour., 32, 936–947, doi:10.1016/j.advwatres.2009.03.008.
de Jong, R., and A. Bootsma (1996), Review of recent developments in soil
water simulation models, Can. J. Soil Sci., 76, 263–273, doi:10.4141/
cjss96-033.
Del Grosso, S. J., W. J. Parton, A. R. Mosier, D. S. Ojima, A. E. Kulmala,and S. Phongpan (2000), General model for N2O and N2 gas emissions
from soils due to denitrification, Global Biogeochem. Cycles, 14,1045–1060, doi:10.1029/1999GB001225.
Del Grosso, S. J., et al. (2001), Simulated interaction of carbon dynamics andnitrogen trace gas fluxes using the DAYCENT model, in Modeling Carbon
and Nitrogen Dynamics for Soil Management, edited by M. Shaffer et al.,pp. 303–332, CRC Press, Boca Raton, Fla.
Del Grosso, S., et al. (2002), Simulated effects of dryland cropping intensification on soil organic matter and greenhouse gas exchanges using theG03040

DAYCENT ecosystem model, Environ. Pollut., 116, S75–S83, doi:10.1016/
S0269-7491(01)00260-3.
Del Grosso, S. J., et al. (2005), DAYCENT model analysis of past and contemporary soil N2O and net greenhouse gas flux for major crops in theUSA, Soil Tillage Res., 83, 9–24, doi:10.1016/j.still.2005.02.007.
Del Grosso, S. J., et al. (2006), DAYCENT national‐scale simulations of
nitrous oxide emissions from cropped soils in the United States, J. Environ.
Qual., 35, 1451–1460, doi:10.2134/jeq2005.0160.
Del Grosso, S. J., et al. (2009), Global scale DAYCENT model analysis ofgreenhouse gas emissions and mitigation strategies for cropped soils,
Global Planet. Change, 67, 44–50, doi:10.1016/j.gloplacha.2008.12.006.
del Prado, A., et al. (2006), N2O and NO emissions from different N sources
and under a range of soil water contents, Nutr. Cycl. Agroecosyst., 74,
229–243, doi:10.1007/s10705-006-9001-6.
Di Luzio, M., et al. (2008), Constructing retrospective gridded daily precipitation and temperature datasets for the conterminous United States, J. Appl.
Meteorol. Climatol., 47, 475–497, doi:10.1175/2007JAMC1356.1.
Dobbie, K. E., and K. A. Smith (2003), Nitrous oxide emission factors for
agricultural soils in Great Britain: The impact of soil water‐filled porespace and other controlling variables, Global Change Biol., 9, 204–218,
doi:10.1046/j.1365-2486.2003.00563.x.
Doltra, J., and P. Muñoz (2010), Simulation of nitrogen leaching from a
fertigated crop rotation in a Mediterranean climate using the EU‐Rotate_N
and Hydrus‐2D models, Agric. Water Manage., 97, 277–285, doi:10.1016/
j.agwat.2009.09.019.
Fenn, M. E., et al. (1996), Evidence for nitrogen saturation in the San Bernardino Mountains in southern California, For. Ecol. Manage., 82, 211–230,
doi:10.1016/0378-1127(95)03668-7.
Fenn, M. E., et al. (2003), Nitrogen emissions, deposition, and monitoring in
the western United States, BioScience, 53, 391–403, doi:10.1641/00063568(2003)053[0391:NEDAMI]2.0.CO;2.
Fenn, M. E., et al. (2008), Empirical and simulated critical loads for nitrogen deposition in California mixed conifer forests, Environ. Pollut., 155,492–511, doi:10.1016/j.envpol.2008.03.019.
Feyen, J., et al. (1998), Modelling water flow and solute transport in heterogeneous soils: A review of recent approaches, J. Agric. Eng. Res., 70,
231–256, doi:10.1006/jaer.1998.0272.
Frolking, S. E., et al. (1998), Comparison of N2O emissions from soils at threetemperate agricultural sites: Simulations of year‐round measurements byfour models, Nutr. Cycl. Agroecosyst., 52, 77–105, doi:10.1023/
A:1009780109748.
Grulke, N. E., et al. (1998), Ozone exposure and nitrogen deposition lowers root biomass of ponderosa pine in the San Bernardino Mountains,
California, Environ. Pollut., 103, 63–73, doi:10.1016/S0269-7491(98)00130-4.
Gu, C., and W. J. Riley (2010), Combined effects of short term rainfall
patterns and soil texture on soil nitrogen cycling–A modeling analysis,
J. Contam. Hydrol., 112, 141–154, doi:10.1016/j.jconhyd.2009.12.003.
Gu, C., F. Maggi, W. J. Riley, G. M. Hornberger, T. Xu, C. M. Oldenburg,
N. Spycher, N. L. Miller, R. T. Venterea, and C. Steefel (2009), Aqueousand gaseous nitrogen losses induced by fertilizer application, J. Geophys.
Res., 114, G01006, doi:10.1029/2008JG000788.
Hendriks, R. F. A., et al. (2007), Predicting soil subsidence and greenhousegas emission in peat soils depending on water management with the
SWAP‐ANIMO model, paper presented at First International Symposium
on Carbon in Peatlands, Wageneingen Univ., Wageningen, Netherlands,
15–18 Apr.
Jones, J. W., et al. (2003), The DSSAT cropping system model, Eur. J.
Agron., 18, 235–265, doi:10.1016/S1161-0301(02)00107-7.
Kitajima, K., K. E. Anderson, and M. F. Allen (2010), Effect of soil temperatureand soil water content on fine root turnover rate in a California mixed conifer
ecosystem, J. Geophys. Res., 115, G04032, doi:10.1029/2009JG001210.
Lee, D. H., and L. M. Abriola (1999), Use of the Richards equation in land
surface parameterizations, J. Geophys. Res., 104, 27,519–27,526, doi:10.1029/
1999JD900951.
Pachepsky, Y., et al. (2006), Information content and complexity of simulated
soil water fluxes, Geoderma, 134, 253–266, doi:10.1016/j.geoderma.
2006.03.003.
Parton, W. J., et al. (1988a), Rates and pathways of nitrous‐oxide production ina shortgrass steppe, Biogeochemistry, 6, 45–58, doi:10.1007/BF00002932.
Parton, W. J., et al. (1988b), Dynamics of C, N, P and S in grassland soils ‐ Amodel, Biogeochemistry, 5, 109–131, doi:10.1007/BF02180320.
Parton, W. J., et al. (1998), DAYCENT and its land surface submodel:
Description and testing, Global Planet. Change, 19, 35–48,
doi:10.1016/S0921-8181(98)00040-X.
Parton, W. J., E. A. Holland, S. J. D. Grosso, M. D. Hartman, R. E. Martin,
A. R. Mosier, D. S. Ojima, and D. S. Schimel (2001), Generalized model forNOx and N2O emissions from soils, J. Geophys. Res., 106, 17,403–17,419,
doi:10.1029/2001JD900101.

14 of 15

G03040

YUAN ET AL.: TRANSIENT SOIL WATER ON N LOSSES

Ranatunga, K., et al. (2008), Review of soil water models and their applications in Australia, Environ. Model. Softw., 23, 1182–1206, doi:10.1016/j.
envsoft.2008.02.003.
Ritchie, J. T. (1998), Soil water balance and plant water stress, in Understanding Options for Agricultural Production, edited by G. Y. Tsuji
et al., pp. 41–54, Kluwer Acad., London.
Rodriguez‐Iturbe, I., P. D’Odorico, F. Laio, L. Ridolfi, and S. Tamea(2007), Challenges in humid land ecohydrology: Interactions of watertable and unsaturated zone with climate, soil, and vegetation, Water
Resour. Res., 43, W09301, doi:10.1029/2007WR006073.
Saxton, K. E., et al. (1986), Estimating generalized soil‐water characteristics from texture1, Soil Sci. Soc. Am. J., 50, 1031–1036, doi:10.2136/
sssaj1986.03615995005000040039x.
Scanlon, B. R., M. Christman, R. C. Reedy, I. Porro, J. Šimůnek, and G. N.
Flerchinger (2002), Intercode comparisons for simulating water balanceof surficial sediments in semiarid regions, Water Resour. Res., 38(12),1323, doi:10.1029/2001WR001233.
Schaap, M. G., et al. (2001), A computer program for estimating soil
hydraulic parameters with hierarchical pedotransfer functions, J. Hydrol.,
251, 163–176, doi:10.1016/S0022-1694(01)00466-8.
Shao, Y., and A. Henderson‐Sellers (1996), Modeling soil moisture: A
Project for intercomparison of land surface parameterization schemes
phase 2(b), J. Geophys. Res., 101, 7227–7250, doi:10.1029/95JD03275.
Simmelsgaard, S. E. (1998), The effect of crop, N‐level, soil type anddrainage on nitrate leaching from Danish soil, Soil Use Manage., 14,
30–36, doi:10.1111/j.1475-2743.1998.tb00607.x.
Šimůnek, J., et al. (2008), Development and applications of the HYDRUS
and STANMOD software packages and related codes, Vadose Zone J., 7,
587–600, doi:10.2136/vzj2007.0077.

G03040Smith, K. A., et al. (1998), Effects of temperature, water content and nitrogen fertilisation on emissions of nitrous oxide by soils, Atmos. Environ.,32, 3301–3309, doi:10.1016/S1352-2310(97)00492-5.
Stolk, P. C., et al. (2011), Simulation of daily nitrous oxide emissions frommanaged peat soils, Vadose Zone J., 10, 156–168, doi:10.2136/
vzj2010.0029.
Twarakavi, N. K. C., et al. (2009), An objective analysis of the dynamic
nature of field capacity, Water Resour. Res., 45, W10410, doi:10.1029/
2009WR007944.
Vargas, R., and M. F. Allen (2008), Dynamics of fine root, fungal rhizomorphs, and soil respiration in a mixed temperate forest: Integrating sensors
and observations, Vadose Zone J., 7, 1055–1064, doi:10.2136/
vzj2007.0138.
Wösten, J. H. M., et al. (2001), Pedotransfer functions: Bridging the gap
between available basic soil data and missing soil hydraulic characteristics, J. Hydrol., 251, 123–150, doi:10.1016/S0022-1694(01)00464-4.
M. E. Fenn, Pacific Southwest Research Station, USDA Forest Service,
4955 Canyon Crest Dr., Riverside, CA 92507, USA.
T. Meixner, Department of Hydrology and Water Resources, Universityof Arizona, 845 N. Park Ave., Tucson, AZ 85710, USA.
J. Šimůnek, Department of Environmental Sciences, University ofCalifornia, 900 University Ave., A135 Bourns Hall, Riverside, CA 92521,
USA
F. Yuan, Institute of Arctic Biology, University of Alaska Fairbanks,
902 Koyukuk Dr., Fairbanks, AK 99775, USA. (fyuan@alaska.edu)

15 of 15