Agricultural Water Management 160 (2015) 64–75

Contents lists available at ScienceDirect

Agricultural Water Management
journal homepage: www.elsevier.com/locate/agwat

Evaluating the impact of groundwater on cotton growth and root zone
water balance using Hydrus-1D coupled with a crop growth model
Ming Han a,b , Chengyi Zhao a,∗ , Jirka Šimůnek c , Gary Feng daState Key Laboratory of Oasis Desert and Ecology, Xinjiang Institute of Ecology and Geography, CAS, Urumqi 830011, China
University of Chinese Academy of Sciences, Beijing 100049, Chinac
Department of Environmental Sciences, University of California Riverside, Riverside, CA 92506, USAdGenetics and Sustainable Agriculture Research Unit, USDA-ARS, P. O. BOX 5367, 810 Highway 12 East, Mississippi State, MS 39762, USAba r t i c l e

i n f oArticle history:
Received 7 December 2014
Received in revised form 27 May 2015
Accepted 30 June 2015
Keywords:
CottonGroundwater
Yield
Evapotranspiration
Hydrus-1D
Crop growtha b s t r a c t
Groundwater is an important factor that needs to be considered when evaluating the water balance ofthe soil-plant-atmosphere system and the sustainable development of arid oases. However, the impact
of shallow groundwater on the root zone water balance and cotton growth is not fully understood. In thisstudy, we have ﬁrst analyzed the inﬂuence of the groundwater table depth on the seasonal maximum
leaf area index of cotton, the average seasonal water stress, cotton yield, actual transpiration, actual
evaporation, and capillary rise using experimental data collected at the Aksu water balance station, in
Xinjiang, northwest of China and the Hydrus-1D variably-saturated soil water ﬂow model coupled with asimpliﬁed crop growth model from SWAT. The coupled model has been ﬁrst calibrated and validated usingﬁeld observations of soil water content, leaf area index, cotton height, the above ground biomass, andcotton yield comparisons between measured and modeled variables have shown a reasonable agreement
for all variables. Additionally, with a validated model, we have carried out numerical experiments fromwhich we have concluded that groundwater is a major water resource for cotton growth in this region.
The capillary rise from groundwater contributes almost 23% of crop transpiration when the averagegroundwater depth is 1.84 m, which is the most suitable groundwater depth for this experimental site. We
have concluded that cotton growth and various components of the soil water balance are highly sensitive
to the groundwater table level. Different positions of the groundwater table showed both positive andnegative effects on cotton growth. Likewise, cotton growth has a signiﬁcant impact on the capillaryrise from groundwater. As a result, groundwater is a crucial factor that needs to be considered whenevaluating agricultural land management in this arid region. The updated Hydrus-1D model developedin this study provides a powerful modeling tool for evaluating the effects of the groundwater table on
local land management.
© 2015 Z. All rights reserved.

1. Introduction
The Aksu oasis, located in Southern Xinjiang, is a major area in
China that produces cotton. Cotton is a leading cash crop in thisregion. The main climatic feature of this arid oasis is the extremely
high ratio of evaporation to precipitation. Water scarcity is oneof the most critical constraints on the sustainable production ofcotton. However, inappropriate land management has resulted inmany problems in this ecologically fragile region. Deteriorationof soils, low water use efﬁciency, and salt build-up in the soil-∗ Corresponding author.
E-mail address: zcy@ms.xjb.ac.cn (C. Zhao).
http://dx.doi.org/10.1016/j.agwat.2015.06.028
0378-3774/© 2015 Z. All rights reserved.

groundwater system threaten the sustainable development of thisarid oasis (Kang et al., 2012).
To solve these problems, many researchers have studied water
ﬂow, solute transport, and root water uptake in this region andobtained a certain conceptual understanding of the system (Houet al., 2009; Hu et al., 2009; Kang et al., 2012; Shen et al., 2012). Characteristics of the cotton root distribution under different amountsof irrigation (Hu et al., 2009; Shen et al., 2012), relationships
between the cotton production and irrigation strategies (Yan et al.,2009), and the spatial and temporal distribution of water and salts
in the soil (Hou et al., 2009; Liu et al., 2012) have been studied indetail using both ﬁeld experiments and numerical modeling. However, no research has been carried out to evaluate the impact ofgroundwater on cotton growth and the root zone water balance,

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

which are important issues that need to be considered to design asustainable land management of this arid region.
It has been shown that groundwater is one of the signiﬁcantdrivers that impact soil water dynamics and vegetation growth.
Shallow groundwater has both positive and negative effects on
agricultural land management. Speciﬁcally, capillary rise fromgroundwater can play an important role in contributing to the cropwater use, which is very important in arid and semi-arid regions(Ayars et al., 2006; Jorenush and Sepaskhah, 2003). Shallow groundwater could be viewed as a potential water resource for crop use.
For example, Soppe and Ayars (2003) found that groundwater contributed up to 40% of daily water used by the safﬂower crop, andon a seasonal basis, 25% of the total crop water use originated fromthe groundwater, when a groundwater table was maintained at a
depth of 1.5 m in weighing lysimeters. On the other hand, excessive
capillary rise from groundwater can also cause a lot of problems,such as an increase in soil evaporation and the risk of soil salinization (Soylu et al., 2014). Shallow groundwater can negatively affectplant physiological functions and can even cause plant mortality(Zaidi et al., 2004). When the water table is persistently too close
to the soil surface during the growing season, the oxygen stresson roots may increase, negatively affecting the crop growth (Soylu
et al., 2014).
There are relatively few process-based models that can consider
groundwater as a part of the soil-vegetation-atmosphere system.
Existing models, such as SWAP (Kroes et al., 2008), DSSAT (Joneset al., 2003), and EPIC (Williams et al., 1989), have been widely used
to quantify production potentials, estimate soil moisture dynamics, and evaluate the performance of irrigation systems and waterproductivity at both ﬁeld and regional scales. However, many ofthese models usually oversimplify the impact of groundwater.
Recently, a new hydrological model has been developed that caneffectively combine processes in the vadose zone and groundwater (Twarakavi et al. 2008). In order to provide a balance between
computational efﬁciency and accuracy, Twarakavi et al. (2008) coupled the Hydrus-1D model (Šimůnek et al., 2008) with MODFLOW(Harbaugh et al., 2000) to simulate one-dimensional (1D) waterﬂow in the vadose zone and three-dimensional (3D) water ﬂow
in groundwater (Twarakavi et al., 2008). However, although waterﬂow is rigorously described in this coupled model, the descriptionof the plant growth and how the plant growth is affected by existing moisture, salinity, and temperature conditions is dramaticallysimpliﬁed.
In this study, we implemented a crop growth model into
Hydrus-1D to study the impact of groundwater on cotton growthand soil water dynamics. The crop growth module is based on
the crop growth model of the SWAT model (Williams et al., 1989;
Williams and Singh, 1995). This model uses a uniﬁed approach tosimulate the growth for more than 80 types of crops and has been
widely used in various eco-hydrological studies (Williams et al.,2006). In this paper, the SWAT crop growth model was simpliﬁed
to simulate (only) cotton growth, including its phenological development based on daily accumulated heat units. The model also
predicts a harvest index for partitioning grain yield and potential
biomass accumulation and can consider water stresses and accordingly adjust the crop growth.
The main objectives of this study thus were to: (1) develop amodeling approach to simultaneously estimate crop production,soil moisture dynamics, evaporation, and transpiration by implementing the SWAT crop growth model into Hydrus-1D, (2) tocalibrate and validate the resulting model using experimental ﬁeld
observations, and (3) to use this new tool to evaluate how thedepth of the groundwater table inﬂuences the cotton root zone
water balance and cotton growth. This new model could be further
incorporated into the Hydrus ﬂow package for MODFLOW, which
would overcome the simpliﬁcation of the description of the soil-

