Groundwater response to leakage of surface water through a thick vadose zone in the middle reaches area of Heihe River Basin , in China

The behavior of groundwater response to leakage of surface water in the middle reaches area of Heihe River Basin is significantly influenced by a thick vadose zone. The groundwater regime is a result of two recharge events due to leakage of Heihe River and irrigation water with different delay time. A nonlinear leakage model is developed to calculate the monthly leakage of Heihe River in considering changes of streamflow, river stage and agricultural water utilization. Numerical modeling of variable saturated flow is carried out to investigate the general behaviors of leakage-recharge conversion through a thick vadose zone. It is found that the recharge pattern can be approximated by simple reservoir models of leakages under a river and under an irrigation district with different delay-time and recession coefficient. A triple-reservoir model of relationship between surface water, vadose zone and groundwater is developed. It reproduces the groundwater regime during 1989–2006 with variable streamflow of Heihe River and agricultural water utilization. The model is applied to interpret changes of groundwater level during 2007–2008 that observed in the Watershed Airborne Telemetry Experimental Research (WATER). Correspondence to: X.-S. Wang (wxsh.cugb@gmail.com)


Introduction
The interaction between groundwater and surface water plays an important role in water balance in most of hydrological basins around the world.An aquifer can obtain recharge from leakage of rivers and irrigation system in a basin.However, in the same basin groundwater can also be the source of base flow for streams, lakes and wetlands.This groundwater and surface water interaction is one of the major components of hydrological process to be observed and investigated in the Watershed Airborne Telemetry Experimental Research (WA-TER) which is carried out in the Heihe River Basin (Li and Cheng, 2008;Li et al., 2008), the second largest inland river basin in the arid region of northwest China.
Interpretation of groundwater regime observed in WATER in the middle reaches area of Heihe River Basin needs a wellestablished model of groundwater recharge and discharge.The typical process of groundwater and surface water interaction is necessary to be considered.Several numerical models of groundwater flow in the middle reaches area of Heihe River Basin developed in the literature can be selected for the purpose.Zhou et al. (1990) established a single-layer two-dimensional finite-difference model (Model-1990).Two numerical models were developed by Zhang et al. (2004) with finite-element method: a single-layer model (Model-2004-1) and a double-layers model (Model-2004-2).Another double-layers model (Model-2006) was also developed by Jia et al. (2006) as a part of the hydrological model WEP-Heihe.The numerical model developed by Hu et al. (2007), which is called Model-2007 in this presentation, is a three-Published by Copernicus Publications on behalf of the European Geosciences Union.
dimensional model with 8 model-layers.In these 5 models, discharge of groundwater into Heihe River is handled by using a constant head boundary (Model-1990, Model-2004-1 and Model-2004-2) or a river-aquifer interaction model (Model-2006 andModel-2007).In Model-2006, a linear formula is applied which assumes that exchange of water between the river and the aquifer is proportional to the difference between the river stage and groundwater level.This linear relationship is according to Rushton and Tomlinson (1979).In Model-2007, horizontal discharge besides vertical flow through the riverbed is also considered while groundwater level is higher than the river stage.The riveraquifer interaction models are successful when surface water and groundwater are directly connected with each other.However, recharge of groundwater through infiltration of surface water influenced by vadose zone is still a problem.In Model-1990, Model-2004-1 and Model-2004-2, infiltration rate is assumed to be equal to recharge rate and is directly input to groundwater model.Therefore, the impact of vadose zone is ignored in the models.In Model-2006 andModel-2007, hypothetical equations are applied to redistribute infiltration rate for recharge rate in considering of the hysteresis between infiltration and recharge.However, the methods are not well expatiated on and these hypothetical equations need to be verified.
Unsaturated flow in the vadose zone is a key factor for infiltration-recharge conversion in the Middile reaches area of Heihe River Basin.This is due to deep buried groundwater table under fluvial fans and irrigation districts where rivers and channels loss surface water through leakage.As released from experiment study with two-dimensional sand box (Vauclin et al., 1979), significant delay time and redistribution of infiltration can be resulted in while surface water is transferred to groundwater through vadose zone.Chen (2007) reports an investigation of stream-aquifer interactions in southcentral Platte River valley, Nebraska.At the site, groundwater level varies seasonally from about 1.0-1.5 m below ground surface.During a rising-falling event of river stage, it is found that the time lag is 0.5 h for the initial response of water table in the aquifer to the rise of stream stage and peak groundwater level came 1.5 h later than the peak crest of streamflow.In a deeper water table condition, Dahan et al. (2007) measured ground water response to flood events in a river beneath a shallow alluvial aquifer.The depth of water table is 3.6-4.0m and rising of groundwater level is sensed 5-6 h after flooding initiation.Indirect recharge from stream leakage is also found for 10-m-depth groundwater table under a sandy alluvial fan (Massuel et al., 2006).The impact of vadose zone in the middle reaches area of Heihe River Basin is more significant than these reported examples because depth of groundwater table in this area can be larger than 100 m.
It is difficult to completely handle the in situ unsaturated flow with moisture measurement in the large scale fluvial fans which is composed of heterogeneous media.How-ever, the characteristics of infiltration-recharge conversion can be recognized through change of groundwater level as response of leakage events of surface water.In study of the upper river Lys basin and the Origny aquifer with an unsaturated zone (depth of water table is 0-50 m), Besbes and De Marsily (1984) suppose a method in recognizing infiltration-recharge conversion with observed rainfall and groundwater level by using a parametric transfer function.The transfer function method is also applied by Zammouri and Feki (2005) in study of artificial groundwater recharge through leakage from reservoirs on ephemeral streams.A model including several linear reservoirs (Nash, 1957) is developed to represent the unsaturated zone which influences the leakage-recharge conversion under streams.The number of the equivalent reservoirs and the recession coefficient of each reservoir are calibrated with observed groundwater level of depth 25-35 m.These previous works suggest that the transfer function is valuable in practice to produce simple analog models of groundwater recharge with a thick vadose zone.However, in the middle reaches area of Heihe River Basin, it is noted in WATER that the change of groundwater level under the fluvial fans is not a simplex response to leakage of surface water, i.e., a transfer function is not enough to configure the whole recharge patterns.
In this presentation, the behavior of groundwater response to leakage of surface water in the middle reaches area of Heihe River Basin is investigated.The purpose is to release the general behaviors of leakage-recharge conversion through a thick vadose zone and check the validation of the transfer function method.A simplified model of relationship between surface water, vadose zone and groundwater is developed to reproduce groundwater regime during a period over 10 years with variable streamflow of Heihe River and agricultural water utilization.The model is applied to interpret changes of groundwater level observed in WATER.

