Improving AVSWAT Stream Flow Simulation by Incorporating Groundwater Recharge Prediction in the Upstream Lesti Watershed , East Java , Indonesia

The upstream Lesti watershed is one of the major watersheds of East Java in Indonesia, covering about 38093 hectares. Although there are enough water resources to meet current demands in the basin, many challenges including high spatial and temporal variability in precipitation from year to year exist. It is essential to understand how the climatic condition affects Lesti River stream flow in each sub basin. This study investigated the applicability of using the Soil and Water Assessment Tool (SWAT) with the incorporation of groundwater recharge prediction in stream flow simulation in the upstream Lesti watershed. Four observation wells in the upstream Lesti watershed were used to evaluate the seasonal and annual variations in the water level and estimate the groundwater recharge in the deep aquifer. The results show that annual water level rise was within the 2800 5700 mm range in 2007, 3900 4700 mm in 2008, 3200 5100 mm in 2009, and 2800 4600 mm in 2010. Based on the specific yield and the measured water level rise, the area-weighted groundwater predictions at the watershed outlet are 736, 820.9, 786.7, 306.4 mm in 2007, 2008, 2009, and 2010, respectively. The consistency test reveals that the R-square statistical value is greater than 0.7, and the DV (%) ranged from 32 55.3% in 2007 2010. Overall, the SWAT model performs better in the wet season flow simulation than the dry season. It is suggested that the SWAT model needs to be improved for stream flow simulation in tropical regions.


INTRODUCTION
The sustainability of water resources could be threatened by over exploitation, contamination and climate change (Principe 2012).Excessive exploitation of aquifers throughout the world has led to groundwater pollution becoming one of the major water resources management concerns to planners and managers.Groundwater is a limited resource and little has been done to map it accurately and regulate it.An understanding of the recharge processes and natural recharge rate quantification are the prerequisites for efficient and sustainable groundwater resources management (Foster 1988;Scanlon et al. 2002;Chand et al. 2005).Groundwater resources evaluation involves several factors, springs (Adelana 2010), rainfall (Liu and Zhang 1993), cli-mate change (Associates 2012) and groundwater recharge (Bredenkamp et al. 1995).When modeling the interactions in aquifers throughout time it is important to understand where and how much recharge enters the system.In arid areas it is often difficult to quantify the recharge and identify the recharge distribution (Henry 2005).Based on the available data and knowledge about the local hydrogeology, a suitable groundwater recharge method can be adopted for groundwater estimation (Kommadath 2000).Therefore, recharge quantification is needed to estimate the sustainable yield of groundwater aquifers, especially in areas like Indonesia where the amount of aquifer production usually greatly outweighs the amount of aquifer replenishment.
Base flow, a stream flow component, is usually associated with water discharged from groundwater storage (Eckhardt 2008).Knowledge about base flow is useful in assessing water quality, forecasting stream flow, and allocating water supply.As hydrological processes at the watershed scale are complex, hydrological models have been developed to simulate such processes.The Soil and Water Assessment Tool (SWAT) is one of the most commonly used hydrological models.It is a semi-distributed model (Raposo et al. 2013), but it cannot accurately reproduce groundwater hydrographs.In the SWAT model water entering a deep aquifer is not considered in the future water budget calculations but as a loss from the system.Therefore, the model focuses only on the groundwater discharge from the shallow aquifer to stream flow.In this study we estimated the groundwater recharge rates in the upstream Lesti watershed, Indonesia to further manage the groundwater resources.In order to improve sustainable water resources management in the upstream Lesti watershed we estimated the groundwater recharge in the deep aquifer using the water table fluctuation (WTF) method (Healy and Cook 2002) and improved the SWAT stream flow simulation by incorporating the groundwater recharge.
The Lesti watershed topography varies at a height between 235 and 3676 m above sea level.Based on satellite imagery and topographic map interpretation at a scale of 1:50000, the entire Lesti watershed is estimated at 59963 hectares.The watershed area is divided into 3 sub-watersheds (Fig. 1).The upstream Lesti watershed is one of the sub watersheds and it is approximately 38093 hectares.The topography is generally flat to gently rolling with a few undulating hills.The soil types are brown andosol, brown red andosol, brown red Mediterranean, brown regosol, and grey regosol.There are 9 types of land use (freshwater, resident, forest, farm, grass, irrigation field, rained field, bush, and field), and the vegetation areas are mainly farms of corn, grasses and trees.
The tropical climate is influenced mainly by monsoon winds.Four rain stations (Dampit, Poncokusumo, Turen, and Tajinan) are spread out in the upstream Lesti watershed.The rainy season is usually from November -April and the dry season is from May -October.The hydrological condition of an area is determined by several factors, such as the state of the river network, topography, soil type and climatic conditions.The Lesti upstream tributaries have 3 branches that merge into one large river, suggesting most of the land is homogeneous.Given the varying topography and homogeneous land, there is less water absorbed by the soil.Thus, floods occur frequently at particular locations.
Groundwater levels are monitored in 4 observation wells equipped with manual water level recorders (data-logger divers) in the study area (Fig. 2).Monitoring the 4 wells is part of a water resource company's (Perum Jasa Tirta 1) information services to provide decision support on integrated water resources management in the upstream Lesti watershed.The water levels were monitored monthly from 2007 -2009 for Sumber Lombok and Sumber Pakem, and from 2007 -2010 for Sumber Suko and Sumber Sari.