65

Fig. 1. A coupling scheme for the Hydrus-1D and crop growth models.

groundwater-crop-atmosphere system in this eco-hydrologicalmodel.
2. Materials and methods
2.1. Model structure
The structure of the coupled model is shown in Fig. 1. The crop
growth and the soil water content dynamics are simulated by the

66

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

crop growth model and Hydrus-1D, respectively, using different
time steps in the two modules. While the crop growth model uses
daily time steps, Hydrus-1D requires smaller time steps in order toproperly solve the highly nonlinear Richards equation. During eachday, Hydrus-1D may perform multiple smaller time steps to simulate soil water dynamics. The two models exchange information
about the potential transpiration, potential evaporation, and actual
transpiration at daily time steps. First, the crop growth model calculates the potential evaporation and potential transpiration eachday based on climate data, crop leaf area index (LAI), and the crop
canopy height. The potential evaporation and potential transpiration are then assigned as an upper atmospheric boundary conditionin Hydrus-1D. The Hydrus-1D model then simulates the soil water
dynamics, actual evaporation, and actual transpiration for one day.
At the end of the day, the crop growth model receives informationfrom Hydrus-1D on actual evaporation and actual transpiration andcalculates the water stress suffered by the crop, the actual LAI increment, the actual biomass increment, the actual crop canopy height,
and the harvest index for the current day. The updated LAI and crop
canopy height are then used to calculate potential evaporation andtranspiration for the next day.
While the entire Hydrus-1D model was used in this modelingframework, the SWAT crop growth model that was implemented
into Hydrus-1D was simpliﬁed. The original crop growth model
considers processes, such as crop interception of solar radiation, theconversion of intercepted light to biomass, the division of biomass
into roots, above ground biomass, and economic yield, root growth,water use, and nutrient uptake. Potential plant growth is simulated using daily time steps and is constrained by the minimum ofﬁve stress factors (water, nitrogen, phosphorus, temperature, andaeration). The water use of the crop and root growth were not considered in the adapted crop growth model since these processes are
directly simulated by Hydrus-1D. Also, the nutrient uptake was not
simulated and the potential plant growth was only constrained byavailable water. In the coupled model, the water stress is calculated
as one minus the ratio of actual and potential transpiration during
each day. Detailed descriptions of water ﬂow and root water uptake
processes can be found in the Hydrus-1D manual (Šimůnek et al.,2005), while a detailed description of the crop growth model canbe found in the SWAT manual (Neitsch et al., 2005). The simpliﬁed
crop growth model implemented into Hydrus-1D and adapted inthis study is described below in Sections 2.2 and 2.3.

i


HUIi =

HUk

k=1

(2)

PHU

where HUIi is the heat unit index for day i and PHU is the number of potential heat units required for maturity. The value of PHU
depends on a particular crop and is usually provided by the user ofthe model.
The potential crop growth is modeled by simulating the development of the leaf area index, light interception, and conversion ofintercepted light into biomass. The interception of solar radiation
is estimated as (Neitsch et al., 2005; Williams et al., 1989):
Hphosyn,i = 0.5 × Hday,i × (1 − exp(−kl × LAIi ))(3)

where Hphosyn,i is the intercepted photosynthetic active radiation
on day I (MJ m−2 ), Hday,i is the incident total solar radiation on dayi (MJ m−2 ), kl is the light extinction coefﬁcient, and LAIi is the leaf
area index on day i.
The potential increase in biomass for a particular day is estimated as (Neitsch et al., 2005; Williams et al., 1989):
bioi = RUE × Hphosyn,i(4)where bioi is the daily potential increase in biomass in (kg ha−1 )on day i and RUE is the crop radiation-use efﬁciency (10−1 g/MJ).
An increase in the potential leaf area on day i is calculated as(Neitsch et al., 2005):

LAIi = frLAImx,i − frLAImx,i−1

× LAImx × (1 − exp (5 × (LAIi−1 − LAImx )))(5)

where LAIi is the potential LAI increment on day i, LAImx is the
maximum leaf area index, and frLAImx,i and frLAImx,i−1 are the fractions of the plant’s maximum leaf area index on days i and i − 1,
respectively:
frLAImx,i =HUIi
HUIi + exp (l1 − l2 × HUIi )(6)

where l1 and l2 are shape coefﬁcients.
The potential harvest index on a particular day is calculated as(Neitsch et al., 2005):