General
Heihe River is the primary stream in Heihe River Basin.As shown in Fig. 1a, it runs down from Qilian Mountain (upper reaches, A) to Zhangye basin on Hexi Corridor(middle reaches, B1), then travels to the Ejina Basin in inner Mongolia (lower reaches, C) and finally stops at Juyan Lakes.The middle segment of Heihe River is controlled by two outlets: Yingluo Gorge (at an altitude of about 1600 m) and Zhengyi Gorge (at an altitude of about 1300 m).As shown in Fig. 2, during 1989-2005, the monthly streamflow of Heihe River at Yingluo Gorge is 10-220 m 3 /s and the corresponding annual streamflow is 11-23×10 8 m 3 /a.
In the middle reaches area of Heihe River Basin shown as B1 and B2 in Fig. 1a, mean annual precipitation is less than 200 mm and mean annual potential evaporation is greater Hydrol.Earth Syst.Sci., 14, 639-650, 2010 www.hydrol-earth-syst-sci.net/14/639/2010/ than 2000 mm in this area.However, oases are strongly developed by introducing streamflow into channels to irrigate farmlands.Rivers come from Qilan Mountain contribute about 30×10 8 m 3 of water resources annually to the middle and lower reaches area.About 2/3 of the water resources is consumed in the middle reaches area.
Groundwater plays an important role in terrestrial hydrology and utilization of water resources in this area.Firstly, a great amount of surface water leaks down as recharge of groundwater on pluvial fans in front of the foot of Qilian Mountain.Secondly, groundwater flows out through springs or discharges into Heihe River in the lower alluvial plain.In addition, groundwater is directly abstracted from aquifers by over 6000 wells for domestic, agricultural and industrial water demand.
The aquifer system in this area is a composite of Quatemary sediments with total thickness of 20-1800 m.Sandy gravel and sand with calyey interbeds are included.Underlying Tertiary sediments of sand-stones and shalestones are normally regarded as the impermeable bed of this groundwater basin.Thick gravel sand sediment exhibits under the pluvial fans.It forms a high permeable aquifer for groundwater flow and percolation of water in the vadose zone from leakage of surface water.

