A distributed hydrological model, the Grid-based Integrated Surface-groundwater MODel (GISMOD), was developed to simulate hydrological processes by considering water interaction among different soil layers. The model integrates six modules. Basic information on catchment, such as the flow direction and the drainage network, can be obtained automatically from digital elevation model (DEM) data by using the preprocessing module. In GISMOD, three methods are available to estimate precipitation and eight to estimate evapotranspiration in each grid. Infiltration excess flow, saturation flow, recharge flow from soil and groundwater are all considered by using a simplified method according to the three-layer structure of GISMOD. A case study in the upper-middle reaches of the Heihe River was presented to evaluate the model performance. The results show that the change tendencies of infiltrated water and recharged water are opposite in the study area; the surface water is mainly infiltrated to the soil layer in the upper streams, and then discharged by groundwater in the middle reaches, which is more consistent with the actual situation in the Heihe River basin. In addition, the simulated runoff of a river grid near the Yingluoxia hydrological station is compared with the observed one and the results also demonstrate that GISMOD has a better performance in runoff estimation on both daily and monthly scales.

INTRODUCTION

The dynamic processes of surface water and groundwater interaction are mainly affected by geological/meteorological characteristics and represent the natural features and hierarchical structure of the watershed (Imagawa et al. 2011). These two kinds of process are not isolated components of the hydrological system in the environment, and the interactions between them are complex and usually difficult to estimate due to the problems of heterogeneity and scale (Halford & Mayer 2000; Kalbus et al. 2006). Therefore, a better understanding of the basic principles of interactions between groundwater and surface water is important for efficient management of water resources (Gardner 1999; Sophocleous 2002). The hydrological model is an important scientific tool, which has focused on the simulation of rainfall-runoff relationships according to a conceptual or physical representation method. Models can be used to reflect groundwater and surface water interaction as well as runoff generation for water resources evaluation, management and utilization on regional or watershed scales (Krause et al. 2009; Khakbaz et al. 2011).

A number of hydrological models have been utilized to quantify stream–aquifer flow in recent years. For instance, Kim et al. (2008) developed the SWAT-MODFLOW model and applied it to estimate groundwater and streamflow reduction by pumping groundwater out. A coupled surface water and groundwater flow model (MODBRANCH) was developed by Swain & Wexler (1996), and it integrated MODFLOW and BRANCH to simulate unsteady flow dynamics in an open channel network. MIKE-SHE and MOD-HMS have been developed to describe the overall exchange and interaction by coupling surface water and groundwater flows (Abbott et al. 1986; Refsgaard & Storm 1995; Graham 2001; Panday & Huyakorn 2004; Graham & Butts 2005). Some of these models were developed to investigate the variations of groundwater level and surface water level for special purposes. For example, Yuan et al. (2011) adopted a coupled model that simulates the interaction between surface water and groundwater flow, and solute transport processes in the coastal wetland. Zhang & Ross (2007) developed a two-layer vadose zone model for surface and groundwater interactions. Imagawa et al. (2011) improved the hydro-environmental watershed model in a water exchange process through surface water-groundwater interaction, and applied it to an agricultural watershed.

However, all of these models have different strengths and limitations. As a matter of fact, hydrological processes are far from simple: watershed is highly heterogeneous in terms of land use, soil characteristics and topography (Pilgrim et al. 1988; Thornes et al. 1996; Kosmas et al. 1997). Conceptually, complex models should perform well by numerical simulation with high quality input data. However, several studies show that simple models can perform as well as complex ones (Michaud & Sorooshian 1994; Refsgaard & Knudsen 1996). Therefore, it is always difficult to select between a simple model and a complex one for hydrological simulation.

In order to describe the interaction between surface water and groundwater in a simple and efficient way, a new kind of hydrological model (a Grid-based Integrated Surface water–groundwater interaction MODel, simplified GISMOD) was developed to estimate the exchange water amount based on a three-layer tank method. The principles of GISMOD are presented first and then a case study in the Heihe River basin is presented to evaluate the performance of this model.

MODEL DESCRIPTION

GISMOD is a physically based distributed hydrological model that was developed to simulate hydrologically relevant processes in river basins, such as water movement, infiltration, interaction and evapotranspiration, by a simplified tank method (Xu et al. 2001).

A grid square method is used for spatial disaggregation in GISMOD, which divides the watershed into hundreds or thousands of grids of equal size, assuming that each grid is homogeneous in structure (Figure 1). In GISMOD, the total number of grids is dependent on the resolution of input data (digital elevation model (DEM), land cover, soil type, etc.) and these input data are needed to keep in line with each other. That is to say that the precision of the simulation is highly dependent upon the quantity and quality of the input data. As with SWIM (Soil and Water Integrated Model) (Krysanova et al. 1998), the soil column is subdivided into three layers (surface layer, soil layer and groundwater layer) in this model to simulate surface flow, interflow and groundwater flow.
Figure 1

Schematic map of the GISMOD.

Figure 1

Schematic map of the GISMOD.

GISMOD is composed of the following six components: a preprocessing module, an interpolation module, an evapotranspiration module, a runoff simulation module, a parameter identification module and a simulation setting module (Figure 2). There are six steps to run the model: (1) DEM data are required at first to extract river-net and simulation order by using the preprocessing module; (2) put the observed precipitation data and meteorological data (such as air temperature, relative humidity, average wind speed) into GISMOD and then spatially resolve into grid cells with the resolution of DEM; (3) calculate potential evapotranspiration of each grid by using the evapotranspiration module; (4) load three kinds of spatial data (land use map, soil map and geology map) and identify the initial value of the model parameters by using the parameter identification module; (5) select the simulation option and set the output grid; and finally (6) run the model and obtain the simulation result for comparison.
Figure 2

Modeling flow chart of the GISMOD.

Figure 2

Modeling flow chart of the GISMOD.

Preprocessing program