2.2. Crop growth model
2.2.1. Potential crop growth
The heat unit theory has proven to be a reliable predictor of thecrop physiological growth. In this theory, the growth of the cropis based on the accumulation of daily heat units, which may becomputed as follows (Neitsch et al., 2005; Williams et al., 1989):
HUi = Tav,i − Tbase when Tav,i > Tbase

A heat unit index (HUI), which ranges from 0 at planting to 1 at
physiological maturity, is computed as follows (Neitsch et al., 2005;
Williams et al., 1989):(1)

HU = 0 when T av,i ≤ Tbase

HIi = HIopt ×(7)

where HIi is the potential harvest index on day i and HIopt is thepotential harvest index for the crop at maturity for ideal growing
conditions.
2.2.2. Actual crop growth
The water stress is the only crop growth constraint consideredby the adapted model. The water stress is equal to 0 for optimal
water contents and increases to 1 as soil water content approaches
the wilting point. The water stress is obtained by comparing actualand potential plant transpiration (Neitsch et al., 2005):
wstrsi = 1 −where HUi is the number of accumulated heat units on day i(heat units), Tav,i is the mean daily temperature on day i (◦ C), andTbase is the crop’s base temperature for growth (◦ C). No growthoccurs when mean daily temperature is at or below base temperature.

100 × HUIi
100 × HUIi + exp(11.1 − 10 × HUIi )Et,i,act
Et,i(8)

where wstrsi is the water stress that crop suffers on day i, Et,iis the potential transpiration on day i (mm), calculated using thePenman–Monteith equation (Allen et al., 1998), and Et,i, act is theactual transpiration on day i (mm), calculated using Hydrus-1D.

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

The potential biomass increment predicted using Eq. (4) isadjusted daily as follows (Neitsch et al., 2005):
bioact,i = bioi × (1 − wstrsi )(9)where bioact, i is the actual biomass accumulated on day i(kg ha−1 ).
The actual biomass on day i is calculated as:
bioi = bioact,i + bioi−1(10)where bioi is the biomass at the end of day i and bioi−1 is thebiomass at the end of day i − 1 (kg ha−1 ).
The above ground crop plant biomass is calculated as (Neitsch
et al., 2005):
bioag,i = (1 − frroot,i ) × bioi

(11)

where bioag, i is the crop above ground biomass on day i (kg ha−1 )and frroot, i is the fraction of total biomass in roots on day i (−).
The fraction of the total biomass in roots on a given day is calculated as (Neitsch et al., 2005):
frroot,i = 0.4 − 0.2 × HUIi

LAIi,act = LAIi ×

1 − wstrsi(13)

where LAIi,act is the actual leaf area index added on day i.
The increase in LAI is simulated as a function of accumulated
heat units. From emergence to the start of leaf decline, LAI is estimated as follows (Neitsch et al., 2005):
LAIi = LAIi−1 + LAIi whenHUIi < dlai
1−HUI
LAIi = LAIms × ( 1−PHUi )whenHUIi > dlaiwhere

LAIi and(14)
LAIi−1 are the leaf area indexes on day i and i − 1, respectively, Anddlai is the fraction of the growing season when the leaf area index
starts declining.
The canopy height on a particular day is calculated as (Neitsch
et al., 2005):
hc,i = hc,mx ×

frLAImx,i(15)

where hc,i is the canopy height on day i (m) and hc,mx is the maximum crop canopy height (m). Once the maximum canopy heightis reached, hc will remain constant until crop is harvested.
The effect of water deﬁcit on the potential harvest index calculated using Eq. (7) is accounted for using the following relationship(Neitsch et al., 2005):
HIac,i = (HI − HImin ) ×

wu,i
wu,i + exp(6.13 − 0.883 × wu,i )(16)

where HIac, i is the actual harvest index on day i, HImin is theminimum harvest index allowed, and  wu,i is the water deﬁciency
factor, which is calculated as (Neitsch et al., 2005):iwu,i = 100 ×

Ea,k

k=1
i


on the day of harvest (kg ha−1 ), and HIac is the actual harvest index
on the day of harvest. Parameters of the simpliﬁed crop growthmodel are listed in Table 1.
2.3. The Hydrus-1D model
2.3.1. Water ﬂow equations
Soil water movement in the experimental ﬁeld was simulatedusing the Hydrus-1D model, which numerically solves the governing ﬂow equation (Richards, 1931):

∂
∂
∂h
=
[k()
]−S
∂t
∂x
∂x

(17)
E0,k

k=1

where Ea,k is the actual evapotranspiration on day k, calculated byHydrus-1D, and E0,k is the potential evapotranspiration on day k,
calculated using the Penman–Monteith equation.
The crop yield is calculated based on the harvest index and theabove ground biomass as (Neitsch et al., 2005):(18) yld = bioag × HIac whenHIac ≤ 1
where yld is the crop yield (kg/ha), bioag is the above groundbiomass on the day of harvest (kg ha−1 ), bio is the total crop biomass(19)

where  is the volumetric soil water content (cm3 cm−3 ), h is
the pressure head (cm), S is a sink term accounting for root water
uptake (cm−1 ), and k() is the unsaturated hydraulic conductivity
function (cm day−1 ). Soil hydraulic properties were described using
the van Genuchten–Mualem analytical functions (van Genuchten,1980):(12)

The potential leaf area added on a particular day in Eq. (5) is also
adjusted daily as (Neitsch et al., 2005):

67 (h) =

⎧
⎨ +
r

⎩

s − r
1 + |˛h|;h < 0n m

 s; h ≥ 0

 

1/m

K  = Ks Sel 1 − 1 − Se

m

2



 − r , m = 1 − 1/n, n > 1
Se = 
s − r(20)

where  s is the saturated water content (cm3 cm−3 ),  r is theresidual water content (cm3 cm−3 ), Ks is the saturated hydraulicconductivity (cm day−1 ), ␣ and n are the shape parameters, and l is
the pore-connectivity and tortuosity parameters.
The Hydrus-1D model solves the Richards equation using theGalerkin ﬁnite-element method. A detailed description of themodel can be found in the Hydrus-1D manual (Šimůnek et al.,2005).
2.3.2. Root water uptake
In Hydrus-1D, the volume of water that is removed from a unit
volume of soil during a unit time due to root water uptake is deﬁned
as (Šimůnek et al., 2005):
Sw (h, x) = ˛(h) × RLT(x) × Et,i(21)