SWAT Model Description and Data Preparation
The SWAT is a semi-distributed, continuous model for predicting the impact of land management practices on water, sediment and agricultural chemical yields in large, complex watersheds with varying soils, land use and management conditions over long periods of time (Du et al. 2013).It simulates the hydrological processes based on weather, soil properties, topography, vegetation, and land management practices occurring in the watershed.In this study, we applied AVSWAT 2000 (ArcView SWAT) program, an extension built in Geographic Information System (GIS) ArcView 3.3 (ESRI).The program was developed by Texas Water Resources Institute, College Station, Texas, USA (Luzio et al. 2002;Neitsch et al. 2002).
For modeling purposes the AVSWAT program allows users to delineate the watershed into numbers of sub basins by specifying a critical source area.The watershed is divided into 39 subbasins, with the watershed outlet located at subbasin 39.The major required input data includes elevation (Digital Elevation Model, DEM), land use, soil and climate data.The Indonesia elevation map (1:25000) is sourced from the Coordinating Agency for Surveys and Mapping (Bakosurtanal agency).All data were digitized in CAD format with the same scale and coordinates of the original maps, and then exported as vector data for GIS program (ArcView GIS 3.3).Daily rainfall (2001 -2010), land use, soil map, stream flow data at Tawangrejeni outlet (2001 -2010), and observed data of wells including discharge and water level fluctuation (2007 -2010) are sourced from the Water Resources Company (Perum Jasa Tirta 1), East Java, Indonesia.Climatological data are sourced from Karangploso Meteorology and Geophysics, East Java, Indonesia.Although there are four climatic stations available in the watershed, only three stations (Dampit, Poncokusumo, and Tajinan stations) are used due to the delineation of sub basins that only takes climatic stations nearby.In this study a total of 9 types of land uses and 5 soil types are classified and simulated for the upstream Lesti watershed.
The SWAT model can simulate the physical processes associated with the movement of water, sediments, plant growth, and nutrient cycles (Neitsch et al. 2002).The SWAT model simulates the hydrological processes in a watershed for two main parts, namely the land and water phases of the hydrological cycles.The hydrological cycle is simulated based on the water balance equation, where precipitation (P) is partitioned into discharge (Q), evapotranspiration (E), and soil storage (ΔS), and the discharge (Q) includes surface runoff (SURQ), lateral flow (LATQ) and return flow (GW).The daily rainfall, maximum and minimum air temperature, solar radiation, wind speed and humidity data are required in this model as it affects the water balance (Neitsch et al. 2002).A watershed is divided into several sub basins and each sub basin is further grouped into various HRUs (Hydrologic Response Units), a unique combination of land use, soil and slope at a sub basin.The amount of water leaving an HRU and entering the main channel is calculated as follows, where TLOSS is transmission loss.
Water Yield (mm) = SURQ + LATQ + GW -TLOSS (1) Many studies have shown that base flow is an important component of the simulated discharge within their study areas (Chekol et al. 2007;Luo et al. 2012).In the SWAT model the lag between the time when water passes through the lowest layer of the soil profile by percolation and the time when water exits the soil profile into the shallow aquifer depends on the depth of the water table and the hydraulic properties of the geologic formations in the vadose and groundwater zones (Yang et al. 2010).The "groundwater table height" is computed to express the thickness of the shallow aquifer in terms of the depth of water stored which will contribute to the stream.The shallow or bypass flow enters and flows through the vadose zone before becoming shallow aquifer recharge.The water in shallow aquifers also revap to the overlaying soil layers by capillary pressure or direct absorption by plant roots.The shallow aquifer contributes base flow only if the amount of water stored in the shallow aquifer exceeds a threshold value specified by the user (Fig. 3).One of the most essential components of an efficient groundwater model is the accuracy of recharge rates within the input data.The conventional groundwater flow analysis performed by the extension MODFLOW program in SWAT often overlooks the accuracy of the recharge rates that are required to be calculated into the model (Kim et al. 2008).A procedure to compute perched groundwater support by DRAINMOD theory has already been used to expand SWAT's capabilities (Vazquez-Amábile and Engel 2005).In SWAT, groundwater prediction is affected by some parameters, such as ALFA_BF (base flow recession constant, an index of groundwater flow response to changes in recharge), GWQMN [threshold depth of water in the shallow aquifer required for return flow to occur (mm)], GW_DE-LAY [Groundwater delay time (days)], and GW_REVAP (the amount of water moving from the shallow aquifer into the overlying unsaturated zone).