Similar to the hydrological tool in ArcGIS, the preprocessing program is used to delineate the river basin and to extract the drainage networks automatically from the DEM data (Li et al. 2013). This program has various functions such as flow direction definition, flow accumulation calculation, and drainage network generation, and the order of runoff simulation can be analyzed conveniently by using a special approach (Xu 2009).

Spatial interpolation

The accuracy of the simulation result depends primarily on how well the variability of the hydro-meteorological information can be specified in a certain area. Thus, the disaggregation of the available hydro-meteorological data from the limited observation station to the entire catchment is essential to estimate the evapotranspiration and precipitation at each grid cell in the GISMOD. Considering the different architectures of the observation networks and climatic characteristics, the Thiessen method, gradient plus inverse distance squared method and inverse distance squared method are developed in the model for spatial interpolation (Nalder & Wein 1998; Teegavarapu & Chandramouli 2005).

Evapotranspiration

In order to measure the ability of water to move from the surface to the atmosphere, and to measure and the processes of evaporation and transpiration without control of water supply, potential evapotranspiration (PET), as an important component of the hydrological cycle and the key variable of the hydrological models, should be estimated as accurately as possible. Calculation methods of PET can be classified into four types (temperature-based, radiation-based, mass-transfer and combination). Most of them need many kinds of input data such as air temperature, air pressure, humidity, wind speed, solar radiation, etc. In order to satisfy the needs of the different available data, GISMOD has integrated eight methods (Turc 1961; Monteith 1965; Priestley & Taylor 1972; Hargreaves & Samani 1985; Allen et al. 1998) to estimate PET (Table 1) and also enable users to import their own PET data.

Table 1

Eight methods for estimating PET in GISMOD

Num Name Equation 
FAO Penman–Monteith  
Penman  
Kimberly–Penman  
Hargreaves–Samani  
Priestley–Taylor  
Makkink  
Turc  
Doorenboss–Pruitt  
Num Name Equation 
FAO Penman–Monteith  
Penman  
Kimberly–Penman  
Hargreaves–Samani  
Priestley–Taylor  
Makkink  
Turc  
Doorenboss–Pruitt  

When ET0 is the potential evapotranspiration, Rn is net radiation flux, G is sensible heat flux into the soil, which can be ignored for daily estimation; T is the temperature, usually taken as the daily mean air temperature (Tmean), Tmax and Tmin is the daily maximum and minimum air temperatures, respectively; es is the vapor pressure of the air at saturation, ea is the actual vapor pressure, Δ is the slope of the saturation vapor pressure temperature relationship curve and γ is the psychometric constant, U2 is the wind speed at 2 m height, λ is the latent heat of vaporization, Rs is the solar or shortwave radiation, and Ra is the extraterrestrial radiation.

Actual evapotranspiration (AET) can be assumed which depletes water from the root zone of the soil (Beven et al. 1995; Kennen et al. 2008), and which is expressed as a simple function of PET and the water content of the surface layer (Xu 2009): 
formula
1
where h1 is the water level of the surface layer and Sf2 is the height of the surface layer.

Interception

The rainfall intercepted by vegetation is calculated according to the relationship of the leaf index (LAI), precipitation and vegetation coverage: 
formula
2
where Cint is vegetation coverage, LAImax is the maximum leaf index and P is precipitation.
LAI is computed by using an equation of effective accumulated temperature (Hu & Zhao 2008): 
formula
3
where δ and ε are the vegetation parameters ranging from 0–8 and 0.5–1, respectively, which are decided based on the vegetation type and soil moisture.
The difference vegetation index (DVI) is calculated by the following formula: 
formula
4
where AT is the total amount of effective accumulated temperature, B is the threshold of daily effective accumulated temperature and Tday is the daily temperature.

Runoff module

Two types of grid are designed in the runoff module of GISMOD: one is the ordinary grid which is vertically divided into three layers: the surface layer, the soil layer and the groundwater layer; the other is the river grid, which looks like a single reservoir for channel routing (Figure 1). In an ordinary grid, each portion of the outflow would pour into its corresponding parts when the downstream grid is also an ordinary grid, otherwise, yielding water will come together as a stream flow into the river grid, then move among these river grids and eventually drain out of the watershed. Figure 3 is a schematic drawing that illustrates the process of water transfer and exchange in the runoff module. It is assumed that evapotranspiration first takes place in the surface layer. The amount of water used for evapotranspiration is mainly dependent upon the water content of the surface layer. If there is no water in the surface layer, water from the soil layer will be supplied to the surface layer.
Figure 3

Runoff generation structure of the GISMOD.

Figure 3

Runoff generation structure of the GISMOD.

Surface layer

The surface layer is a thin top soil layer covered by vegetation and regarded as a kind of reservoir with an opening bottom in GISMOD. Rainwater, the discharge from upstream grids and the water supplied by the soil layer are regarded as the inflow to the surface layer. Apart from evapotranspiration and infiltration, water in the surface layer is generated by two kinds of overland flow (saturation excess overland flow and infiltration excess overland flow), and by the lateral and infiltration flow.

The continuity water balance equation in the surface layer is expressed as: 
formula
5
where h1 is the water level of the surface layer, t is the time step, R is the daily precipitation, Qr0 is the inflow from the upstream neighboring surface layer, Qb0 is the recharge flow from the soil layer, E is evapotranspiration, Qs1 is the saturation excess flow, Qs2 is the infiltration excess flow, Qc0 is the lateral outflow (surface flow) and Qx0 is the infiltration flow.

Checking for water balance in each layer is carried out independently in GISMOD. First of all, the amount of available water in the surface layer should be estimated, including outflow from upstream, rain water and existing water in the layer. Then, available outflow from the surface layer will be calculated by considering different situations. When the available water is not enough to meet the demand for release, all of the water from the surface layer should be pouring out on the basis of proportional distribution. Afterwards, the water level in the surface layer will be decreased to 0. Otherwise, each part of the outflow should be estimated by using Equations (6)–(9), and the water level in the surface layer will be re-computed according to the continuity equation of water balance. In order to take account of the recharge flow from the soil layer, the water level in the surface layer will be finally adjusted.