where Et,i is the potential transpiration of crop on day i, which is
obtained from the crop growth model, and ␣(h) is the soil water
stress response function (dimensionless) of Feddes et al. (1978).
The parameters of the water stress response function in the Feddes
model for cotton are: h1 = −10 cm, h2 = −25 cm, h3max = −200 cm,
h3min = −600 cm, and h4 = −14000 cm (Forkutsa et al., 2009; Wanget al., 2014). For h < h4 (the wilting point pressure head) and h > h1
(the anaerobiosis point pressure head), water uptake is assumed
to be zero. Water uptake is considered optimal between pressure heads h2 and h3 , whereas for pressure heads between h3 andh4 (or h1 and h2 ), water uptake decreases (or increases) linearlywith h (Šimůnek and Hopmans, 2009). Note that the h3 value isinterpolated between h3min and h3max depending on the potential
transpiration rate Tp (Feddes et al., 1978). RLT(x) is the normalized
root water uptake spatial distribution, which is deﬁned as:b(x)RLT(x) = 
b(x)

(22)

68

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

Table 1
Parameters of the simpliﬁed crop growth model.
Parameter

Units

Description

RUE
HIopt
LAImx
dlai
frgrw1(kg/ha)/(J/m2 )(kg/ha)/(kg/ha)
m2 /m2––laimx1–

frgrw2–laimx2

–

hc,mx
Tbase
HImin
kl
PHU

m
◦
C
–
–
–

Biomass-energy ratio (Eq. 4)
Potential harvest index for the crop at maturity (Eq. 7)
Maximum leaf area index (Eq. 5)
A fraction of the growing season when leaf area declines (Eq. 14)
A fraction of the growing season corresponding to the 1st point on the optimal leaf area development curve, used to
calculatel(Eq. 6)
A fraction of the maximum leaf area index corresponding to the 1st point on the optimal leaf area development curve,
used to calculatelin Eq. 6
A fraction of the growing season corresponding to the 2nd point on the optimal leaf area development curve, used to
calculatel(Eq. 6)
A fraction of the maximum leaf area index corresponding to the 2nd point on the optimal leaf area development
curve, used to calculatel(Eq. 6)
Crop maximum canopy height (Eq. 15)
Crop’s base temperature (Eq. 1)
The minimum harvest index allowed for the crop (Eq. 16)
A light extinction coefﬁcient (Eq. 3)
Potential heat units required for maturity (Eq. 2)

Table 2
Soil bulk density and particle size distribution.
Depth (cm)

0–10
10–20
20–30
30–40
40–60
60–80
80–100

Bulk density (g/cm3 )

1.52
1.58
1.59
1.58
1.56
1.56
1.54

Table 4Irrigation schedule in years 2007 and 2008.

Soil particle size distribution (%)

Irrigation date

Irrigation amount (mm)<0.002 mm

0.002–0.02 mm>0.02 mm

5
5
5
11
16
6
4

46
48
49
74
80
48
39

49
47
46
15
4
46
57

02/16/2007
07/10/2007
07/23/2007
08/16/2007
11/07/2007
02/02/2008
06/01/2008
06/23/2008
07/09/2008
08/06/2008
12/03/2008

200
150
150
150
200
200
150
150
150
150
200

where b(x) is the root water uptake distribution function. There aremany ways how the b(x) function can be expressed (e.g., constant,linear, exponential with depth), and Hydrus-1D does not limit its
users in any way. The actual transpiration on day i is then calculated
as:
Et,i,act =

N
i 


Sw,k (h, x) × tk(23)

k=i−1 x=1

where k is the kth time step to solve the Richard equation between
days i − 1 and i, tk is the duration of the time step for the kth time
step, and N is the nodal number.
2.4. Field experiment and model parameterization
2.4.1. Site location
The present study was conducted at the Aksu National Water
Balance Station (latitude: 40.37◦ N, 80.45◦ E, 1028.0 m a.s.l.), locatedin the Aksu Oasis, northwest of the Tarim Basin in the Xinjiang
province of northwestern China. The study site was located in a
typical cotton cultivation area, which has inland arid climate conditions. Mean annual precipitation is approximately 45.7 mm. Most
precipitation occurs from June to October, accounting for 65% of
the annual precipitation. Mean annual evaporation from the freewater surface is 2500 mm. The average annual temperature in the
study area is approximately 8 ◦ C, and the annual accumulated temperature higher than 10 ◦ C is 4428 ◦ C. The average depth of thegroundwater table is approximately 2 m. Soil texture at the experi-

ment site is medium loam. The soil texture and the soil bulk densityof the experimental site is shown in Table 2.
2.4.2. Field experiment design
The experiment, conducted in years of 2007 and 2008, wasdesigned to observe interactions between the soil, soil water,
groundwater, and crop. The experimental plot size was 150 × 90 m.
ZhongMian 49, a cotton species widely grown in China,was usedin the ﬁeld trial. The sowing date was on April 4th and the harvesting date in early November. The phenological phases of thecotton growth were closely monitored (Table 3). The ﬁeld was
ﬂood-irrigated during the growing season 3 times in 2007 and 4
times in 2008. The irrigation schedule was the same as that usedby most farmers in the region (Table 4), and the irrigation depth
was 150 mm. In order to leach accumulated salts, 200 mm of water
was additionally applied twice each year, before seeding in spring
and after harvest in winter.
2.4.3. Measurements
The soil water content was measured every 5 days using a neutron probe (CNC503DR) with three replications. The depths of thesoil moisture measurements were 10, 20, 30, 40, 50, 70, 90, 110,
130, and 150 cm. At the same time, in order to calibrate the neutron probe, soil samples were taken each month to measure thegravimetric soil water content in the laboratory. The position of
the groundwater table was measured every 5 days using a con-Table 3
Phenological phases of cotton growth in 2007 and 2008.
Year

Seeding

Emergence stage

Squaring stage

Flowering stage

Topping stage