Water Table Fluctuation Method
The WTF method is one of the most widely used techniques for estimating groundwater recharge (qR) over a wide variety of climatic conditions.Due to the availability of water level data and the simplicity of estimating recharge rates from temporal fluctuations or spatial patterns of water levels, the groundwater recharge can be calculated as follows (Hall and Risser 1993;Healy and Cook 2002;Scanlon et al. 2002).
Where S y (-) is the specific yield of the groundwater aquifer material, h (L) is the water table height and t is time (T).
The WTF method is based on the premise that a rise in groundwater level results from recharge arriving at the water table.Favorable aspects of the WTF method include its simplicity and capability for any well that taps the water table (Sharma 1989).Compared to the potential recharge given by other methods, the actual recharge estimated by the WTF method is more reliable (Ordens et al. 2012).
The water level rise (Δh) is generally computed as the difference between the peak water level rise and the value of the extrapolated antecedent recession curve at the time of the peak.The recession curve is the trace that the well hydrograph would have followed had there not been any recharge (Delin et al. 2007).There are two major approaches for estimating water level rise: the master recession curve (MRC) and the graphical extrapolation approach.The MRC for a given site is a characteristic water table recession hydrograph which represents the average behavior of a declining water table for the site and can be used to predict the decline in water table in the absence of recharge (Jie et al. 2011).We applied the graphical extrapolation method because it is simple and less time-consuming.The groundwater hydrograph of each well was plotted and the antecedent recession curves were manually extrapolated to obtain the water level rise during the recharge period.However, the graphical extrapolation approach has more subjectivity compared to the MRC, since different persons could obtain slightly different recession curves and subsequently different values for the water level rise.
The specific yield of a rock or soil is defined as the ratio of water volume that drains from a saturated geomaterial due to gravity to the total volume of the geomaterial (Meinzer 1923;Healy and Cook 2002).It can also be seen as a fraction of the porosity of an aquifer that can be drained by gravity.The value depends on the grain size, shape and distribution of pores and compaction of the strata (Gupta and Gupta 1999), lithology, temperature (Meizer 1923) and depth to water table.In theory, specific yield is treated as a storage term that does not depend on time, accounting for the instantaneous release of water from storage.In practice, however, the release of water is often not instantaneous but time-dependent (Lerner et al. 1990;Healy and Cook 2002;Nachabe 2002).This is more evident in situations when the water table declines relatively fast and the drainage from the unsaturated zone may lag behind depending on the soil properties (Storm 1991).Lerner et al. (1990) ascribed standard specific yield values to be taken from literature rather than field test measurements when values from laboratories are unavailable.We did not determine the exact specific yield values for the aquifer material, but selected the specific yield values based on the literature (Prickett 1965, cited in Lerner et al. 1990).The specific values (0.12 -0.18) of the same soil (red brown latosol) and texture (sandy loam) in India were used to represent the characteristics of the upstream Lesti watershed.
Four observation wells are located at sub basins 5, 6, and 33, spreading from upstream to downstream.It is assumed that the water level rise from these four wells can represent the groundwater recharge in all sub basins.Thus, the whole upstream Lesti watershed was divided into two parts that cover different groups of sub basins: groundwater 1 (sub basins 1 -18) and groundwater 2 (sub basins 19 -39) (Fig. 2).Separation was conducted to determine groundwater recharge from upstream to downstream by adding the recharge mechanism of the upstream (groundwater 1) and downstream (groundwater 2) to represent the recharge at the watershed outlet.The area-weighted average groundwater prediction for the entire watershed was calculated for each year.