Saturation excess flow will be generated while the water level (h1) exceeds the height of the surface layer (sf2), which is estimated by an equation derived from Manning's formula: 
formula
6
where i is the slope of grid, L is the length of the grid and n is Manning's roughness coefficient.
In contrast to the saturation excess flow, the infiltration excess flow only occurs when rainfall (R) is greater than the infiltration capacity (f0) of the surface layer. 
formula
7
where f0 is the stable infiltration rate of the surface layer, A is the area of the grid and tha is the discharge coefficient of the surface layer.
Lateral flow and infiltration flow are decided by the water level (h1) of the surface layer, which should be greater than 0 for keeping these two kinds of flow. The equations are written as: 
formula
8
 
formula
9
where α is the storage coefficient of the surface layer, so that the surface flow decreases as the storage coefficient α increases.

Soil layer

A thin layer of soil over underlying rock (called the soil layer), which is used in place of the leaching layer and subsoil layer in GISMOD, is designed to simply describe the dynamic changes of soil moisture. The supplements from the groundwater layer and surface layer are largely determined by the water content of the soil layer. The transmission and storage capacities of this layer are decided according to the soil type. Similar to the surface layer, the water level of the soil layer will be adjusted after considering the recharge water from the groundwater layer. The water balance equation of the soil layer can be written as: 
formula
10
where h2 is the water level of the soil layer, t is the simulation time step, Qr1 is the inflow of the soil layer, Qb1 is the recharge flow from the groundwater layer, Qc1 is the outflow of the soil layer and Qx1 is the infiltration flow of the soil layer.

Recharge flow will be generated if the available water exceeds the capacity of the soil layer. Another possibility is that the moisture content of the surface layer falls to fill the needs of evapotranspiration and infiltration. Meanwhile, free water from the soil layer could be provided for the requirements of the surface layer. There are three possible scenarios that can be described by the following equations.

  1. Less water can be supplied from the soil layer during the periods when there is a water shortage in the surface layer: 
    formula
    11
where sf0 is the particular water level of the surface layer, which means that the surface layer needs a water supply from the soil layer if the water is lower than this value.

  • 2. The available water of the soil layer is sufficient to meet the needs of surface layer: 
    formula
    12
where s1 is the particular water level of the soil layer, which means that the soil layer can afford water to the surface layer once the water level is higher than this value.

  • 3. Additional water is unnecessary in the surface layer, while there is an abundance of water in the soil layer: 
    formula
    13

The only limiting condition of interflow and infiltration flow of the soil layer is the water content, which should be greater than 0 for those two kinds of flow. The equations are written as: 
formula
14
 
formula
15
where s2 is the height of the soil layer, kx is the horizontal hydraulic conductivity of the soil layer, and kz is the vertical hydraulic conductivity of the soil layer.

Groundwater layer

The groundwater layer is adapted for distinguishing unconfined and confined flow by using a special tank which has one inlet and two outlets. The groundwater flow is more stable and changes slowly; rainfall and human activities have little influence. As with a semi-closed or closed system, only unconfined water can flow into the river grid in GISMOD. Regarded as the base flow with little water exchange, therefore, confined flow is not considered in flow routing. The water balance equation of the groundwater layer is shown as: 
formula
16
where h3 is the water level of the groundwater layer, t is the simulation time step, Qr2 is the inflow of the groundwater layer, Qc2 is the unconfined groundwater flow, and Qc3 is the confined groundwater flow.

A designated water level is used to distinguish unconfined flow and confined flow in this model, below which only confined water can be generated from the groundwater layer. Specifically, it is worth noting that the flow pattern of this layer is not only controlled by the configuration of the water table but also by the distribution of hydraulic conductivity of the rocks.

Similar to the soil layer, three kinds of water supply from the groundwater layer are shown as follows.

  1. The water level of the soil layer is below the particular value (s0) and the water level of the groundwater layer is above the designated value (ss1): 
    formula
    17
where s0 is the particular value of the water level of the soil layer, which means that the soil layer needs water from the groundwater layer. If the water level is lower than this value, ss2 is the height of the groundwater layer, and ss1 is the designated water level of the groundwater layer.

  • 2. The amount of available water held by the groundwater layer can satisfy the requirements of the soil layer: 
    formula
    18
  • 3. The water level of the soil layer is above the particular value (s0) and the total water is simultaneously beyond the limit of water-holding capacity in the groundwater layer. 
    formula
    19
Referring to Darcy's law, the two kinds of outflow (unconfined and confined flow) from the groundwater layer are calculated by equations shown as: 
formula
20
 
formula
21
where Au is the flow coefficient of unconfined water, Ag is the flow coefficient of confined water.

River grid

The river grid is designed in the form of a general single layer tank to represent the stream channel in the model. Based on the flow direction determined by a preprocessing program, water from ordinary grids will transfer between each other and then flow into the river grid, and eventually run off to the outlet of the catchment. In order to resolve the spatially limited extent of the river handled in this grid-based model (Clark et al. 2011), a basic equation from Wang & Xu (2011) is used as follows: 
formula
22
where L is the channel length, B is the channel width, Qr3 is the inflow of the river grid, is the discharge of the river grid.
The channel flow is simulated by coupling with the Manning equation: 
formula
23
B is simulated as a function of the upper catchment area Aup: 
formula
24
where k and s are constants (k = 7.0, s = 0.5). For the segments where the above equation is not suitable, the value of B may be input manually.

Parameter setting

A number of parameters are used in GISMOD, which should be categorized into three groups according to the grid structure and relevance with the properties of the soil, the geology and the type of land cover.

There are six parameters responsible for surface layer simulation (Table 2), the values of which are evaluated according to the land use map. In order to utilize the land use data effectively, land use type is classified into six categories on the basis of the Technical Rules of Land Use Survey (TCADC 1984). Initial values of the parameters in different land use type are listed in Table 2.