Boll opening stage

Harvest2007
2008

04/26/2007
04/26/2008

05/07/2007
05/07/2008

06/18/2007
06/08/2008

07/19/2007
06/28/2008

07/23/2007
07/14/2008

09/05/2007
08/25/2008

11/05/2007
11/01/2008

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

69

Fig. 2. Irrigation and precipitation rates (A), and volumetric soil water contents at a depth of 10 cm and groundwater depths (B) in 2007 and 2008.

ductance method in two observation wells located inside of theexperiment ﬁeld. A fully equipped meteorological station, located
50 m from the experiment site, was used to collect meteorological
data such as solar radiation, air temperature, rainfall, wind speed,
and relative humidity. Temporal changes in soil water contents andgroundwater depths during 2007 and 2008, together with irrigationand precipitation amounts, are shown in Fig. 2.
The leaf area index (LAI), crop height, and aboveground biomass
were measured by destructive sampling with three replications ineach crop phenological phase. The root samples were taken with
six replications at the end of the growing season each year. The root
auger (125.6 cm3 ) was used to take root samples at depths of 0–10,
10–20, 20–30, 30–40, 40–60, 60–80, and 80–100 cm. The roots inthe auger were carefully cleaned and then the dry weight of roots ineach sample depth was measured. The root density was deﬁned as
the total root dry weight divided by the sampling volume in eachsampling depth. The measured root densities in years 2007 and2008 are shown in Fig. 3. The cotton yield was measured before
harvest by destructive sampling with three replications.
Additional undisturbed soil samples (100 cm3 ) were collected at
depths of 0–10, 10–20, 20–30, 30–40, 40–60, 60–80, and 80–270 cm
to measure the saturated hydraulic conductivity, the saturated
water content, and the dry bulk density. Values of the bulk density, saturated water content, and saturated hydraulic conductivity
are presented in Table 5.
2.4.4. Numerical modeling
The soil water dynamics in the experiment ﬁeld were evaluatedas a one-dimensional problem. The depth of the simulation proﬁle was 2.7 m, which was discretized into a grid size of 0.01 m. The

Fig. 3. Observed root densities in 2007 and 2008. Error bars represent standard
deviations of the measurements.

time-variable pressure head boundary condition (BC) was speciﬁed at the bottom boundary using the measured groundwater table
depths. The atmospheric boundary condition was used at the upperboundary. The potential transpiration and evaporation, required as
input by the atmospheric BC, were obtained each day from the crop
growth model. Irrigation was considered using the atmospheric BC
without surface runoff, since surface runoff was not observed during the ﬁeld experiment. Initial conditions were set according to soil

70

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

Table 5
Soil hydraulic parameters of the van Genuchten–Mualem model for different soil
horizons at the experiment site. Parameters Ks and  s are measured, while parameters  r , ˛, and n are estimated using inverse simulations.
Depth (cm) r (−)

 s (−)

␣ (m−1 )n

Ks (m day−1 )

0–10
10–20
20–30
30–40
40–60
60–80
80–270

0.04
0.04
0.04
0.04
0.04
0.04
0.04

0.45
0.48
0.48
0.45
0.48
0.47
0.42

0.92
0.65
0.90
0.37
0.37
0.50
0.23

1.30
1.36
1.30
1.37
1.44
1.31
1.17

0.27
0.13
0.09
0.11
0.09
0.05
0.05

Table 6
Optimized parameters for the crop growth model.
Parameter

Optimized Value

Units

RUE
HIopt
LAImx
dlai
frgrw1
laimx1
frgrw2laimx2
hc,mx
Tbase
HImin
kl
PHU

25.0
0.60
4.00
0.85
0.28
0.11
0.38
0.48
0.65
25.0
0.33
0.1
2400

(kg/ha)/(J/m2 )(kg/ha)/(kg/ha)
m2 /m2
–
–
–
–
–
M
◦
C
–
–
–

water content measurements. The time period from April 7, 2007 to
December 5, 2007 was used to calibrate the coupled model, while
the time period from April 5, 2008 to December 1, 2008 was used
to validate the model.
Seven soil horizons were deﬁned in the soil proﬁle, reﬂecting
the sampling protocol (0–10, 10–20, 20–30, 30–40, 40–60, 60–80,
and 80–270 cm) (Table 2). Measured values of the bulk density,
saturated soil water content, and saturated hydraulic conductivity(Table 5) were used as model inputs and remained constant during
the model calibration and validation. Pore-connectivity parameter (l) was kept ﬁxed and equal to 0.5. Additional soil hydraulicparameters ( r , ˛, and n) were optimized using the model independent Parameter Estimation Tool (PEST, Watermark Computing,2005) and observed volumetric soil water contents. The minimized
objective function was deﬁned as the sum of squared deviations
between observed and simulated values, while the same weights
were assigned to all data points. The optimized parameters are
given in Table 5.
The measurements of LAI, cotton height, above ground biomass
and yield in 2007 and 2008 were used to calibrate and validate theparameters of the crop growth model. The calibrated and validated
parameters of the crop growth model are shown in Table 6.
2.5. Criteria of model evaluation
Agreement between simulated results and observed data wasevaluated using the correlation coefﬁcient R2 (−) and the root mean
standard error (RMSE) for each treatment:


RMSE =

1
(Csi − Cob )2
n
n

1/2

R2 =