Model Performance Evaluation
Simulated stream flow values are compared with observed data to test the accuracy.Three statistical criteria were used to test accuracy: R-square value, deviation of runoff volumes (D v ) and Nash-Sutcliffe coefficient (E NS ).The R-square value describes the degree of colinearity between the observed and predicted time series.Higher R square analysis values indicate better agreement between the observations and simulations.A perfect model has a correlation coefficient equal to 1.According to American Society of Civil Engineers (ASCE Task Committee 1993), E NS and D v , also known as the percentage bias, are simple goodness-fit criteria.Smaller D v values correspond to better model calibration, with a perfect match generating D v = 0 (MacLean 2005).The smaller the D v value the better the model performance.It should be noted that D v can take any values but in this study the smaller values of D v were satisfied (ASCE Task Committee 1993).In addition, the model performance with an E NS value greater than 0.5 is regarded as satisfactory as recommended by Santhi et al. (2001).

Streamflow Simulation
Before model calibration the simulated streamflow was first compared with the observed streamflow at the gauge Tawangrejeni (subbasin 39) (Fig. 4a).It was found that the model underestimated the streamflow in terms of the annual simulated and observed stream flow are 1937 and 2610 mm year -1 , respectively.Moreover, during the dry season the simulated stream flow is much lower than the observed stream flow.It is because the Lesti River, located in extensive and highly productive aquifers, is regarded as a perennial river along with many water resources.
Water yield simulated at HRU level in the SWAT model is the net amount of water that leaves the sub basin and contributes to stream flow in the reach during the time step.The difference between the simulated water yield sum at sub basins and the simulated stream flow at the watershed outlet is the water loss from the reach by transmission through the streambed (Fig. 4b).In addition, the simulated stream flow was examined using the water balance equation as the sum of difference between precipitation and evapotranspiration [ ( ) P E -

/
] is the watershed discharge and soil water storage during a time step (Table 1).Therefore, the difference between ( ) P E -

/
and simulated water yield is the monthly soil water storage in the watershed (Fig. 4c).The similarities among the flow simulation results, water yield and the difference between precipitation and evapotranspiration [ ( ) P E -

/
] indicate that the SWAT model can reasonably simulate the stream flow and reflect the impact of rainy and dry seasons on the stream flow.In order to improve the stream flow simulation, parameters that affect the stream flow were selected for model calibration.Suggested by Alansi et al. (2009), we selected and calibrated 4 parameters (ALPHA_BF, GWQMN, GW_DELAY, and GW_REVAP) (Table 2).However, the simulated stream flow was not much improved after calibration with little improvement in R-square about 0.4.

Groundwater Recharge Estimation
The WTF method was applied to estimate the groundwater recharge in the upstream Lesti watershed.The monitoring of water levels in Sumber Lombok and Sumber Pakem started from 2007 -2009, for Sumber Suko and Sumber Sari started from 2007 -2010 (Table 3).The highest and lowest annual area-weighted recharge estimations were found at Sumber Suko and Sumber Lombok, 820.9 and 306.4 mm in 2008 and 2010, respectively (Table 4).The groundwater recharge was then added into the simulated stream flow at the watershed outlet to improve the stream flow simulation.
The monthly rainfall for the study area was found higher during November -April (rainfall season), and lower during May -October (dry season) (Fig. 5).Although the rainfall season in the study area starts in November, the water level in all wells does not start to rise until May/April.This month lag can be described as a period of soil refilling due to moisture deficit inherited from the past dry season.The    lag shows that there is a threshold water table, soil texture, and non-linear relationship between rainfalls which fill the underground cavity space.Recharge prediction was calculated by multiplying the annual maximum recharge rate in every well that the sub basin covers.The highest and lowest predicted groundwater depths were found in 2008 and 2010 (Table 4).The monthly groundwater prediction can be assumed as a contribution from the deep aquifer (> 20 m) (Rao and Yang 2010), which SWAT does not take into account for the base flow simulation.Therefore, we averaged the annual groundwater prediction using the number of dry months in each year (Table 4) and the monthly groundwater predictions were then directly added to the simulated stream flow, resulting in an improved flow simulation that approximates the stream flow observations.For example, compared to the flow simulation before adding groundwater prediction in 2007 (Fig. 4a), the improved flow simulation was closer to the flow observation, especially during the dry season (Fig. 6).
Moreover, the improved model showed that the simulated flow during 2007 -2010 is positively correlated with the observed flow with the R 2 values ranging between 0.7 -0.81, while D v (%) with the range of 32 -55.3% and Ens values greater than 0.7 were satisfied in this study (Table 5).