Groundwater observation
On the west of Zhangye City, as shown in Fig. 1b, there are 8 observation wells allocated near Heihe River with data of groundwater level during 1989-2005.Variations of groundwater level at five typical wells are shown in Fig. 3. Depth of water table in the aquifer observed in the wells ranges between 4-81 m, as shown in Fig. 1c.However, in the south of the site where Heihe River runs out of Yingluo Gorge, groundwater level is 100-300 m below the riverbed.
It is found that variation patterns of streamflow and groundwater level are significantly different.As shown in Fig. 4, the highest streamflow occurs during midsummer but the highest groundwater level occurs during autumn every year.In the winter and the following spring, a general descending of groundwater level follows decrease of streamflow in the autumn.However, during December and following January and February when streamflow is very small, drop of groundwater level is significantly slowed down (in 1993(in -1994, Fig. 4, Fig. 4) or broken by a rising-falling event (in 1994(in -1995, Fig. 4, Fig. 4).It is a temporary disturb and the drop of groundwater level will continues until June.
This special relationship between streamflow and groundwater level is due to the special river-aquifer interactions at the site.The involved hydrological processes are leakage of streams, leakage of irrigation water, unsaturated flow in the vadose zone and groundwater discharge.

Leakage of Heihe River
In the south area of Zhangye-Jiuquan Basin as shown in Fig. 1a, more than 20 rivers provide surface water that originates form Qilian Mountain.Most of the rivers loss their total flow by leakage and irrigation before they arrive in the lower reaches area.Heihe River is the largest one which provides the most leakage of surface water in Zhangye basin.
Heihe River mainly leaks through the riverbed along segment between Yingluo Gorge and Zhangye Bridge, as shown in Fig. 3.This leakage is dependent on streamflow, river stage and properties of sediment on the riverbed.Investigation of leakage of Heihe River has been carried out in 1967, 1985 and 2001-2003  when channels are closed.Total leakage along a measured segment of the river is calculated by difference between inflow and outflow.It is found that the total leakage is a function of inflow and length of the segment.

by measuring change of streamflow
A simplified model is presented in this study to simulate the leakage of Heihe River, as shown in Fig. 5.In this model, decreasing of streamflow due to leakage is described by a nonlinear relationship: where Q is streamflow, L is the distance along the river, h is the river stage, ξ 1 and η 1 are two factors.The effect of evaporation over water table is ignored in Eq. ( 1).The maximum annual evaporation loss is calculated as 0.042×10 8 m 3 /a according to length of the investigated river segment.It is about 0.38% of the minimum annual streamflow at Yongluo Gorge (11.0×10 8 m 3 /a).Therefore, the relative error of Eq. ( 1) by ignoring evaporation is less than 0.4% which is acceptable in this study.
The relationship between the river stage and the streamflow can be generally described by an empirical formula as follows where ξ 2 and η 2 are two factors.Substituting Eq. ( 2) into Eq.
In accounting for monthly behavior of the streamflow, a quasi-steady state flow is assumed in each month so that ∂Q/∂L in Eq. ( 3) can be approximated by dQ/dL.For constant ξ and η, solution of Eq. ( 3) as outflow, Q B in Fig. 5, is obtained as where L 0 is the threshold distance dependent on inflow, Q A , and can be calculated by It means that the river comes to be dry after leakage along the threshold distance.Total leakage along segment A-B, F AB , is estimated by For leakage between Yingluo Gorge and Zhangye Bridge, L AB =30 km.If irrigation channels are closed, the actual leakage of Heihe River can be calculated as difference in Hydrol.Earth Syst.Sci., 14, 639-650, 2010 www.hydrol-earth-syst-sci.net/14/639/2010/ streamflow between the two measurement sites.It is found that ξ =0.09 km −1 (for unit of streamflow is m 3 /s) and b=0.54 leads to a best agreement between estimation and observation in this situation.
In considering of irrigation water introduced from Heihe River, streamflow after export to channels at site-1 and site-2 (Fig. 5) are reduced to where Q y is the measured streamflow at Yingluo Gorge, I 1 is the total flow diverted to CW1 and CE1, I 2 is the total flow diverted to CW2 and CE2, F 12 is the leakage of Heihe River between site-1 and site-2.Eqs.(4-7) are applied to estimate F 12 by replacing A to 1 and replacing B to 2. Distance between site-1 and site-2 is L 12 =14 km.Streamflow at site-3 (Fig. 5) is determined by where F 23 is the leakage of Heihe River between site-2 and site-3.Equations (4-7) are applied to estimate F 23 by replacing A to 2 and replacing B to 3. Distance between 2 and 3 is L 23 =16 km.The total leakage of Heihe River is given by Export of water from Heihe River to channels is input to I 1 and I 2 in Eqs. ( 8) and ( 9) to calculate the total leakage.As an example, the calculated monthly leakage of Heihe River according to measured streamflow at Yingluo Gorge and export flow for irrigation in 1996 is shown in Table 1.It indicates that in this year 27.0% of total flow in Heihe River lost with leakage.As shown in Fig. 6, this monthly leakage of Heihe River is a function of excess streamflow and can be approximated by an empirical equation as follows where F 0 =0.164×10 8 m 3 /month.F 0 is a reference leakage rate introduced herein to form a dimensionless power function in Eq. ( 12).It is found that annual leakage of Heihe River is 15%-30% of the annual streamflow at Yingluo Gorge during 1989-2008.

Numerical modeling of leakage-recharge conversion through a thick vadose zone
In the middle reaches area of Heihe River Basin, recharge of groundwater is mainly contributed by leakage of streams and leakage of irrigation water.Before moving into groundwater, the leaked surface water must pass through the vadose zone under pluvial fans and irrigation districts.Thickness of the vadose zone can be larger than 100 m.Unsaturated flow in the vadose zone significantly influences the timelag behavior of leakage-recharge conversion and results in a special pattern of groundwater response as shown in Fig. 4.However, few in situ observations on the unsaturated flow in the thick vadose zone are available for analysis on leakagerecharge conversion.In this study, numerical models are applied to investigate the impact of vadose zone on leakage-recharge conversion.Due to complexity in spatial heterogeneity of the vadose zone and lack of hydraulic data of the soils, it is difficult to carry out a real world model of unsaturated flow in this area.However, it is possible to find general behaviors of the leakagerecharge conversion by using well-defined models with simple boundary conditions, uniform isotropic porous media and empirical parameters.The focus in this section is the timelag feature of recharge following leakage under a river bed or an irrigation district.

Equations of flow in the vadose zone
Richards equation based on hydraulic head is applied to describe variable saturated flow in the vadose zone, as follows where H hydraulic head given by H = h − z, m of water; h soil water pressure, m of water; for unsaturated condition, h <0; where θ is the volumetric water content, θ s and θ r are saturated and residual moisture content, respectively; α, n, m are factors and there is m=(n−1)/n.θ s =0.25, θ r =0.05, α=1.2 m −1 and n=1.5 are applied in this study.Hydraulic conductivity is determined with an empirical equation of Gardner (1958) where K s is the saturated hydraulic conductivity, c k is a factor.K s =10 m/d and c k =0.5 m −1 is applied in this study.
In this study, Eq. ( 13) is numerically solved with an open source code, AquiferFlow-2007(Wang, 2008).AquiferFlow-2007 deals with a three-dimensional numerical model for variable saturated flow based on block-centered finite difference method.The parameters are specified according to normal empirical conditions of moderate sand soil as an approximation of the sediments under Heihe River.However, the exact parameter values are not available currently.It is not a trouble because the model is built to analyze the general behaviors of leakage-recharge processing.

Leakage-recharge conversion under a river
In this section, a model is built up to investigate the patterns of groundwater recharge from leakage of a river through a thick vadose zone.It is a two-dimensional model on profile.The river is 2-m-width, extends infinitely along a beeline with the same aquifer condition on both sides.In considering the symmetry of the problem, only the right part of cross section originated from the center of the river is treated in the model.As shown in Fig. 7a-b, the model is composed of a given head cell on the left-top in representing the river, the top boundary with constant uniform infiltration of precipitation and the 100-m-depth bottom boundary with constant head that equal to −100 m.Initial head of the given head cell is 0 m which means the river stage is equal to the top of the model.Infiltration of precipitation is 0.001 m/d.The bottom boundary approximately represents a ground water table to accept recharge from leakage of the river.
An initial steady state flow is firstly simulated while river stage is steady on 0 m. Figure 7b shows the initial distribution of volumetric water content in the vadose zone.Under the river, it can be seen that a thick unsaturated wet zone is developed between the top saturated zone and the bottom saturated zone.Moisture content in the wet zone is 0.13-0.25,larger than that far away the river in horizontal where x > 18 m.Figure 7c shows the distributions of initial flow velocity along horizontal axis at different depth.They approximately represent a group of gauss curve with different highest velocity below center of the river.With increasing of depth, both the highest velocity and deviation of the curve decrease.It indicates a diffusion of infiltration with depth.
A rising-falling event of river stage is specified on the given head cell to investigate the recharge response on variable leakage.The river stage rises up from 0 m to 2 m in five days and then falls down from 2 m to 0 m in the next five days.This transient flow is simulated with a time step of 0.2 d.Redistribution of volumetric water content on the fifth day is shown in Fig. 7d.Comparing with the initial condition as shown in Fig. 7b, this rising of river stage enlarges the saturated zone and the wet zone (θ=0.13-0.25)under the river.Range of the wet zone above z=90 m increases to x=25 m.It is found that this wet zone is diminished gradually during falling of river stage.Variations of total infiltration at different depth response to the rising-falling event can be seen in Fig. 8.The total flow rate at one depth is an accumulation of Darcy velocity at the depth along the horizontal axis.The curve of z=1 m gives the variable leakage of the river.Change of groundwater recharge can be approximately represented by variation of the total flow rate at z=90 m, as shown in Fig. 8. Obviously, a significant delay exhibits in response of the recharge to the variable leakage.The peak of recharge rate occurs at t=6.5 d, which is 1.5 days after the peak of leakage rate.As proposed by Besbes and De Marsily (1984) and Zammouri and Feki (2005), the effect of vadose zone to groundwater recharge can approximately described by a transfer function.For a thick vadose zone, the transfer function of several linear reservoirs can be considered.However, the number of the reservoirs is not a well-defined parameter though it is explained that one reservoir represents a vadose zone layer of several meters.To satisfy the nature of infiltration-recharge conversion, fractional value of the number of the reservoirs is allowed.In this study, a different way in description of the time-lag behavior of groundwater recharge is proposed.It is found from the numerical results that the leakage-recharge conversion under a river can be simply simulated by a modified linear reservoir (tank) model described as follows where F (t) is the leakage/infiltration rate on the top, R is the recharge rate of groundwater, a is the recession coefficient, t and τ are the current time and the delay-time, respectively.The delay-time indicates how quickly the recharge over wa- ter table begins to answer the change of leakage.It is subject to the initial moisture and the thickness of the vadose zone.Equation ( 16) can be numerically solved for arbitrary F (t) with finite difference method as follows where F t is the discrete leakage rate at time t, R t+τ and R t+τ + t are the discrete recharge rate at time t+τ and t+τ + t, respectively, t is time step.In the numerical scheme, the arbitrary function, F (t), is approximated by discrete value at time t 0 , t 0 + t, t 0 +2 t, . . ., and R(t) is calculated with Eq. ( 17) for time t 0 +τ , t 0 +τ + t, t 0 +τ +2 t, . . ., where t 0 is a beginning time.
The result of recharge in the previous numerical modeling is analyzed with the simple reservoir ( tank

Leakage-recharge conversion under an irrigation district
In this section, a vertical one-dimensional model is built up to investigate the patterns of groundwater recharge from leakage of irrigation water through a thick vadose zone.This problem is a little different from leakage of a stream.Leakage under a river leads to a concentrated infiltration (a linear source in space or a point source in profile) and then horizontally spreads in the vadose zone with increasing of depth.However, leakage of irrigation water occurs in a larger scale on the ground surface as an areal source.It is an excess of irrigation water upon evaporation and transpiration within the root zone.Vertical flow in the vadose zone is the dominative process of leakage-recharge conversion under these areal sources.Therefore, a vertical flow model is considered.
As shown in Fig. 9a, the model is composed of a top boundary of given infiltration rate and a bottom boundary of a draining cell.The draining cell is assumed to be connected with a constant head body of water table at z=80 m.Drainage of groundwater out off the cell is given by where q is the discharge rate, H cell is the hydraulic head in the draining cell, H 0 is the elevation of the constant head body, C d is the conductance.Transport of water in the vadose zone is simulated in the numerical model.Figure 9b shows the redistribution of volumetric water content in leakage-recharge conversion.Before t = 10 d, a wetting front exhibits and moves gradually toward the groundwater table without any response of water table.After t = 16 d, the wetting front disappears quickly following with a rising of groundwater table and a decreasing of water content near ground surface.
Recharge of groundwater in this model can be estimated by balance of water in the saturated zone: where z g is the elevation of water table where h = 0, it can be approximated by H cell in this model.Substituting Eq. ( 18) into Eq.( 19) we have Therefore, the recharge rate can be numerically calculated by where j is the time step.The result of recharge as a response of the infiltration event is shown in Fig. 10.It is obviously that the recharge process is delayed and is smoother than the infiltration process.Similar to the leakage-recharge conversion under a river, response of recharge to the infiltration event can also be approximated by a simple reservoir model described in Eq. ( 16) and solved with Eq. ( 17).For this vertical flow problem, τ = 17.9 d and a = 0.059 d −1 are identified.As shown in Fig. 10, the time of peak recharge rate calculated by the reservoir model is about 3 days later than that simulated by the numerical model.However, shapes of the two recharge curves are quite similar to each other.It indicates that the reservoir model is valuable to simplify analysis on leakage-recharge conversion through a thick vadose zone under an irrigation district.

A triple-reservoir model of groundwater response to leakage of surface water
Previously, with a two-dimensional model on profile and a one-dimensional vertical model, the leakage-recharge conversions through a thick vadose zone under a river and an irrigation district are investigated with finite difference method.The numerical modeling is undertaken with assumptions of isotropic and homogenous porous media, simple boundary conditions and empirical parameters.However, they provide the knowledge of the impacts of a thick vadose zone on leakage-recharge conversion.Similarity has been found between leakage-recharge conversion through a vadose zone and input-output transfer within a simple reservoir model.Accordingly, an alternative method besides a complex real world numerical model is suggested in this study to investigate groundwater response to leakage of surface water.The method applies reservoir models in simulation of water transfer through the vadose zone.A triple-reservoir model is established to describe interactions among the river, the irrigation districts, the vadose zone and the aquifer system.This model can deal with the essential characteristics of hydrological process at the studied site.

Description of the triple-reservoir model
The investigation is concentrated in the Ganzhou district around Zhangye City as shown in Fig. 1b.Heihe River between Yingluo Gorge and Zhangye Bridge (Fig. 5) is included.As shown in Fig. 11, the triple-reservoir model is composed of three reservoirs representing the vadose zones under Heihe River (S 1 ) and under the irrigation district (S 2 ) and the aquifer system (S g ), respectively.
In the Ganzhou district, the total stream flow, Q y , at Yingluo Gorge is divided into three parts: export flow to channels, leakage through river bed and excess flow for the lower reaches.Export flow to channels includes Q i and Q f where Q i is export flow to farmlands in the Ganzhou district and Q f is export flow to other irrigation districts 40 km away from Heihe River in the west.There is where I 1 and I 2 are defined in Eqs.(8-9).Leakage of Heihe River is calculated with Eq. ( 12).Both leakage of Heihe River and leakage of trunk channels are given to the first reservoir (S 1 ).Leakage of small channels which densely distribute in the area of farmlands is included in the leakage of irrigation water and given to the second reservoir (S 2 ).
Governing equations of the triple-reservoir model are F s is the total leakage of streams and trunk channels; F i is the total leakage of irrigation water including leakage under small channels; Q w is the total pumping rate of groundwater from wells; Q d is the total discharge of groundwater to springs and the lower part of Heihe River; R 1 and R 2 are recharge of groundwater under streams and irrigation districts, respectively; S g and S g0 are current storage and basic storage of groundwater; λ c and λ i are leakage coefficient of trunk channels and irrigation water, respectively; τ 1 and τ 2 are delay-times of recharge items R 1 and R 2 , respectively; a 1 , a 2 and a g are recession coefficient of reservoirs for S 1 , S 2 and S g , respectively.
It is assumed that storage of groundwater, S g , is proportional to groundwater level measured in observation wells at the site.This relationship can be described as follow where H k is the groundwater level measured in well-k, H k0 is a reference head related to the basic storage of groundwater, ϕ k is a factor to address the relationship of S g and H k .Substituting Eqs.(27-28) into Eq.( 26) we have Solution of Eq. ( 29) gives a simulated response of groundwater level to leakage of surface water.
A problem arisen in the model is the probable relationship between groundwater use and surface water use for irrigation.If the conjunctive use of surface water and groundwater is a policy in water resources management, Q w would perhaps depend on how diverted surface water can satisfy the demand of agriculture.In the last two decades, the pumping rate of groundwater in this studied area is less than 1.0×10 8 m 3 /a and gradually increases with time.However, no strong correlation is found between Q w and Q i that indicates a relatively poor integrated practice in water resources management.As a result in the model both Q w andQ i are treated as known data.

Estimation of the parameters
There are nine parameters, λ c , λ i , τ 1 , τ 2 , a 1 , a 2 , a g , ϕ k , H k0 , need to be identified in the triple-reservoir model for the studied area.Among them, λ c , λ i can be determined with observation of irrigation system.It is reported that during 1990-2005, efficiency of water transport in the trunk channels and smaller channels have been increased form 85% to 91% and form 73% to 82% respectively due to preventing of leakage.In addition, 20%-35% (23% in average) of water diverted into the farmlands by channels can penetrate through the root zone and comes to be the source of groundwater recharge.Therefore, λ c =0.15-0.09and λ i =0.44-0.37 are applied and change gradually during 1990-2005 in this study.The other parameters are identified with a simple trial-err method to make the simulation result matches the observation.
Observations of groundwater level in G5 and G79 during 1989-2006 are selected in this study for the purpose.Simulation of the process during 1986-2006 is undertaken with data of streamflow, irrigation water demand and total pumping rate of groundwater.The model is run on monthly time steps.Backward finite difference method is applied to numerically solve Eq. ( 29).
The four parameters, τ 1 , τ 2 , a 1 , a 2 , are identified according to time-lag properties of groundwater level as shown in Fig. 4. It is implied in variation of groundwater level that there are two significant recharge events every year.During July and August, groundwater level rapidly rises up to the highest value in a year.To reproduce the groundwater regime as response of two recharge events in the model, τ 1 = 0 month, τ 2 = 5 month, a 1 = 5 month −1 and a 2 = 0.25 month −1 are available.It means that the delaytime of leakage-recharge conversion under rivers is less than a month.However, the delay-time of leakage-recharge conversion under the irrigation districts is about 150 d.The peak of recharge rate increases with increasing of a 1 and a 2 .However, it is not sensitive to change of a 1 -value while a 1 > 5 month −1 .This result indicates that the first recharge event is almost the same of the leakage event under rivers in modeling of monthly behaviors.
The three parameters, a g , ϕ k , H k0 , are identified according to long-term groundwater level movement.Finally, a g = 0.39 month −1 is obtained which can lead to a best agreement between simulation and observation.The other two parameters are identified for G5 and G79 as shown in Fig. 12.
To evaluate the efficiency of the model, Nash-Sutcliffe coefficient, E k , is applied for well-k that defined as where H p k and H p kw are the model calculated and observed groundwater level, respectively, at time index of p for well-k, H ka is the mean groundwater level observed at well-k, M is the total number of data.For groundwater level during 1989-2005.The Nash-Sutcliffe coefficients of the triple-reservoir model are 0.80 and 0.65 for G5 and G79, respectively, as a result of the optimization of parameters.

Verification of the model in WATER and discussions
In 2008, the watershed scientific experiment WATER was undertaken.The observation well G5 is selected to obtain on time information of groundwater in the vicinity of Heihe River.While the new data of groundwater level measured in G5 is released in WATER, it is expected that the triplereservoir model can provide an on time simulation of groundwater response to change of streamflow.
With observed streamflow and previously estimated parameters, the triple-reservoir model gives a simulation result of groundwater level in G5 during 2006-2008, as shown in Fig. 13.This result approximately fits the observation in the period between January and September in 2008.However, it is not successful to approach the highest groundwater level in October.It is found that the triple-reservoir model briefly underestimate groundwater level in G5 during 2003-2008.This underestimation leads to a low Nash-Sutcliffe coefficients: E G5 = 0.54 and E G5 = 0.63 with respect to observation data during 2005 and 2008, respectively.
The most probable reason of the underestimation is the changes of hydrological condition after 2001 are not fully captured.In order to retrieve the ecological environment in the lower reaches area, Chinese government started a project (Dispatching Water Project) in 2001 to control export of surface water in Heihe River at the middle reaches area.This project has significantly increased the excess streamflow in Heihe River, especially after 2003.It enhances the leakage of Heihe River but reduce the leakage of irrigation water.As a response, fluctuate of groundwater level comes to be stronger which can be also found in the modeling result.However, why increasing of water table in the summer is higher than that ontained with the triple-reservoir model is still a problem.In the subsequent studies, the unaware changes of hydrological condition should be carefully investigated and captured in the model.In particular, a systematic observation of streamflow along Heiher River during Dispatching Water Project is suggested as well as measurement of river stage.It can provide essential direct data of leakage to check the leakage model and the triple-reservoir model.

Conclusions
The behavior of groundwater response to leakage of surface water in the middle reaches area of Heihe River Basin is significantly influenced by the thick vadose zone.At the studied site, it is found that the groundwater regime is a result of two recharge events due to leakage of Heihe River and irrigation water with different delay-time.
Leakage of Heihe River depends on the streamflow.A nonlinear leakage model is developed to calculate the monthly leakage of Heihe River in considering changes of streamflow, river stage and agricultural water utilization.The model is calibrated with measured data.The annual leakage of Heihe River is 15%-30% of the annual streamflow at Yingluo Gorge during 1989-2008.
Numerical modeling of variable saturated flow is carried out to investigate the general behaviors of leakage-recharge conversion through a thick vadose zone.It is found that the variable recharge from both leakage under a river and leakage under an irrigation district can be approximated by simple

Fig. 1 .
Fig. 1.A schematic representation of the studied area: (a) Heihe River Basin; (a) observation wells near Heihe River; (c) comparison between ground water level and ground surface along Heihe River.

Fig. 5 .
Fig. 5.A schematic representation of streamflow, river stage and leakage along a segment of Heihe River.

Fig. 7 .
Fig. 7.A simplified model for water flow in a thick vadose driven by leakage of a river: (a) finite-difference grid near the river (represented by a given head cell); (b) initial distribution of water content; (c) initial flow rate of water at different distance and depth; (d) redistribution of water content on the fifth day during a rising-falling event of river stage.

Fig. 8 .
Fig. 8. Change of total flow rate in the vadose zone at different depth during a rising-falling event of river stage.TM: the result of simulation with a simple tank model.
) model.During the rising-falling event of river stage, F is obtained though curve of z=1 m as shown in Fig. 8.It starts from 28.15 m 2 /d, rises up to 46.73 m 2 /d linearly with time in five days and then falls down to 28.15 m 2 /d in the next five days.The recharge for z=90 m, R, is reproduced through Eq. (20) with a time step t=0.2 d and suitable parameters.It is found that τ =0.77 d and a=2.3 d −1 lead to a best agreement between the reservoir model and the numerical model as shown in Fig. 8.This agreement confirms the validation of a simple reservoir model in recognizing the behavior of leakage-recharge conversion through a thick vadose zone under a river.

Fig. 9 .
Fig. 9. Numerical modeling of vertical transient flow in the vadose zone induce by an infiltration event: (a) finite difference grid of the model; (b) redistribution of water content at different time.

Fig. 10 .
Fig. 10.Change of groundwater recharge in the vertical model as a response to infiltration event on the top.

Fig. 11 .
Fig. 11.A schematic representation of the triple-reservoir model.S 1 and S 2 are the vadose zones under Heihe River and the Ganzhou irrigation district, respectively.