n i=1i=1¯ )(Csi − C¯si )(Cob − Cobn(Csi − C¯si )i=1

2.6. Model sensitivity experiments
Additional simulation experiments were carried out to examine the effects of the groundwater depth on plant physiological
functions and the root zone water balance. The effects of thegroundwater table depth on the maximum LAI, cotton yield, above
ground biomass, potential transpiration, actual transpiration, andcapillary rise (deﬁned as ﬂow from the bottom of the soil proﬁleinto the root zone) were evaluated. In these additional simulations,the root zone was considered to extend from the soil surface to adepth of 60 cm, while additional input parameters from 2007 wereused. Different groundwater table depths were considered in 16
different modeling scenarios. For each scenario, groundwater table
time series were generated by adding a constant value (ranging
from −2.5 to 2.0 m) to the groundwater table depth measured in2007. One additional simulation was carried out, in which a freedrainage BC was considered at the bottom boundary. In this ﬁnal
control simulation, there is no groundwater present in the soil proﬁle and it thus can have no effect on processes such as the root zone
water balance and/or crop growth.
3. Results and discussion
3.1. Model calibration and validation
3.1.1. Soil water content
The correspondence between measured and simulated volumetric soil water contents during both calibration and validation
years is shown in the scatter plot displayed in Fig. 4. Good agreement between simulated and measured volumetric soil watercontents was found for both calibration (Fig. 4A) and validation(Fig. 4B) years. The model produced slightly lower values of watercontents than those observed during the calibration process, while
the opposite occurred during the validation process. Various statistical tests were carried out to investigate the performance of thecoupled model. The R2 and RMSE values for soil water contents arepresented in Table 7. As shown in Table 7, the R2 values were in therange of 0.70–0.81 during both calibration and validation periods.
The RMSE values were 0.032–0.027 for the calibration and validation periods, respectively. These statistical measures indicate a highconsistency between simulated and measured values during both
calibration and validation periods.
Simulated and measured soil water contents at depths of 20 and
150 cm are presented in Fig. 5, which shows that simulated valuesof the volumetric soil water content are in close agreement withobserved values during both calibration (Fig. 5A) and validation(Fig. 5B) periods. Furthermore, the pattern of ﬂuctuating soil watercontents consistently reﬂected the irrigation events, especially inthe shallow soil depth of 20 cm. Apparently, the soil water content
in the upper soil layers produced more dramatic changes than in
the deeper soil layers. Based on Table 7 and Figs. 4 and 5, it canbe concluded that the coupled model performed well in simulating
volumetric soil water contents in the root zone.

(24)

i=1

n

¯ are the mean values of simulated andtion process, and C¯si and Cob
observed data points, respectively. The units of the RMSE are thesame as those of a compared variable.

¯ )(Cob − Cob

2
(25)

where Csi is a simulated value, Cob is an observed value, n is thetotal number of observed values used in the calibration and valida-

3.1.2. Cotton growth
The measured and simulated LAIs, above ground biomass, andcotton heights are shown in Fig. 6. Good agreement was again foundbetween measured and simulated values during both calibration(Fig. 6A) and validation (Fig. 6B) periods. The dynamics of LAI, above
ground biomass, and the cotton height during the cotton growthseason were captured well by the coupled model during both cali-M. Han et al. / Agricultural Water Management 160 (2015) 64–75

71

Fig. 4. A comparison of simulated and observed volumetric soil water contents during the calibration (2007; A) and validation (2008; B) periods.

Table 7
Statistical tests for modeling results.
PeriodYear

Item

LAI (−)

Aboveground biomass (kg/ha)

Cotton height (cm)Soil water content (−)

Calibration process

2007

Validation process

2008

R2
RMSE
R2
RMSE

0.91
0.54
0.91
0.62

0.97
1019.80.95
1761.5

0.94
8.48
0.96
4.38

0.81
0.032
0.70
0.027

Table 8Simulated and observed cotton yields.
Year

Observed cotton yield (kg/ha)

Simulated cotton yield (kg/ha)

2007
2008

5060.0 ± 694.8
5063.0 ± 373.9

5022.95034.1bration and validation periods. Only the simulated crop height wasslightly overestimated in 2007 and underestimated in 2006, likely
due to inadequate accounting for the water stress in Eq. (15). The
statistical tests, carried out to evaluate the performance of the coupled model, such as R2 and RMSE for LAI, above ground biomass, andcotton height are presented in Table 7. The R2 values for LAI, above
ground biomass, and cotton height were all higher than 0.9. The
RMSE values for LAI, above ground biomass, and cotton height were
all close to their corresponding maximum standard measurement
errors, which were 0.7, 2405 kg/ha, 5.1 cm, respectively.
The model performance was also evaluated using the measured
and simulated cotton yield. As shown in Table 8, the simulated yield
was similar to the observed yield for both calibration and validation
periods. Standard deviations of the observed yields were larger thanthe difference between the observed and simulated yields. Basedon Tables 7 and 8 and Fig. 6, it can be concluded that the coupledmodel performed well in simulating cotton growth.
3.2. Model sensitivity experiments
3.2.1. Impact of the groundwater table on the root zone waterbalance
The impacts of the groundwater table level on total actual
transpiration, total actual evaporation, total root zone soil watercontent, and capillary rise from groundwater are shown in Fig. 7.

The total root zone soil water content was deﬁned here as thewater content in the top 60 cm of the soil proﬁle and the groundwater table change was the difference between the groundwatertable depth observed in 2007 and groundwater table depths considered in the sensitivity simulations (see also Section 2.5). As canbe observed, the groundwater table level has a signiﬁcant impact on
the root zone water balance. As the groundwater table approaches
the soil surface, the response of various components of the rootzone water balance can be divided into two parts. As the change
in the groundwater table level increased from −2.27 to 0.72 m, thetotal transpiration of cotton showed an increasing trend (Fig. 7A).
This occurred because as the groundwater table level increases, so
increases the capillary rise from groundwater into the root zone
(Fig. 7B). An increase in capillary rise subsequently leads to anincrease in the total soil water content in the root zone and cotton
transpiration (Fig. 7B). Similarly, the average water stress index ofcotton (1 – actual transpiration/potential transpiration) decreased
to values close to 0, indicating that increasing the groundwatertable reduces the cotton water stress and plays a positive role incotton growth.
Note that the reversal point is slightly different for the water
stress index and actual transpiration. The reversal point would
be the same if potential transpiration was the same in all simulations. However, due to different crop growth and different LAIs
for different positions of the groundwater table, different potential
transpirations were obtained for different simulations.
The response of various components of the root zone water balance was different when the groundwater table change further
increased between 0.72 and 1.72 m. The total cotton transpiration decreased and the cotton water stress increased (Fig. 7A).
This occurred because as the total root zone water content further increased, it induced anaerobic conditions, which produced an

72

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

Fig. 5. Simulated (lines) and measured (dots) volumetric soil water contents at depths of 20 (top) and 150 (bottom) cm during the calibration (2007; left) and validation(2008; right) periods. Error bars represent standard deviations of the measurements.

oxygen stress on crop roots. Since cotton transpiration decreased,
the capillary rise from groundwater showed a decreasing trend forthese groundwater table levels. A further increase of the groundwater table level thus has a negative impact on cotton growth. On theother hand, soil evaporation further increased because more water
was available in the root zone with the increasing groundwatertable level.
3.2.2. Impact of the groundwater table on cotton growth
The impact of the groundwater table on the cotton maximum
LAI and yield is shown in Fig. 8. Similarly, as the response of various components of the root zone water balance, crop transpirationin particular, changes in LAI and cotton yield can also be divided
into two parts. First, as the groundwater table level increases (from−2.27 to 0.52 m), the maximum cotton LAI and yield are bothincreasing. This is due to an increase in the total soil water content
in the root zone and a decrease in the cotton water stress. Second, as
the groundwater table level increases further (from 0.52 to 1.72 m),
the maximum cotton LAI and yield both decrease. This is caused byanaerobic conditions in the root zone.
Evaluated scenarios indicate that the optimal change in the
groundwater table level is about 0.52 m. For these conditions, theaverage groundwater table depth during the cotton growth season
is 1.84 m, total transpiration is about 514 mm, the maximum LAI is3.53, and the cotton yield is 5118 kg/ha. Compared to the control

simulation without groundwater (a scenario with the free drainage
boundary condition), total transpiration increased by 117.6 mm,
which represents an increase of 23%. This indicates that when the
groundwater table is raised during the season by 0.52 m, 23% ofcrop transpiration is supplied by the capillary rise from groundwater. The cotton yield increased by 1067 kg/ha, which represents
an increase of 20%. This means that groundwater is an importantsource of water in this region, which has a signiﬁcant impact on thesoil water balance and cotton growth.
The sensitivity analysis also revealed the interactions between
cotton growth and hydrological processes. The impact of cotton growth on hydrological processes is more signiﬁcant when
the groundwater table is shallow. A decrease in cotton actual
transpiration reduces the capillary rise from groundwater. Whenthe groundwater table is at an optimal depth, an increase incotton LAI produces an increase in the cotton potential transpiration, which is another driver for a capillary rise fromgroundwater. An increase in the groundwater table can signiﬁcantly increase cotton transpiration, reduce the cotton water stress,
and increase cotton yield. However, when the groundwater table
is too close to the soil surface, the oxygen stress on cotton roots
increases, inducing anaerobic conditions in the root zone and dramatically reducing the cotton yield. It is important to carefully
control the shallow groundwater table when groundwater is to be
used as a water source for crops.

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

73

Fig. 6. Simulated (lines) and measured (dots) leaf area index (top), cotton height (middle), and aboveground biomass (bottom) during the calibration (2007; left) andvalidation (2008; right) periods. Error bars represent standard deviations of the measurements.

Fig. 7. An impact of the groundwater table depth change on components of the root water balance obtained from simulated results. The groundwater table change wascompared with groundwater table depth measurements in 2007.

74

M. Han et al. / Agricultural Water Management 160 (2015) 64–75

Acknowledgements
The study was supported by the Key National Natural Science Foundation (U14003281), (4117095), National Project (No.
2013BAC10B01) of China and Project SuMaRiO of BMBF (German
Federal Ministry of Science and Technology) – Sustainable Management of River Oases along the Tarim River. The authors are gratefulto Dr. Markus Disse and Keilholz Patrick of Technical University
Munich, Germany, for their help.

References

Fig. 8. An impact of the groundwater table depth change on simulated maximum
LAI and cotton yield.

4. Conclusions
Groundwater is an important factor that needs to be considered when designing sustainable land management for arid oases.
However, its impact on plant functioning and the root zone waterbalance has not been sufﬁciently examined, partly because numerical models have not fully considered all the interactions between
the soil water balance and plant growth. In this study, simulations
were performed to examine the inﬂuence of groundwater on cotton
growth and the root zone water balance using a numerical model,
which coupled Hydrus-1D with a simpliﬁed crop growth modelfrom SWAT.
The simulation results of the coupled model were comparedwith experimental data obtained from cotton ﬁeld experiments.
Results suggest that volumetric soil water contents, LAIs, above
ground biomass and cotton yields simulated by the coupled model
were in good agreement with the measurements.
Additional model simulations showed that groundwater is amajor source of water for cotton growth. Compared to a control
simulation that had no groundwater, 23% of crop transpiration
is supplied by a capillary rise from groundwater, producing an
increase in cotton yield by 20%.
The model simulations also showed that the cotton growth androot zone water balance are very sensitive to the depth of thegroundwater table and that cotton growth, in turn, affects subsurface water ﬂuxes, such as capillary rise. When the groundwatertable depth is raised from its positions in 2007 by less than 0.52 m, it
has a positive effect on cotton growth by making more water available in the root zone and reducing cotton water stress and thereby
enhancing the maximum cotton LAI and yield in this region. However, if the groundwater table is raised by more than 0.52 m, it thenhas a negative effect on cotton growth, mainly by creating anaerobic conditions in the root zone. Similarly, when the groundwatertable has a positive effect on cotton growth, cotton growth wouldin turn increase capillary rise from groundwater by increasing LAI,
and potential transpiration. On the other hand, when the groundwater table is raised more than 0.52 m, decreasing cotton LAI andactual transpiration decrease capillary rise from groundwater.
It can be concluded that groundwater is a crucial factor thatneeds to be taken into consideration when evaluating agriculturalland management in this arid region. Our work, presented in thismanuscript, provides a useful modeling tool for evaluating local
land management and for designing sustainable conditions of arid
oases.

Allen, R.G., Pereira, L.S., Raes, D., Smith, M., 1998. Cropevapotranspiration—guidelines for computing crop water requirements. FAO
Irrigation and Drainage Paper 56.
Ayars, J.E., Christen, E.W., Soppe, R.W., Meyer, W.S., 2006. The resource potential of
in-situ shallow ground water use in irrigated agriculture: a review. Irrigation
Sci. 24, 147–160.
Computing, W., 2005. Model-Independent Parameter Estimation User Manual, 5thed. Watermark Computing, Brisbane Australia.
Feddes, R.A., Kowalik, P.J., Zaradny, H., 1978. Simulation of Field Water Use andCrop Yield. John Wiley & Sons, New York.
Forkutsa, I., Sommer, R., Shirokova, Y.I., Lamers, J.P.A., Kienzler, K., Tischbein, B.,
Martius, C., Vlek, P.L.G., 2009. Modeling irrigated cotton with shallow
groundwater in the Aral Sea Basin of Uzbekistan: II. Soil salinity dynamics.
Irrigation Sci. 27, 319–330.
Harbaugh, A.W., Banta, E.R., Hill, M.C., McDonald, M.G., 2000. MODFLOW-2000, theUs Geological Survey Modular Ground-water Model: User Guide to
Modularization Concepts and the Ground-water Flow Process. US Geological
Survey Reston, VA, USA.
Hou, Z., Chen, W., Li, X., Xiu, L., Wu, L., 2009. Effects of salinity and fertigation
practice on cotton yield and 15 N recovery. Agric. Water Manage. 96,
1483–1489.
Hu, X.T., Chen, H., Wang, J., Meng, X.B., Chen, F.H., 2009. Effects of soil watercontent on cotton root growth and distribution under mulched drip irrigation.
Agric. Sci. China 8, 709–716.
Jones, J.W., Hoogenboom, G., Porter, C.H., Boote, K.J., Batchelor, W.D., Hunt, L.,
Wilkens, P.W., Singh, U., Gijsman, A.J., Ritchie, J.T., 2003. The DSSAT cropping
system model. Eur. J. Agron. 18, 235–265.
Jorenush, M.H., Sepaskhah, A.R., 2003. Modelling capillary rise and soil salinity forshallow saline water table under irrigated and non-irrigated conditions. Agric.
Water Manage. 61, 125–141.
Kang, Y., Wang, R., Wan, S., Hu, W., Jiang, S., Liu, S., 2012. Effects of different waterlevels on cotton growth and water use through drip irrigation in an arid region
with saline ground water of Northwest China. Agric. Water Manage. 109,
117–126.
Kroes, J., Van Dam, J., Groenendijk, P., Hendriks, R., Jacobs, C., 2008. SWAP Version
3.2: Theory Description and User Manual. Alterra Wageningen, The
Netherlands.
Liu, M., Yang, J., Li, X., Yu, M., Wang, J., 2012. Effects of irrigation water quality anddrip tape arrangement on soil salinity, soil moisture distribution, and cottonyield (Gossypium hirsutum L.) under mulched drip irrigation in Xinjiang, China.
J. Integr. Agric. 11, 502–511.
Neitsch, S., Arnold, J., Kiniry, J., Williams, J., King, K., 2005. SWAT: Soil and Water
Assessment Tool: Theoretical documentation, version 2005. Texas, USA.
Richards, L.A., 1931. Capillary conduction of liquids through porous mediums.
Physics 1, 318–333.
Shen, X., Zhang, J., Liy, Z., Li, M., Yang, G., Sun, J., 2012. Effects of different water
treatments on yield and water use efﬁciency of cotton with drip irrigation
under ﬁlm mulch. Agric. Res. Arid Areas 30, 118–124.
Šimůnek, J., van Genuchten, M.T., Sejna, M., 2005. The HYDRUS-1D software
package for simulating the one-dimensional movement of water, heat, andmultiple solutes in variably-saturated media. In: Version 3.0, HYDRUS
Software Series 1. Department of Environmental Sciences, University ofCalifornia Riverside, Riverside, CA, USA, pp. 1–270.
Šimůnek, J., van Genuchten, M.T., Šejna, M., 2008. Development and applications of
the HYDRUS and STANMOD software packages and related codes. Vadose Zone
J. 7, 587–600.
Šimůnek, J., Hopmans, J.W., 2009. Modeling compensated root water and nutrient
uptake. Ecol. Model. 220, 505–521.
Soppe, R.W.O., Ayars, J.E., 2003. Characterizing ground water use by safﬂower
using weighing lysimeters. Agric. Water Manage. 60, 59–71.
Soylu, M.E., Kucharik, C.J., Loheide II, S.P., 2014. Inﬂuence of groundwater on plantwater use and productivity: development of an integrated ecosystem –variably saturated soil water ﬂow model. Agric. For. Meteorol. 189, 198–210.
Twarakavi, N.K.C., Šimůnek, J., Seo, S., 2008. Evaluating interactions between
groundwater and vadose zone using the HYDRUS-based ﬂow package for
MODFLOW. Vadose Zone J. 7, 757–768.
van Genuchten, M.T., 1980. A closed-form equation for predicting the hydraulicconductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44, 892–898.
Wang, Z., Jin, M., Šimůnek, J., van Genuchten, M., 2014. Evaluation of mulched drip
irrigation for cotton in arid Northwest China. Irrigation Sci. 32, 15–27.

M. Han et al. / Agricultural Water Management 160 (2015) 64–75
Williams, J., Wang, E., Meinardus, A., Harman, W., Siemers, M., Atwood, J.D., 2006.
EPIC users guide v. 0509. Blackland Research and Extension Center, Temple,
Texas.
Williams, J.R., Jones, C.A., Kiniry, J.R., Spanel, D.A., 1989. The epic crop
growth-model. Trans. Asae 32, 497–511.
Williams, J.R., Singh, V., 1995. The EPIC model. Computer Models of Watershed
Hydrology, 909-1000.

75

Yan, Y.Y., Zhao, C.Y., Sheng, Y., Li, J.Y., Peng, D.M., Li, Z.L., Feng, S.L., 2009. Effects ofdrip irrigation under mulching on cotton root and shoot biomass and yield. J.
Appl. Ecol. 20, 970–976.
Zaidi, P.H., Raﬁque, S., Rai, P.K., Singh, N.N., Srinivasan, G., 2004. Tolerance toexcess moisture in maize (Zea mays L.): susceptible crop stages andidentiﬁcation of tolerant genotypes. Field Crops Res. 90, 189–202.