DISCUSSION
Model calibration is performed through a combination of trial and error adjustments and limited parameters optimization (Arnold and Allen 1996).The calibration process requires precision and accuracy in evaluating the effect of any change in each parameter on model simulations.However, there are limits to changing groundwater parameters in the SWAT model.Peterson and Hamlett (1988) found the SWAT model was not able to simulate base flow due to the presence of soil fragipans.Moreover, the SWAT model has difficulties expressing the spatial distribution of groundwater levels and recharge (Jha 2009).The SWAT's groundwater module lumps the distributed parameters, such as hydraulic conductivity distribution, thus the groundwater distribution could not be well presented (Arnold et al. 1993).These parameters are applicable to all types of land use and soils, whereas other parameters, USLE_C, SOL_AWC, CN2, ESCO, allowed in accordance with the default program of SWAT without any change.In this study, we adjusted some parameter values in January until December in the whole year as indicated that in dry months the results are unpredictable, but the simulated stream flow still has little or no change compared to noncalibrated simulated stream flow.The underestimated stream flow is similar to the results conducted in other tropical watersheds (Reungsang et al. 2010;Phomcha et al. 2011).
The groundwater recharge is controlled by factors such as the rate and duration of precipitation, the antecedent moisture condition of the soil profile, geology, soil properties, the depth to water  and Barry 2010).The water that recharges the groundwater reservoir may further re-emerge as stream flow.In the upstream Lesti watershed, precipitation plays an important role in groundwater recharge.Critical examination of the groundwater hydrographs and water level data of the observation wells suggests that groundwater recharge in the upstream Lesti watershed is almost entirely from the seasonal rainfall, since water level rise occurred mostly in the dry season.Though there was some accumulation of recharge in the rain season possibly due to the regional groundwater, the amount is very small.Therefore, it can be reasonably concluded that, the contributions from aquifers outside the study area to groundwater recharge in the upstream Lesti watershed is insignificant.
The incorporation of groundwater recharge estimates by this study into SWAT stream flow simulations greatly improved the stream flow estimation, especially in dry season.The result revealed that base flow contributed by deep aquifers plays a big part in the water balance system.The base flow is usually greater than surface flow in some areas.Sathian and Syamala (2010) simulated water balance components of the basin and found that the base flow contribution is the major part in the stream flow, followed by lateral flow and surface runoff.Chekol et al. (2007) found that surface runoff and base flow contributed 16.8 and 53%, respectively of the stream flow.
The possible reasons for the discrepancy between the simulated and observed currents can be discussed a number of ways.For example, the overestimated flow in dry season in 2009 (Fig. 6) may be due to the lack of water usage information.Seventy percent of the area is dominated by agriculture powered by small scale irrigation systems.However, information about the water used by the system, such as pumping for irrigation and water diverted from the river is not available.
Percolation to the deep aquifer was analyzed to prove the assumption that the rising water is released from soil as discharge to the river.Percolation has a relatively linear relationship with rainfall in the rainy season, while percolation in the dry season approaches zero (Fig. 7).The annual percolation is estimated using the area-weighted of average water level rise in different groundwater zones.The small percolation value equals the value without the addition of groundwater.After the addition of the value of groundwater in the dry months in 2007 until 2010 (Table 6), simulated stream flow approached observed stream flow (Fig. 7).Therefore, groundwater is very important for predicting how much water from the aquifer supplies the river.

CONCLUSION
We investigated the impact of groundwater recharge on SWAT streamflow simulation, which did not perfom well during the dry season.Due to the inability of the SWAT model to take into acount the groundwater discharge from the deep aquifers, manual calibration was needed to improve groundwater recharge estimation in the upstream Lesti watershed.Moreover, quauntification of groundwater recharge from a deep aquifer using the WTF method was prooven to be useful in improving SWAT simulated streamflow at the watershed outlet (Tawangrejeni).The fluctuations in the water table were analyzed for better understanding of the change in water level rise in the watershed.However, the limitations of this study are the data avalibility and accuracy.Therefore, it is suggested that related government agencies improve the completeness of the inventory data [i.e., installation of Outlet Automatic Water Level Recorder (AWLR)] for future studies in sustainable water resources management.

Fig. 2 .
Fig. 2. Soil distribution and locations of the rain gauges and observation wells in the upstream Lesti watershed.

Fig. 4 .
Fig. 4. Comparison of observed stream flow, simulated stream flow, water yield and watershed discharge (a) simulated and observed flow (b) simulated flow and water yield (c) simulated water yield and watershed discharge.

Table 1 .
Description of value output in this research.

Table 5 .
table and aquifer properties, vegetation and land use, topography and landform (Obuobie The model evaluation for the goodness of stream flow.