Table 2

Parameters for surface layer in the GISMOD

  Land use type
 
Parameters Woodland Grassland Farmland Urban Water Unused land 
Sf2(m) 0.8 0.5 0.5 0.5 
Sf0(m) 0.2 0.2 0.1 0.01 0.05 0.1 
h1(m) 0.5 0.4 0.3 0.01 0.3 0.3 
f0(m/d) 0.004 0.002 0.001 0.0004 0.001 0.004 
tha 0.24 0.28 0.28 0.32 0.24 0.32 
n(m–1/3d–10.7 2.0 1.5 0.3 0.03 0.7 
  Land use type
 
Parameters Woodland Grassland Farmland Urban Water Unused land 
Sf2(m) 0.8 0.5 0.5 0.5 
Sf0(m) 0.2 0.2 0.1 0.01 0.05 0.1 
h1(m) 0.5 0.4 0.3 0.01 0.3 0.3 
f0(m/d) 0.004 0.002 0.001 0.0004 0.001 0.004 
tha 0.24 0.28 0.28 0.32 0.24 0.32 
n(m–1/3d–10.7 2.0 1.5 0.3 0.03 0.7 

Sf2, height of surface layer; Sf0, height of outlet; h1, initial water level; f0, stable infiltration rate; tha, discharge coefficient; n, Manning's roughness coefficient.

The parameters responsible for the soil layer are adjusted according to the land use map and generalized into three sorts by the soil infiltration capacity (Table 3). Initial values of parameters of the three soil types are given in Table 3.

Table 3

Parameters for soil layer in the GISMOD

  Soil infiltration capacity
 
Parameters Low Medium High 
S2(m) 30 40 50 
S1(m) 10 15 
S0(m) 10 
h2(m) 15 20 25 
Kz(m/d) 0.0001 0.0002 0.0005 
Kx(m/d) 0.0003 0.0015 0.008 
  Soil infiltration capacity
 
Parameters Low Medium High 
S2(m) 30 40 50 
S1(m) 10 15 
S0(m) 10 
h2(m) 15 20 25 
Kz(m/d) 0.0001 0.0002 0.0005 
Kx(m/d) 0.0003 0.0015 0.008 

S2, height of soil layer; S1, exchange water level; S0, height of outlet; h2, initial water level; Kx, horizontal hydraulic conductivity; Kz, vertical hydraulic conductivity.

As with the soil layer, the parameters of the groundwater layer are also classified into three categories depending on the different rock permeabilities in the model. Although there is no classification standard for rock permeability, the rock porosity can be assumed instead. In general, greater grain size and smaller porosity tend to decrease the penetrability, while a decrease in range of particle size and greater porosity tends to increase penetrability. Table 4 shows the reference value of parameters of the groundwater layer used in GISMOD.

Table 4

Parameters for groundwater layer in the GISMOD

  Rock penetrability
 
Parameters Low Medium High 
Ss2(m) 30 40 50 
Ss1(m) 10 15 
h3(m) 15 20 25 
Au (m–1/2d–1/20.0001 0.0002 0.0003 
Ag(1/d) 0.00002 0.00004 0.00008 
  Rock penetrability
 
Parameters Low Medium High 
Ss2(m) 30 40 50 
Ss1(m) 10 15 
h3(m) 15 20 25 
Au (m–1/2d–1/20.0001 0.0002 0.0003 
Ag(1/d) 0.00002 0.00004 0.00008 

Ss2, height of groundwater layer; Ss1, height of outlet; h3, initial water level; Au, flow coefficient of unconfined water; Ag, flow coefficient of confined water.

Simulation setting

While the output grid is determined, calibration and validation can be conducted on a daily or monthly time step according to the spatially interpolating hydro-meteorological data. The simulation result can be compared on the basis of model performance with regards to objective criteria representing the water balance.

Two objective functions are adopted in order to optimize the parameters, including relative error (RE) and Nash–Sutcliffe efficiency (NSE) (Nash & Suthcliff 1970). RE represents systematic water balance error, and the positive values mean overestimation while negative values mean underestimation. NSE measures the fraction of the variance of observed values explained by the model, which can vary from − ∞ to 1 and should be as near as possible to 1. RE and NSE are estimated as follows: 
formula
25
 
formula
26
where and are mean values of the observed and simulated discharge, and are the observed and simulated discharge series.

CASE STUDY

Study area description

The Heihe River is the second largest inland river in Northwest China. The river originates from the north foot of Qilian Mountain and passes nearly 821 km north through Qinghai and Gansu provinces, and finally reaches Lake Juyan in the Inner Mongolia Autonomous Region (Tang et al. 1992; Lu et al. 2003). Yingluoxia (Yingluoxia Gorge) and Zhengyixia (Zhengyixia Gorge) divide the drainage area of the Heihe River into upper, middle and lower reaches (Qin et al. 2011). The climate in the upper reaches is damp and cold; 70% of the area is mountainous. The average annual temperature and precipitation in the upper reaches are about 2 °C and 350 mm, respectively. In the middle reaches, annual precipitation is only 140 mm because of the dry weather, and the average annual temperature is between 14 and 18 °C.

As a particular hydro-geological characteristic, the transformation processes between the surface water and the groundwater are frequent and complex in the upper and middle reaches of the Heihe River. The local hydrological cycle starts with surface water going underground at the foot of the Qilian mountains, then re-emerging at the lower portion of the plain to form surface water again (Qiu 1987). The surface water–groundwater–surface water cycle is due to the different geological structures and often repeats two or three times in the alluvial plain (Wang & Cheng 2000). In order to investigate the characteristics and trends of spatial alternation of water exchange, the upper-middle reaches of the Heihe River basin are selected as the study area with an area of 72,200 km2 (Figure 4).
Figure 4

Location of the study area and hydrological and meteorological stations.

Figure 4

Location of the study area and hydrological and meteorological stations.

Data preparation

For parameterization, topography data, land use data, soil data, geological data, and observed metrological and hydrological data are required in GISMOD provided by the DIGITAL HEIHE website (http://heihenew.westgis.ac.cn).

DEM data

The size of the grid is determined by DEM data with a horizontal resolution of 1 × 1 km in this study. DEM data should be converted from other formats to ASCII to satisfy the input requirement of GISMOD, and other related data such as soil type and land use should be transformed using the same method.

Land use data

Raster data of land use were obtained with a spatial resolution of 1 km from the thematic mapper images in 1990. Land use data have been reclassified into six types referring to the sub-classification system of land cover (Figure 5).
Figure 5

Land use map.

Figure 5

Land use map.

Soil data

Spatial soil data, with a resolution of 1 × 1 km, were reclassified into three texture classes according to the soil database (ISSAS 1990) (Figure 6), and the detailed information of soil type reclassification is shown in Table 5.
Table 5

Reclassification of soil typea

Soil type Class Soil type Class 
Luvisol Primitive soil 
Semi-luvisol Saline-alkali soil 
Pedacal soil Hydromorphic soils 
Xerosol Alpine soil 
Ferrallitic soil Man-made soil 
Irrigation-silted soil Aeolian sandy soil 
Soil type Class Soil type Class 
Luvisol Primitive soil 
Semi-luvisol Saline-alkali soil 
Pedacal soil Hydromorphic soils 
Xerosol Alpine soil 
Ferrallitic soil Man-made soil 
Irrigation-silted soil Aeolian sandy soil 

aThe value of class is defined as 1, 2 and 3, which means high, medium and low infiltration capacity of soil, respectively.

Figure 6

Reclassified soil map.

Figure 6

Reclassified soil map.

Geological data

To estimate the initial values of parameters of the groundwater layer, geological data were also grouped into three types considering the permeability and porosity of the different kinds of rock (Figure 7). The details of geological type reclassification are shown in Table 6.
Table 6

Reclassification of geological typea

Geological type Class Geological type Class 
Black schist Conglomerate 
Diabase amphibolite Rhyolite 
Green schist Porphyry 
Silicalite Granite 
Sandstone Serpentinite 
Gravel Impurity 
Mudstone Andesite 
Geological type Class Geological type Class 
Black schist Conglomerate 
Diabase amphibolite Rhyolite 
Green schist Porphyry 
Silicalite Granite 
Sandstone Serpentinite 
Gravel Impurity 
Mudstone Andesite 

aThe value of class is defined as 1, 2 and 3, which means high, medium and low infiltration capacity of rock, respectively.

Figure 7

Reclassified geological map.

Figure 7

Reclassified geological map.

Observed meteorological and hydrological data

Observed daily data series from nine meteorological stations located in or around the study area were selected (Figure 4), which include daily atmospheric pressure, mean, maximum and minimum air temperature, sunshine duration, relative humidity and average wind speed data over the period 1961–2006. All these data were provided by NCC (National Climatic Centre of China). In addition, 18 precipitation stations and three hydrological stations were selected with the daily data series of precipitation and discharge provided by CAREER (Cold and Arid Regions Environmental and Engineering Research Institute, Chinese Academy of Sciences). Unfortunately, some of the daily precipitation data were only available for several non-consecutive series, which included daily series for 1964, 1971, 1976–1979, 1981, 1983–1984, 1990–1993 and 2001–2002, and therefore the period of 1990–1993 is selected in this study.

Result analysis and discussion

In order to present the model performance in the upper-middle reaches of Heihe River basin, which is divided by the Yingluoxia and Zhengyixia hydrological station, grid A (100.18, 38.81) near the Yingluoxia hydrologic station in the upper reaches and grid B (100.41, 39.13) and C (99.46, 39.81) around Zhangye City in the middle reaches were selected in this study. Simulated results of the three ordinary grids (see Figure 4) were exported for better understanding the water exchange characteristics of Heihe River basin. In addition, the simulated runoff of a river gird near Yingluoxia hydrological station was used to compare with the measured data for model validation.

Estimation of runoff

In order to evaluate the model performance, the simulated daily and monthly results of a river grid near Yingluoxia hydrological station were used to compare with the measured one in this study. The calibration period was set from 1990–1992, and 1993 was selected for validation. It can be seen that GISMOD reproduces the monthly runoff series with NSE values of 0.86 during the whole period (Figure 8). The best efficiency coefficient was 0.88 in 1993 with small differences between 1990–1992 (NSE 0.81–0.88), and the corresponding deviation in water balance was −5%. The daily runoff series with NSE value of 0.71 shows that some of the peak flows were underestimated in 1993 (NSE 0.70) and several base flows were overestimated from 1990 to 1992 (NSE 0.63–0.70). This may have resulted from the intensive seasonal variation of snow cover in mountainous regions. According to previous studies, performances of GISMOD at Yingluoxia station are acceptable in both the calibration and validation periods (He & Ge 1986; Wang & Cheng 2000; Zhu et al. 2004).
Figure 8

Comparison between the simulated and observed daily/monthly streamflow.

Figure 8

Comparison between the simulated and observed daily/monthly streamflow.

AET

It can be seen from Figure 9 that the changing trends of daily AET were almost the same among the grids. The maximum AET usually appeared in July, ranging from 1.8 to 3.9 mm, and the minimum value was close to 0.5 mm which mainly occurred in January. Shan (2007) also gave the same varying trend by combing MODIS data and the SEBS model together, the monthly AET reached a peak at July and August and then stared to decrease and reached the lowest value in December and January. Moreover, there is a significant difference in the total amount of AET among these grids. Figure 10 shows that the annual AET changed slightly during 1990–1993 at each grid, the average of which was about 422 mm, 338 mm and 395 mm at grid A, B and C, respectively. It demonstrated that AET was decreasing from the upper to the middle reaches in the Heihe River basin. Yang (2010) calculated the daily AET of Heihe River and also found that AET higher in the upper stream where it is mainly covered by forest and has plenty of rain, and it was lower in the middle and downstream. Luo et al. (2012) reported the same result based on a revised three-temperature model, which indicated that the highest AET values were from the forest or alpine meadow areas in the upper reaches, and the lowest values were from the Gobi desert in the lower reaches.
Figure 9

Simulated daily AET.

Figure 9

Simulated daily AET.

Figure 10

Simulated annual AET.

Figure 10

Simulated annual AET.

Water dynamics of different soil layers

The water content of each layer from GISMOD was exported for a better understanding of the variation trend of water resources in the different soil layers. Figure 11 shows that the water content of the surface layer was very sensitive to precipitation, especially at grid A, which was quite stable with precipitation of less than 5 mm/day and increased suddenly when the precipitation was more than 20 mm/day. The water content of the surface layer mainly fluctuated from 11 to 18% at grid B and changed from 15 to 23% at grid C.
Figure 11

Variation trend of water content in each soil layer.

Figure 11

Variation trend of water content in each soil layer.

In the soil layer, the water content varied more gently, the peak value of which had a half-month delay with the surface layer during the whole period. In general, the soil water content showed an increasing tendency at grid A and C and decreased at grid B, the average of which was about 40%, 33% and 32%, respectively.

Groundwater increased rapidly in 1990 and then remained stable from 1991 to 1993 at grid A. However, groundwater descended linearly at grid B and C, which demonstrates that the groundwater resources increased in the upper streams and is then exhausting in the middle reaches. Field monitoring data by former studies have also found that the groundwater level of the middle reach is decreasing in recent decades. For instance, Hu et al. (2007) has shown the groundwater level from four observation wells around Zhangye City, which has fallen rapidly from 1987 to 2000. Zhu et al. (2004) also found that the water-table had dropped the most in the alluvial-diluvial fan zone in the middle reaches of the Heihe River basin.

Water exchange analysis

The infiltrated and recharged water from each layer were exported to analyze the water exchange rule of the Heihe River basin. From Figure 12 it can be seen that the infiltrated and recharged water exhibits a regular trend from the upper reach to the middle reaches.
Figure 12

Pie charts of the percentage of the infiltrated and recharged water between different soil layers.

Figure 12

Pie charts of the percentage of the infiltrated and recharged water between different soil layers.

From grids A–C, the water infiltrated from the surface layer is about half of the total amount of infiltrated water, which increased from 54.27 to 65.46%, with the soil infiltration decreasing from 45.73 to 34.54%. Recharged water from the soil layer occupied the main part of recharge water, the ratio of which ranged from 65.19 to 93.8%. According to grid A, the ratio of groundwater recharge increased at grids B and C (from 6.20 to 34.81%). These results show that the groundwater layer is mainly supplied by the soil layer in the upper streams and then recharging recharged water to the soil layer in the middle streams.

For a more detailed understanding of the water exchange process, the water infiltrated from the surface layer and soil layer were summed up as the infiltrated water to compare with the recharged water, which was considered the total volume of water supplied from the groundwater layer and the soil layer. The statistical results from these three sites are shown in Table 7.

Table 7

Statistic data of the infiltrated and recharged water (monthly) (unit: mm)

Grid Item Max Min Mean Annual average 
Infiltrated 13.47 7.02 13.47 102.29 
Recharged 7.11 1.92 5.52 66.21 
Infiltrated 9.61 6.49 7.68 92.15 
Recharged 10.83 5.39 8.15 97.76 
Infiltrated 15.11 11.23 14.09 169.17 
Recharged 19.39 13.08 15.31 183.71 
Grid Item Max Min Mean Annual average 
Infiltrated 13.47 7.02 13.47 102.29 
Recharged 7.11 1.92 5.52 66.21 
Infiltrated 9.61 6.49 7.68 92.15 
Recharged 10.83 5.39 8.15 97.76 
Infiltrated 15.11 11.23 14.09 169.17 
Recharged 19.39 13.08 15.31 183.71 

The annual average of infiltrated water is nearly 1× larger than the recharged water at grid A, where few water resources were supplied from the groundwater (ranging from 1.92 to 7.11 mm/month). As opposed to grid A, the recharged water was larger than the infiltrated water at grid B and C, which illustrates that the infiltrated water was reduced, with the recharge water increasing from 66.21 to 183.71 mm/year in the middle streams of the Heihe River basin. As proposed by Gao & Li (1990) and Cheng & Qu (1992), the runoff was mainly generated by rainfall, infiltration, glaciers and snowmelt in the upper stream and produced a great amount of surface water leakage down at the foot of Qilian Mountain with less groundwater recharge (Figure 13). Along with the stream, the aquifer system began to recharge at the front of Zhangye City and then mainly discharged on the pluvial fans of the middle stream (Wang & Cheng 2000; Hu et al. 2007). Moreover, a case study on groundwater and river water exchange using a tracer test also demonstrated that the groundwater recharge increased from Yingluoxia to Gaotai (Zhang et al. 2005; Song et al. 2006). Therefore, the simulated results of GISMOD are generally consistent with previous studies and show the special water exchange characteristics in the upper-middle reaches of Heihe River.
Figure 13

Transformation of precipitation, surface water and groundwater in the Heihe river basin (Zhao et al. 2011).

Figure 13

Transformation of precipitation, surface water and groundwater in the Heihe river basin (Zhao et al. 2011).

In addition, monthly results are given for analysis of the seasonal dynamics of water exchange (Figure 14). The change tendency of infiltrated water and recharged water is quite different at all three sites. When infiltrated water decreases, recharged water increases, and vice versa. In summer (from June to August), a large amount of water leaks down from the surface layer while little water could be supplied from the groundwater with the ratio of recharged to infiltrated ranging from 0.38 to 1.06 (Table 8). The most likely reason is that the soil water content becomes saturated after continuous rainfall in the upper reaches, which may increase the infiltration and decrease the replenishment. From Figure 14 it can be seen that water recharge took place mainly in winter and spring at all these sites. Wang et al. (2009) indicated that there is a general descent of groundwater in winter and the following spring in the Heihe River basin. Through a seasonal observation of δ18O, Akiyama et al. (2007) found that the river water was supplied by groundwater discharge in the middle reaches of Heihe River during the non-irrigation period with little precipitation. Furthermore, the results of Zhao et al. (2011) show the same seasonal dynamic of water exchange where the recharge of surface runoff occurred mainly from June to mid-September, whereas the supply of surface runoff in winter was from the base flow. Based on the above analysis, it can be concluded that the interaction between surface water and groundwater was significantly different by season and location in the Heihe River basin. Groundwater is mainly recharged in the upper reaches and then discharged in the middle reaches with the maximum of recharged water appearing in winter and the minimum value appearing in summer.
Table 8

Comparison between the infiltrated and recharged water in different seasons (unit: mm)

  Spring Summer Autumn Winter 
Grid Variables (Mar–May) (Jun–Aug) (Sep–Nov) (Dec–Feb) 
Recharged/infiltrated (ratio) 6.29/7.88 3.86/10.17 5.68/8.29 6.24/7.75 
79.76% 38% 68.45% 80.52% 
Recharged/infiltrated (ratio) 8.32/7.62 7.58/8.13 8.37/7.47 8.31/7.49 
109.19% 93.26% 112.07% 110.95% 
Recharged/infiltrated (ratio) 15.34/14.03 15.21/14.27 15.46/13.99 15.22/14.10 
109.31% 106.63% 110.53% 107.95% 
  Spring Summer Autumn Winter 
Grid Variables (Mar–May) (Jun–Aug) (Sep–Nov) (Dec–Feb) 
Recharged/infiltrated (ratio) 6.29/7.88 3.86/10.17 5.68/8.29 6.24/7.75 
79.76% 38% 68.45% 80.52% 
Recharged/infiltrated (ratio) 8.32/7.62 7.58/8.13 8.37/7.47 8.31/7.49 
109.19% 93.26% 112.07% 110.95% 
Recharged/infiltrated (ratio) 15.34/14.03 15.21/14.27 15.46/13.99 15.22/14.10 
109.31% 106.63% 110.53% 107.95% 
Figure 14

Amount of water exchanged at the three grids (monthly).

Figure 14

Amount of water exchanged at the three grids (monthly).

CONCLUSIONS

GISMOD is a grid-based distributed hydrological model which has been designed to furnish a description of hydrological processes with the emphasis on surface-groundwater exchange in a catchment. A specific method based on a water balance equation is used to simulate the hydrological processes of different soil layers by taking into account precipitation, evapotranspiration and infiltration in GISMOD. The advantage of this model is that water interaction between different soil layers can be simply described by a generalized physical process under various situations.

A case study was performed in the upper-middle reaches of the Heihe River basin to demonstrate the applicability of the GISMOD. Compared with the observed runoff data, we found that both calibration and validation of GISMOD were performed with reasonable accuracy at the daily and monthly scales. The hydrological processes in the different soil layers are generally reflected by GISMOD, which also gives the seasonally changing laws of water exchange in the Heihe River basin. Generally, the GISMOD can estimate hydrological processes and model the spatial and temporal variations of the infiltrated water and recharged water in the study area.

GISMOD is a modular modeling system developed independently and primarily for a program task, and it needs to be improved in future studies to add the module of human water use, automatic parameter calibration, snowmelt simulation and so on.

ACKNOWLEDGEMENTS

This study was supported by the Major Research Plan of the National Natural Science of China (No. 91125015). Thanks are given to the Environmental and Ecological Science Data Center for West China (http://westdc.westgis.ac.cn) and the Digital River Basin (http://heihe.westgis.ac.cn) for providing the necessary data.

REFERENCES

REFERENCES
Akiyama
T.
Sakai
A.
Yamazaki
Y.
Wang
G. X.
2007
Surface water-groundwater interaction in the Heihe River basin, Northwestern China
.
Bull. Glaciol. Res.
24
,
87
94
.
Allen
R. G.
Pereira
L. S.
Raes
D.
Smith
M.
1998
Crop Evapotranspiration: Guidelines for Computing Crop Water Requirements
.
FAO Irrig. Drain. Paper No. 56. FAO
,
Rome
,
Italy
, p.
300
.
Beven
K. J.
Lamb
R.
Quinn
P. F.
Romanowicz
R.
Freer
J.
1995
TOPMODEL
. In:
Computer Models of Watershed Hydrology
(
Singh
V. P.
, ed.).
Water Resources Publications
,
Highlands Ranch, Colorado
, pp.
627
668
.
Cheng
L. H.
Qu
Y. G.
1992
Water and Land Resources and their Rational Development and Utilization in the Hexi Region
.
Science Press
,
Beijing
.
Clark
M. P.
Kavetski
D.
Fenicia
F.
2011
Pursuing the method of multiple working hypotheses for hydrological modeling
.
Water Resour. Res.
47
,
W09301
.
Gao
Q. Z.
Li
F. X.
1990
The Rational Utilization of Water and Land Resources in Hei River Basin
.
Gansu Science & Technology Press
,
Lanzhou
.
Gardner
K. M.
1999
The Importance of Surface Water/Groundwater Interactions. US Environmental Protection Agency Issue Paper: EPA-910-R-999-013
.
Graham
N. R. A.
2001
MIKE SHE: a distributed, physically based modeling system for surface water/groundwater interactions
. In:
Proceedings of MODFLOW 2001 and Other Modelling Odysseys. Golden, Colorado
.
Graham
D. N.
Butts
M. B.
2005
Flexible, integrated watershed modelling with MIKE SHE
. In:
Watershed Models
(
Singh
V. P.
Frevert
D. K.
, eds).
CRC Press
,
Boca Raton, FL
,
USA
, pp.
245
272
.
Hargreaves
G. H.
Samani
Z. A.
1985
Reference crop evapotranspiration from temperature
.
Appl. Eng. Agric.
2
,
96
99
.
He
S.
Ge
S.
1986
Groundwater dynamically changes and its influence on the environment in Jinta irrigation area (in Chinese)
.
J. Desert Res.
6
,
38
49
.
Hu
L.-T.
Chen
C.-X.
Jiao
J. J.
Wang
Z.-J.
2007
Simulated groundwater interaction with rivers and springs in the Heihe river basin
.
Hydrol. Process.
21
,
2794
2806
.
Imagawa
C.
Takeuchi
J.
Kawachi
T.
Ishida
K.
Chono
S.
Buma
N.
2011
A hydro-environmental watershed model improved in canal-aquifer water exchange process
.
Paddy Water Environ.
9
,
425
439
.
Institute of Soil Science, Academia Scinica (ISSAS)
1990
Soils of China
.
Science Press
,
Beijing
.
Kalbus
E.
Reinstorf
F.
Schirmer
M.
2006
Measuring methods for groundwater, surface water and their interactions: a review
.
Hydrol. Earth Syst. Sci. Discuss.
3
,
1809
1850
.
Khakbaz
B.
Imam
B.
Sorooshian
S.
Koren
V. I.
Cui
Z.
Smith
M. B.
Restrepo
P.
2011
Modification of the National Weather Service Distributed Hydrologic Model for subsurface water exchanges between grids
.
Water Resour. Res.
47
,
W06524
.
Kim
N. W.
Chung
I. M.
Won
Y. S.
Arnold
J. G.
2008
Development and application of the integrated SWAT–MODFLOW model
.
J. Hydrol.
356
,
1
16
.
Kosmas
C.
Danalatos
N.
Cammeraat
L. H.
Chabart
M.
Diamantopoulos
J.
Farand
R.
Gutierrez
L.
Jacob
A.
Marques
H.
Martinez-Fernandez
J.
Mizara
A.
Moustakas
N.
Nicolau
J. M.
Oliveros
C.
Pinna
G.
Puddu
R.
Puigdefabregas
J.
Roxo
M.
Simao
A.
Stamou
G.
Tomasi
N.
Usai
D.
Vacca
A.
1997
The effect of land use on runoff and soil erosion rates under Mediterranean conditions
.
Catena
29
,
45
59
.
Krause
S.
Hannah
D. M.
Fleckenstein
J. H.
2009
Hyporheic hydrology: interactions at the groundwater-surface water interface
.
Hydrol. Process.
23
,
2103
2107
.
Li
L.
Xu
Z. X.
Li
Y. L.
Liu
W. F.
Zhang
L. Y.
2013
A preprocessing program for a distribute hydrological model: Development and application
.
J. Hydroinform.
15
,
1258
1275
.
Michaud
J.
Sorooshian
S.
1994
Comparison of simple versus complex distributed runoff models on a midsized semiarid watershed
.
Water Resour. Res.
30
,
593
605
.
Monteith
J. L.
1965
Evaporation and environment
. In:
Symposium of the Society for Experimental Biology, The State and Movement of Water in Living Organisms
(
Fogg
G. E.
, ed.).
Academic Press, Inc.
,
New York
, pp.
205
234
.
Nash
I. E.
Suthcliff
I. V.
1970
River flow forecasting through conceptual models
.
J. Hydrol.
10
,
282
290
.
Qiu
Y.
1987
The transformation between surface water and groundwater and the calculation of the total water resources in Hexi region (in Chinese)
.
Nat. Resour.
2
,
7
15
.
Refsgaard
J. C.
Storm
B.
1995
MIKE SHE
. In:
Computer Models of Watershed Hydrology
(
Singh
V. P.
, ed.).
Water Resources Publications
,
Highlands Ranch, Colorado
, pp.
809
846
.
Shan
X.
2007
Regional EVAPTRANSPIRATIOn over Arid Inland Heihe River Basin in Northwest China
.
ITC
,
Enschede
,
The Netherlands
.
Swain
E. D.
Wexler
E. J.
1996
A coupled surface-water and ground-water flow model (MODBRANCH) for simulation of stream-aquifer interaction
. In:
U.S. Geological Survey Techniques of Water Resources Investigations (USGS, ed.), Book 6, Government Printing Office, Washington, DC, p. 125.
Tang
Q.
Qiu
Y.
Zhou
L.
1992
Hydrology and Water Resource Use in Arid Regions of China
.
Science Press
,
Beijing
,
China
,
195
pp.
TCADC
1984
Technical Rules of Land Use Survey
.
The Committee Agricultural Divisions of China
,
Beijing
.
Thornes
J. B.
Shao
J. X.
Diaz
E.
Roldan
A.
McMahon
M.
Hawkes
J. C.
1996
Testing the MEDALUS hillslope model
.
Catena
26
,
137
160
.
Turc
L.
1961
Estimation of irrigation water requirements, potential evapotranspiration: a simple climatic formula evolved up to date
.
Ann. Agron.
12
,
13
49
.
Xu
Z. X.
2009
Hydrological Models
.
Science Press
,
Beijing
,
China
.
Xu
Z. X.
Ito
K.
Schultz
G. A.
Li
J. Y.
2001
Integrated hydrologic modeling and GIS in water resources management
.
J. Comput. Civil Eng.
15
,
217
223
.
Yang
Y. M.
2010
Estimation of Evapotranspiration in Heihe River Basin Based on Remote Sensing (in Chinese)
.
Lanzhou University
,
Lanzhou
,
China
.
Yuan
L.-R.
Xin
P.
Kong
J.
Li
L.
Lockington
D.
2011
A coupled model for simulating surface water and groundwater interactions in coastal wetlands
.
Hydrol. Process.
25
,
3533
3546
.
Zhang
Y.
Wu
Y.
Ding
J.
Wen
X.
2005
Exchange of groundwater and river water in a basin of the Middle Heihe River by using δ18O (in Chinese)
.
J. Glaciol. Geocryol.
27
,
106
110
.