Development and application of a conceptual hydrologic model to predict soil salinity within modern Tunisian oases

Journal of Hydrology 380 (2010) 45–61

Contents lists available at ScienceDirect

Journal of Hydrology journal homepage: www.elsevier.com/locate/jhydrol

Development and application of a conceptual hydrologic model to predict soil salinity within modern Tunisian oases Brahim Askri a,*, Rachida Bouhlila b, Jean Olivier Job c a

Ecole Nationale d’Ingénieurs de Gabes, 6029 Gabes, Tunisia Ecole Nationale d’Ingénieurs de Tunis, B.P. 37 Le Belvédère, 1002 Tunis, Tunisia c Centre Régional de l’Eau et de l’Environnement, ESIB, Université St. Joseph, BP 1514, Beyrouth, Lebanon b

a r t i c l e

i n f o

Article history: Received 30 April 2009 Received in revised form 7 October 2009 Accepted 14 October 2009 This manuscript was handled by L. Charlet, Editor-in-Chief, with the assistance of Eddy Y. Zeng, Associate Editor Keywords: Oasis Soil Water table Salinity Conceptual model

s u m m a r y In modern oases situated in the south of Tunisia, secondary salination of irrigated lands is a crucial problem. The visible salt deposits and soil salination processes are the consequence of several factors including the excessive use of saline irrigation water, seepage from earthen canal systems, inefficient irrigation practices and inadequate drainage. Understanding the mechanism of the secondary salination is of interest in order to maintain existing oases, and thus ensure the sustainability of date production in this part of the country. Therefore, a conceptual, daily, semi-distributed hydrologic model (OASIS_MOD) was developed to analyse the impact of irrigation management on the water table fluctuation, soil salinity and drain discharge, and to evaluate measures to control salinity within an oasis ecosystem. The basic processes incorporated in the model are irrigation, infiltration, percolation to the shallow groundwater, soil evaporation, crop transpiration, groundwater flow, capillary rise flux, and drain discharge. OASIS_MOD was tested with data collected in a parcel of farmland situated in the Segdoud oasis, in the south-west of Tunisia. The calibration results showed that groundwater levels were simulated with acceptable accuracy, since the differences between the simulated and measured values are less than 0.22 m. However, the model under-predicted some water table peaks when irrigation occurs due to inconsistencies in the irrigation water data. The validation results showed that deviations between observed and simulated groundwater levels have increased to about 0.5 m due to under-estimation of groundwater inflow from an upstream palm plantation. A long-term simulation scenario revealed that the soil salinity and groundwater level have three types of variability in time: a daily variability due to irrigation practices, seasonal fluctuation due to climatic conditions and annual variability explained by the increase in cultivated areas. The irrigation interval was found to be important with irrigating once each ten days leading to soil salinity increase during the dry summer season and to a rising water table during the autumn–winter period. The annual increase in the irrigated area caused a decrease of the irrigation water depths, and thus an augmentation of the soil and groundwater salinities. The surface area affected by a soil salinity concentration above 15 g/L has increased from 2% of the study parcel area in June 1992 to about 50% four years later due to the abandonment of several cultivated basins. Ó 2009 Elsevier B.V. All rights reserved.

Introduction The southern domain of Tunisia extends from the low plains situated south of Gafsa Mountains up to the Sahara boundaries (Mtimet, 1983). In this arid zone, oasis agriculture using spring water has been performed since ancient times (Lasram, 1990). The traditional oasis is an ecosystem with particular properties (Sellami and

* Corresponding author. Tel.: +216 97 278 492; fax: +216 75 392 190. E-mail address: [email protected] (B. Askri). 0022-1694/$ - see front matter Ó 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2009.10.022

Sifaoui, 2008). Typically, several crops such date palms, fruit trees and market gardening are supported on the same field. The cultivation of these three levels together makes a microclimate, commonly called ‘‘oasis microclimate”. For many centuries, oases were used by authorities to stabilise indigenous semi-nomadic people in areas where water is scarce (Kassah, 1996). A change of the agricultural policies in Tunisia took place in the 1960s creating

46

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

an intensive modern agriculture in the southern part of the country, in order to maximise the production and the exportation of dates and to promote economic development (Sghaier, 1984). The exploitation of two deep groundwater aquifers and the use of water for irrigation have allowed the creation of a new type of oases with new irrigated schemas. Irrigation water allocation, equipment renovation and infrastructure implementation were sponsored by the state to encourage farm development. However, the dry climate and the lack of irrigation knowledge by young farmers induced an excessive use of saline irrigation water. As soils are predominantly sandy and drainage is almost non-existent, the intense irrigation has led to the creation of a perched and saline water table. The shallow groundwater has contributed to salt accumulation through capillary rise and has resulted in the accumulation of salts at the soil surface. Degradation of soil had strong negative impacts on the growth and vigour of palms, resulting in a decrease of date yield and deterioration of date quality. As a solution, since 1990, the Tunisian government has started encouraging farmers to replace the earthen irrigation canals using PVC pipes, subsidising 40–60% of the total investment. This policy focused on the importance to water saving. Therefore, in 2001, a project called ‘‘Amelioration of Irrigated Areas in South Oases (APIOS)” was initiated to install concrete canals for irrigation and subsurface drainage system, the objectives being to improve irrigation frequency, to increase water productivity and to save water resources since they are limited (Ben Mohamed, 2002). In order to prevent continuing soil salination within modern Tunisian oases and thus degradation of date production, it is necessary to clarify the mechanism of secondary salination in relation to water management including irrigation and drainage. Numerical models play an important role in the studies of land salination in areas with irrigation and drainage system (Xu and Shao, 2002). In particular, they help to understand hydrologic processes and to derive strategies for salinity control. Over the last two decades, two different approaches have emerged to analyse water flow and solute transport in structured soils within irrigated lands. One consists of a physically based or mechanistic approach based on Richard’s equation for the movement of water in unsaturated soil, in combination with a differential salinity dispersion equation. Typical examples are the SWAP (Singh, 2004; Ninghu et al., 2005) and HYDRUS-1 models (Tedeschi and Menenti, 2002; Rassam and Cook, 2002). Most of these models solve the Richards equation in 1-D (vertical direction) and regionalisation is carried out based on land use, soil, or topography, or two or three of these combined (Renschler et al., 2001). This approach can become limiting quickly because of the over-parameterisation and the lack of readily available data. In a review study, Moran et al. (2004) reported that these models are generally of a limited practical use because of the difficulties of specifying parameters and the initial and boundary conditions. The other approach uses capacity storage-type models called ‘‘Simple models” which attempt to characterise the average response of each cell rather than to capture in detail the physics and the variations that may occur within the cell (Bouraoui et al., 1997). Initially, this approach was used to develop hydrologic models that relate runoff to rainfall at catchment or sub-catchment scales (Vardavas, 1988). Since the last decade, several capacity models were extended to understand hydrologic processes occurring within agricultural lands in relation to water management including irrigation and drainage (Connell et al., 1999; Lamsal et al., 1999; Singh et al., 1999). A typical example is the SALTMOD model, which uses seasonal water balance components as input data (Srinivasulu et al., 2004; Bahceci et al., 2006). This model characterises the agricultural land in a lumped manner, assuming uniform distribution of the cropping, irrigation and drainage characteristics over the entire area. This assumption is not realistic for oasis ecosystems because of the

large-scale spatial variability in crop levels, drainage condition, and water availability and management. Roost et al. (2008) developed the OASIS model which is a planning tool for medium to large-scale canal irrigation systems. The model relies in a conceptual representation based on three levels: system, irrigation unit and field. It integrates in a lumped manner surface and groundwater flows and contains a crop growth module. The actual conceptualisation of the OASIS model does not include processes of salt transport. Thus, it is not useful to simulate scenarios of land salinisation. This paper describes the development of a conceptual, continuous, daily hydrologic model called OASIS_MOD and its usefulness in investigating the mechanism of secondary salination within Tunisian oases. Because of the high heterogeneity of water management within these ecosystems, OASIS_MOD is a spatially semi-distributed model dividing the field in homogeneous volume units. It has been developed to highlight the possibility of improving irrigation management on the basis of long-term predictions of the impact of irrigation practices on water table fluctuation, soil salinity and drain discharge. The lack of available input data at the oasis scale has required the application of the model to a 1.5 ha parcel of farmland situated in the Segdoud oasis in the southwest of Tunisia. We describe the calibration and the validation of the model and then performed a simulation scenario to understand the mechanism of the secondary salination occurring in an oasis ecosystem.

Water management within Tunisian oases The main source of irrigation water in the south of Tunisia is the Northwest Sahara Aquifer System (NWSAS). It consists of two main aquifers: the complex terminal (CT) and the continental intercalary (CI) (Horchani, 2007). Apart from these two aquifers, water can also be drawn from the shallow groundwater. This aquifer was often observed in irrigated areas in the south of Tunisia due to the hydraulic anisotropy caused by gypsum born crust which could be present at medium depths of 40–60 cm (Mtimet, 1983). Irrigation management within Tunisian oases is based on longstanding traditions and water allocation. Irrigation interval and time are based on historic water rights (Belloumi, 2006). The irrigation network operates on rotation delivery, and irrigation duration is imposed on each parcel of farmland, depending on the served area. Water allocation is defined by an organisation charged with water management called AIC (Association of Collective Interest). An operator, who is responsible for opening and closing the irrigation outlets, handles the distribution of water. The only guidance given to the operator is the time required per hectare and the sequence of gate opening. Irrigation water distribution to parcels is unequal due to several factors including the socio-economic situation of farms, the parcel elevation, the soil salinity and water table depth. Modern oases are divided into 1 to 2 ha individually parcel of farmland. Each parcel is a rectangular area of land, typically 150 m long and 100 m wide. The plant density is about 100 palms/ha, and the spacing between them is about 10 m (Fig. 1). The parcel is subdivided into square basins of 10–20 m2 called huwdh, limited by ridges that surround one palm tree to retain ponded water on the soil surface. Lands around basins are uncropped, and farmers increase annually the basin area into the surrounding area in response to the development of the palm roots. The irrigation system in the parcel consists of an irrigation outlet, an open concrete canal and several earthen canals called suqyyus (Fig. 1). The flood irrigation is the main irrigation method practiced within oasis ecosystem. This method induces seepage to the groundwater and a high spatial variability of the irrigation

47

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

inflow from the head of the suqyyu to the tail. Experiments conducted in several modern oases showed that about 30% of the amount of irrigation water was lost to infiltration in a 100 m stretch of earthen canal (SANYU Consulting, 1996). In modern oases, the intense irrigation and the poor drainage condition have led to a rising water table, mainly during the monsoon months (January, February and March). Drainage was therefore required to ensure that the high saline water table does not encroach into the root zone being initially accomplished by open-ditches. This system was later abandoned due to problems related to sediment deposition and vegetation development. Within the framework of the APIOS project, a subsurface drainage system was implemented in several oases during 1991–1997.

Model formulation

Irrigation outlet

Uncropped area

Saqyya

basins

Irrigation concrete canal

Fig. 1. Topology of irrigation system in the parcel level within Tunisia modern oases.

Conceptualisation and modeling approach Tunisian oases are typically highly heterogeneous, with several crop levels, various soil types, a range of irrigation schemas and different drainage systems (Mtimet and Hachicha, 1995). The oasis ecosystem is divided into parcels of farmland, and each parcel is connected to a network of irrigation and drainage canals. Each parcel is subdivided in square surface elements of uniform parameters, called cultivated units CUs (Fig. 2). The number of CUs depends on the surface properties, the number of cultivated basins and availability of data related to the soil and groundwater properties. Within a CU, the soil properties are assumed to be homogeneous, and the groundwater level is considered uniform. CUs thus define the actual spatial resolution of the model. As shown in Fig. 2, the unsaturated zone within a single CU has two parts, the irrigated land inside basins and the non-irrigated land in the surrounding area. The shallow groundwater provides the pathway for the movement of excess water from the irrigated land to the non-irrigated land (Konukcu et al., 2006). Because of the high infiltrating properties of soils within Tunisian oases (Mtimet, 1983), it is assumed that: – there is no lateral subsurface flow from one CU to another in the unsaturated zone; – the unsaturated zone is vertically homogeneous; – the water movement among CUs occurs only in the shallow groundwater; – the solute transport of salts is considered to be an advective transport of conservative matter; – a complete mixing of the irrigation water with the soil solution. A simple schematisation of a single CU is the reservoir model as illustrated in Fig. 3. The first reservoir corresponds to the irrigated land inside basins. Here, the input variables are evapotranspiration, irrigation and rainfall. The second reservoir is the nonirrigated land in the surrounding area. Rainfall and evapotranspiration are the main driving forces for this reservoir. The third reservoir is the shallow groundwater. The last reservoir is the lower boundary of the model and represents the observed groundwater level. The basic processes incorporated in the model at the CU scale are irrigation, rainfall, infiltration, percolation to groundwater, evapotranspiration, groundwater flow, capillary rise and drainage discharge. Fig. 3 illustrates these components within a single CU during a specified period of time, which is usually selected as 1 day. Two hydrologic phases operate sequentially in the OASIS_MOD model. Initially, there is a slow daily infiltration phase and then an even slower evaporation phase. These two phases occur each day, but at vastly different rates (Boughton, 1968).

During the infiltration phase, the movement of water is assumed to occur only in the liquid phase. The irrigated land is supplied with water by irrigation, rainfall and capillary rise fluxes. The last two processes are the main source of water for the non-irrigated land. When the soil water content reaches field capacity (defined at 5 kPa), the water on excess drains through the bottom of the soil zone to become potential recharge to the shallow groundwater. The groundwater aquifer receives percolation from the irrigated and non-irrigated lands. It can contribute via capillary fluxes to soil evaporation and crop transpiration. In Tunisian oases, the capillary rise is the main cause of salt accumulation in the root zone whenever the water table is close to the ground surface (Mtimet and Hachicha, 1995). The capillary rise induces the soil water movement in the liquid phase. Thus, it was assumed to occur during the infiltration phase. The model simulates two groundwater flow processes in addition to capillary rise. The first represents the effect of a drainage system (ditches or pipe drainage). The second aims to capture the impact of hydraulic heads on groundwater flow among CUs. During the daily evaporation phase, crop transpiration and soil evaporation dry out water stored within the irrigated and non-irrigated lands. These two processes induce the transfer of water in the form of vapour from the land surface to the atmosphere. Infiltration phase The unsaturated zone In Tunisian oases, soils are predominantly sandy with high infiltration capacity, and ponding will be of limited duration. The infiltration rate, F (m/day), is computed as follows:

 F¼

IR þ Ra for the irrigated land Ra

for the non-irrigated land

ð1Þ

where IR is the irrigation rate (m/day), Ra is the rainfall intensity (m/day). The percolation rate through the bottom of the irrigated land, Pe (m/day), is computed according to (modified Vardavas, 1988):

( Pe ¼

P max



W irr W FC W SAT W FC

0



zS hSUB



if

W irr  W FC and zS  0

ð2Þ

otherwise

where Pmax is the maximum percolation rate (m/day); zS is the depth of the water table (m); hSUB is the substratum depth of the shallow groundwater (m); WSAT and WFC are the soil water storage at saturation and field capacity (m); Wirr is the actual water storage in the irrigated land (m), which can be written as follows:

W irr ¼ hirr  zS

ð3Þ

48

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

(2006) propose the WAVE model to simulate the soil water fluxes at 1.0 m depth for different groundwater depths and evapotranspiration rates. The obtained results showed that upward flux is less influenced by evapotranspiration rates but is highly dependent upon the groundwater depth. Empirical and semi-empirical approaches have been used to estimate the upward flux from the bottom of the root zone as a term of the soil balance (Jorenush and Sepaskhah, 2003; Raes and Deproost, 2003; Konukcu et al., 2006). The methodology proposed by Doorenbos and Pruitt (1977) to estimate the rate of this process, CR (m/day), is given by:

(a)

Crop transpiration

8 if W  W WP CR > >  < max  W p W if W WP 6 W  W p CR ¼ CRmax W p W WP > > : 0 if W  W p

Rainfall

Basin

Evaporation Irrigation

(b) Capillary rise Lateral inflow

where W is the actual soil water storage in the root zone (m), WWP is the soil–water storage at permanent wilting point (m), Wp the root zone soil water storage corresponding to the depletion fraction for no stress, p; and CRmax is the potential rate of the capillary flux (m/day). CRmax depends upon the soil characteristics (Doorenbos and Pruitt, 1977), and the depletion fraction p is crop specific (Allen et al., 1998). The approach described by Eq. (4) is relatively easy to implement but could be affected by inaccuracies in estimating CRmax, since it depends not only upon the soil characteristics but also on the groundwater depth and evapotranspiration demand (Liu et al., 2006). In addition, this approach is based on the assumption that the upward flux occurs only when the soil water stored in the root zone is lower than the storage at the non-stress threshold (W  Wp). However, this flux occurs also when W  Wp if the groundwater table is high, the soil moisture is relatively low and the crop climatic demand is high. In MODFLOW model, the upward flux from a static water table is assumed negligible when the soil water content in the topsoil is above field capacity (Raes and Deproost, 2003). In the present model, Eq. (4) was modified as follows:

Infiltration Lateral outflow Percolation

Shallow groundwater Substratm

Drainage

Fig. 2. Schematic representation of OASIS_MOD model structure. (a) Discretisation of the model domain into cultivated units CU. (b) Processes simulated in a single CU.

where hirr is the mean volumetric water content in the irrigated land (m3/m3). In a similar way, the percolation rate through the bottom of the non-irrigated land, Punc (m/day), is calculated by Using Eq. (2). In this case, the actual water storage in this reservoir, Wunc (m), is used instead of Wirr. Computations of Pe and Punc with Eq. (2) assume that the percolation process ceases when either the soil moisture is lower than the moisture at field capacity or the water table rises to the ground surface. Capillary rise is defined as the upward flux leaving a static water table due to crop transpiration and soil evaporation (Jorenush and Sepaskhah, 2003). The rate of upward flux depends on several factors including textural and structural soil characteristics, evapotranspiration demand, available soil water storage and groundwater depth (Mermoud and Seytoux, 1989). Liu et al. ETirr

8 CR if W irr  W WP >  < max W FC W irr if W WP 6 W irr  W FC CR ¼ CRmax W FC W WP > : 0 if W irr  W FC

ETunc

Ra

Irrigated land

Non-irrigated land

Upstream CU

Soil zone

CR

ð5Þ

A steady-state condition is assumed to occur in the soil profile when the actual soil water storage is lower than the storage at permanent

IR

Downstream CU

ð4Þ

Pe

CR

Shallow groundwater

Groundwater inflow

Drainage outflow Groundwater ouflow CU boundary Fig. 3. Flow diagram of OASIS_MOD model.

49

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

wilting point. Under this condition, the upward fluxes are assumed to occur at the potential rate (CRmax), which may be approximated by the approach proposed by Gardner (1958). This approach proposed the following empirical relationship between the hydraulic conductivity, K, and the matric potential, wm, that has been largely used (Gardner, 1958):

 Kðwm Þ ¼

a b þ wnm

 ð6Þ

where a (mn+1  s1), b (mn) and n () are parameters determined in the laboratory experiments. Replenishing this expression in the generalised Darcy equation and integrating to infinity lead to the following analytical solution (Gardner, 1958):

CRmax ¼

b znS

ð7Þ

where b (mn+1  s1) and n (the same as in Eq. (6)) depend on the soil under consideration. In the approach of Gardner (1958), b = 3.77a for n = 3/2, b = 2.46a for n = 2, b = 1.76a for n = 3 and b = 1.52a for n = 4. In the OASIS_MOD model, the potential rate of the upward flux described by Eq. (7) is incorporated in Eq. (5). The continuity equation for water stored in the irrigated land is expressed as follows:

dW irr ¼ IR þ Ra þ CR  Pe ¼ L dt

ð8Þ

where L represents the time rate of change of the irrigated land storage during the infiltration phase (m/day), dt is the model time step (day). Eq. (8) does not include the evapotranspiration rate, because this process is assumed to occur during the evaporation phase. In the general case, there is no analytical solution satisfying Eq. (8) and, therefore, a numerical method is employed to obtain an approximate solution. Based on Euler-trapezoidal integration over time step Dt, Eq. (8) can be re-written for day i in the form:

W irr;i ¼ W irr;i1 þ

ðLi þ Li1 ÞDt 2

ð9Þ

with Dt = 1 day. The continuity equation for water stored in the non-irrigated land is calculated as follows:

dW unc ¼ Ra þ CR  Punc ¼ B dt

ð10Þ

where Wunc is the non-irrigated land water store (m), B represents the time rate of change of the irrigated land water store (m/day) and so:

W unc;i ¼ W unc;i1 þ

ðBi þ Bi1 ÞDt 2

ð11Þ

The saturated zone Several Tunisian oases are close to sabkhas and chotts that are ephemeral salt lakes. In these ecosystems, where groundwater surface tends to follow the ground surface, topographic gradients cause significant horizontal movements of the shallow groundwater (Kadri and Van Ranst, 2002). In OASIS_MOD model, it is assumed that the groundwater flow occurs from higher topographic area to lower topographic area. Furthermore, the groundwater flow is assumed bidirectional occurring from the north to the south and from the east to the west of the model domain. Fig. 4 demonstrates the interaction between CUs, focusing on the lateral groundwater flow among CUs. The groundwater flow out of one CU becomes the groundwater inflow to the adjacent downstream CUs located in the south and in

the west of the given CU. The groundwater flow among CUs is modeled by a simple linear reservoir. This approach assumed a linear relationship between the groundwater flow and the height of water table above a user defined, reference level for downstream groundwater flow. This reference level is assumed invariable in time and may correspond to the groundwater level before the oasis creation. In OASIS model, Roost et al. (2008) applied the following equation to estimate groundwater flow among irrigation units:

qdowns ¼



cðhGW  hdowns Þ if hGW  hdowns 0

otherwise

ð12Þ

where qdowns is the groundwater flow (m/day), hGW is the groundwater level (m) (hGW is equal to zero at the sea level and is defined positively upwards), hdowns is the reference level for downstream groundwater flow (m), c is an input parameter to be calibrated (day1). The approach described by Eq. (12) is adopted in the OASIS_ MOD model to estimate groundwater flow. When the groundwater levels in the upstream areas outside the active model domain are known, the groundwater inflow to this domain can be estimated using Eq. (12). If these groundwater levels data are not available, the groundwater inflow to the model domain can be considered as an input parameter to be calibrated. Drainage discharge within a CU occurs whenever the groundwater level rises above the topographic level of subsurface drainage system (ditch or pipe drains). The groundwater flow to the drainage system, qdrain (m/day), may be expressed as linearly increasing with the height of water table above the topographic level of the drainage system as follows (Roost et al., 2008):

qdrain ¼



kðhGW  hdrain Þ if hGW  hdrain 0

otherwise

ð13Þ

where hdrain is the topographic level of the subsurface drainage system (m), k is the inverse of the drainage resistance (day1). The OASIS_MOD model uses the continuity equation for water stored in the shallow groundwater within each CU. The groundwater level within a CU varies according to: dhGW dt  ¼

A Atot



    P e þ 1  AAtot  P unc þ qNups þ qEups  2  qdowns  CR  qdrain

l

¼D

ð14Þ

where D represents the time rate of change of the groundwater level (m/day), A* and Atot are the areas of the basin and of the CU, respectively, (m2), qNups and qEups are the groundwater inflows from adjacent CUs located, respectively, in the north and in the east of the given CU (m/day), the last two flows are calculated by using Eq. (12), l (a constant for homogeneous aquifers) is the specific yield of the shallow groundwater, and so:

hGW;i ¼ hGW;i1 þ

ðDi þ Di1 ÞDt 2

ð15Þ

The Gauss–Seidel method is used to resolve Eqs. 9, 11 and 15. Salt balances within the cultivated unit The water balance for each reservoir of the CU is extended with salt balance, taking account of the amount of water and solute. The salinity concentration of the outflow (either from one reservoir into the other or by subsurface drainage) is computed by multiplying this outflow by the salinity concentration in the incoming reservoir solution. Thus, the salt balance equation for the irrigated land within a CU is taken as follows:

50

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

Drainage system of CU

qdrain

Shallow groundwater system of downstream CU

hdrain

hGW

qdowns

qups hdowns

Shallow groundwater

hsub Substratum

N

Upstream area

Downstream area

Groundwater flow direction Fig. 4. Shallow groundwater flow modelling (adapted from Roost et al., 2008).

dðW irr  C irr Þ ¼ IR  C IR þ CR  C GW  P e  C irr dt

ð16Þ

where Cirr and CGW are the salinity concentrations in the irrigated land and groundwater solutions, respectively, during the infiltration phase, CIR is the salinity concentration in the irrigation water (g/L). The combination of Eqs. (8) and (16) leads to the following equation:

dC irr C IR  IR  C irr  Pe þ ðC GW  C irr Þ  CR ¼ W irr dt

C GW;i ¼

Gi C GW;i1 þ Dt

h

A Atot



ð17Þ

In a similar way, the salinity concentration in the non-irrigated land solution, Cunc (g/L), can be expressed as follows:

C unc;i ¼

W unc;i  C unc;i1 þ Dt  C GW;i  RCi W unc;i þ DtðRai þ RCi Þ

The salinity concentration in the groundwater solution, CGW (g/L), is computed as follows:

  i   Pe;i  C irr;i þ 1  AAtot  C unc;i  Punc;i þ C upsN ;i  qupsN ;i þ C upsE ;i  qupsE ;i Gi þ DtðPe ;i þ qupsN ;i þ qupsE ;i Þ

The left-hand side of the above equation represents the time rate of change of the salinity concentration in the irrigated land solution during the infiltration phase. Using the classical differential method with implicit schema, this equation can be written for a given day i as follows:

ð19Þ

ð20Þ

where C Eups andC Nups are the salinity concentration in the groundwater solution (g/L) within CUs located, respectively, in the east and in the north of the given CU. The Gauss–Seidel method is used to resolve Eqs. (18)–(20). Evaporation phase

C irr;i ¼

W irr;i  C irr;i1 þ DtðC IR;i  IRi þ C GW;i  RCi Þ W irr;i þ DtðIRi þ Rai þ RCi Þ

ð18Þ

During the evaporation phase, water stored in the irrigated land Wirr and non-irrigated land Wunc are taken to be depleted by crop

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

transpiration and soil evaporation. These two processes are integrated in a single variable, which is estimated using the simplified version of the Penman–Monteith (FAO56) approach (Allen et al., 1998). FAO56 calculates a reference crop evapotranspiration ET0 (m) for a hypothetical grass cover. To calculate the potential evapotranspiration rate for other surfaces, the ET0 is multiplied by a coefficient (Kc) that is crop specific.

ET c ¼ K c  ET 0

ð21Þ

where ETc is the potential crop evapotranspiration under standard conditions (m), Kc () is the crop coefficient. When the soil is sufficiently wet, the soil supplies water fast enough to meet the atmospheric demand of the crop. Under this condition, the water uptake occurs at its potential rate ETc that is controlled by the atmospheric demand (Allen et al., 1998). In oasis ecosystem, this stage occurs few days after irrigation or rainfall events. As the soil water content decreases, the crop comes under stress and the transpiration rate decreases. Transpiration ceases at the permanent wilting point (Allen et al., 1998). Eilers et al. (2007) consider that moisture can be extracted by the roots without water stress when the actual soil moisture is higher than the moisture at field capacity. In the OASIS_MOD model, the crop transpiration is assumed to occur at the potential rate ETc when the actual soil water storage is higher than the storage at field capacity. Under this condition, the depth of soil water which can be used by the crop (Wirr–WWP) is higher than the potential crop evapotranspiration. If the actual soil water storage is lower than the storage at field capacity, the actual evapotranspiration is assumed to decrease linearly in proportion to the soil water storage in excess to the storage at permanent wilting point (W–WWP) when W  WWP. The actual evapotranspiration from the irrigated land, ETirr,i (m), is computed as follows (Wang et al., 2007; Rodriguez et al., 2008):

( ET irr;i ¼

ET c;i



W W WP ET c;i Wirr;i FC W WP



if W irr;i  W FC if W WP 6 W irr;i  W FC

ð22Þ

The water stored in the irrigated land after the daily evaporation phase, W 0irr;i (m), is thus given by:

W 0irr;i

8 W irr;i  ET c;i > > <   ET ¼ W irr;i  exp W FC Wc;iWP > > : W WP

if W irr;i  W FC if W WP  W irr;i  W FC

ð23Þ

if W irr;i 6 W WP

The soil evaporation and crop transpiration may induce the loss of water within the irrigated and non-irrigated lands. However, the amounts of salts in these compartments are assumed invariant after the evaporation phase. Thus, the salinity concentration in the irrigated land solution is written as follows:

C 0irr;i ¼

W 0irr;i  C irr;i W irr;i

ð24Þ

where C 0irr;i is the salinity concentration in the irrigated land solution after the daily evaporation phase (g/L). Summary of the model results Operation of the model requires basic input data for each CU, which consists of (a) meteorological data including the daily values of reference evapotranspiration and rainfall; (b) irrigation data including discharge, duration and salinity concentration; (c) initial water storages and salinity concentrations in various reservoirs of the CU and (d) groundwater inflow from upstream areas and its salinity concentration. The model output can be classified into three main groups: (i) water stored in various reservoirs of the CU including groundwater level; (ii) flow processes including infil-

51

tration, percolation to the shallow groundwater, capillary rise flux, drain discharge, groundwater flow and actual evapotranspiration and (iii) salinity concentrations in solutions within various reservoirs of the CU and the salinity of the drain discharge solution. The model requires both morphological features and physical properties. The morphological parameters are listed in Table 1 and encompass for each CU: the CU and the basin areas, the substratum depth of the groundwater aquifer and the reference levels for downstream groundwater flow and drainage system. The basin areas can be obtained from field survey. The substratum depth of the shallow groundwater can be obtained using field measurements or deduced from geological maps. The reference level for downstream groundwater flow can be estimated from the historic of the groundwater level. The physical proprieties of the soil and the shallow groundwater are listed in Table 2; most entries can be deduced from physical considerations, either literature reviews or field and laboratory measurements (Rodriguez et al., 2008).

Materials and methods Site description Data used for evaluating OASIS_MOD model have been collected in a modern oasis located in the Segdoud region in the southwest of Tunisia (geographic coordinates are latitude 34°140 N and longitude 8°100 E) at about 30 km, north of the town of Tozeur (Fig. 5). The Segdoud region belongs hydrographically to Chott el Gharsa basin and geographically includes three units: isolated palm plantations, hydrologic discharges area (Chott el Gharsa) and steppe. Geological formations are mainly Quaternary deposits. The region shows an inclined landscape with the lower altitude at the center of the Chott (20 m below mean sea level) but limited to the north by the Negueb and Chouabine mountains (700 m). Except for palm plantations, the natural steppe vegetation is very ephemeral and scarce, and for most of the year, the soil surface is a sandy and stony desert (Brunel et al., 2006). The groundwater flow occurs from the north to the south, and the predominant hydraulic gradients are in the range 0.7–1.5%, with an average of about 1.1% (Askri, 2002). The general soil type is alluvial sand with an average content above 70%. The soil infiltration rate at saturation ranges from 0.1 to 7.2 m/day, with an average of about 2.2 m/day. Organic matter content is of about 2% in the superficial soil horizon (SANYU Consulting, 1996). Weather in the Segdoud region is hot and dry in the summer and cold in the winter. The annual mean temperature is about 21.6 °C, while the absolute extreme temperatures are 44 °C in August and 2.5 °C in January (SANYU Consulting, 1996). Annual rainfall is about 96 mm, while the pan evaporation is about 2500 mm. The reference evapotranspiration measured inside Tozeur oasis is about 1643 mm per year (Richter, 1986). The Segdoud oasis covers about 166 ha and is divided into 1.5 ha parcels of farmland. The main crop is palm tree. Some fruit trees, such as pome-granate, apricot and fig, are cultivated inside the earthen canal systems. Occasionally, farms use basins located close to the irrigation outlets to cultivate market gardening crops to satisfy their food needs. Irrigation started in the oasis in July 1989, and the available cultivated area had increased from 10 ha at that time to about 20 ha in May 1995. The oasis irrigation water is pumped up from a depth of about 1000 m in the CI aquifer using two wells, which provide a total inflow of about 76 L/s. The model calibration and validation are carried out at the parcel scale, because problems related to irrigation management in the oasis level are very complex with water distribution to parcels being unequal. In May 1995, the oasis included about 10,000 individual basins or huwdhs. Thus, it was very difficult to know with

52

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

accuracy the volume of irrigation water applied at the basin scale on daily basis. The test parcel selected for this study is located at the northeast edge of the oasis at a high elevation, in order to reduce the influence of surrounding parcels on the water and salt dynamics. In 1995, the test parcel received a subsurface drainage system, consisting of three pipes 250 m long, buried at a depth of 1.5 m and spaced about 70 m apart (Fig. 6).

Available data The rainfall, pan evaporation, wind speed and temperature data are obtained at a daily time step from a meteorological station, located in the town of Tozeur. The reference evapotranspiration is calculated using the Penman–Monteith equation. In 1995, about 0.2 ha area was irrigated within the test parcel. Dates and quantities of water applications were reported for each irrigation event from 19th May to 30th September, 1995. Mean value of the irrigation water salinity concentration is of about 4 g/L. OASIS_MOD model is calibrated and validated against the groundwater level

Table 1 Morphological features of one CU. Feature

Unit

Description

A* Atot hSUB hdrain hdowns

m2 m2 m m m

Area of cultivated basin Area of cultivated unit (CU) Substratum depth of the shallow groundwater Topographic level of subsurface drainage system Reference level for downstream groundwater flow

Table 2 Model parameters of one CU. Parameter

Unit

Description

Pmax Wsat WFC WWP Kc b n k

m/day m m m – mn+1/day – day1 day1 –

Maximum percolation rate Soil water storage at saturation Soil water storage at field capacity Soil water storage at permanent wilting point Crop coefficient Parameter related to the soil texture Parameter related to the soil texture Inverse of the drainage resistance Groundwater flow coefficient Specific yield of the shallow groundwater

c l

and drain discharge measurements. A piezometer survey is set up in the test parcel to assess water table fluctuations and to characterise the movement of the shallow groundwater along two small-scale toposequences. Focusing on periods with intense irrigation groundwater levels and drain discharges were monitored, in the summer of 1995, at a daily time step. The drain discharges were measured manually by using a bucket and stop watch in each junction box. This system consists of a cement structure with a diameter of 1.2 m and a depth of 2 m. It is made to connect the subsurface pipe drains to a collector pipe for discharge to an open earthen channel that discharges the drainage water to the Chott el Gharsa. The shallow groundwater in the Segdoud region is composed of alternate sand and sandy clay layers (Askri, 2002). The depth of the gypsum crust was measured in the whole oasis using a soil sampling tube. The depth varied between 0 and 4 m and was about 3 m in the test parcel. We assumed that palm trees are the only crop cultivated inside the oasis. The crop-specific parameter (Kc) is derived from the literature (Allen et al., 1998) and is adjusted for local conditions on the base of model calibration. Saturated hydraulic conductivity (Ks) and volumetric water content at saturation (hSAT), field capacity (hFC) and permanent wilting point (hWP) are obtained in laboratory from soil samples taken in the field at depths of 10 and 60 cm (Fig. 6). Table 3 presents the values of these properties for each soil layer. The soil water storages at saturation (WSAT), field capacity (WFC) and permanent wilting point (WWP) are obtained using the values of the soil thickness and water contents. Soil water holding properties are found to be uniform over the first meter of depth, with a typical plant water availability of 80–100 mm/m. This range of water availability is higher than the main daily potential evapotranspiration which is about 4.5 mm/day. Sand with water stored in the soil zone at saturation (438 mm/m), at field capacity (140 mm/m) and at wilting point (54 mm/m) is assumed representative of the test parcel for simulation. The maximum percolation rate (Pmax) that is derived through measurement of the saturated hydraulic conductivity varies between 0.96 and 3.36 m/day. An average of 2 m/day is deemed appropriate for the soil in the test parcel. This ecosystem is bordered in the north by a palm plantation covering an area of about 4 ha (Fig. 6). The palm plantation is located in the upstream area. The groundwater inflow from this ecosystem varies according to the hydraulic gradient between the palm plantation and the test parcel. In the present study, this inflow is assumed constant in time for all simulations. Its value was estimated by using the Darcy equation (0.024 m/day).

Fig. 5. Geographical location of the Segdoud oasis.

53

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

N Abandoned plot

Palm plantation

LEGEND:

P3

Irrigation earthen canal P5

P7

P4

P2 Uncropped area

Abandoned plot

Irrigation outlet Irrigation concrete canal

P1

P6

Pipe drain Collector drain Junction boxes Piezometer Soil sam pling

D3

D1

D2

0

20 m

Fig. 6. Layout of the test parcel and instrumentation.

Model calibration and validation The base model unit is fixed to 10 m  10 m, which is compatible with palm tree spacing to control input parameters such as climatic and irrigation. In the model development, six calibrated parameters are retained in defining the key hydrologic processes within each CU. These parameters can be classified into four main groups: (i) the crop-specific parameter, Kc, (ii) soil parameters including b and n (Eq. (9)), (iii) shallow groundwater parameters including the specific yield, l, and the parameter related to the groundwater flow, c, and (iv) the inverse of the drainage resistance, k. The values of parameters a and n in Eq. (8) are deduced from literature. For a sandy loam soil, Ashraf (1997, 2000) obtained a = 1.20  107 mn+1/day and n = 4. Thus b = 1.52a = 1.82  107 mn+1/day. These values of parameters b and n are used in all simulations. The shallow groundwater is assumed homogeneous within the test parcel. The model was calibrated in two steps. First, the values for Kc, l and c parameters were selected automatically by using the Fibonacci search technique to minimise the differences between the observed and the simulated groundwater levels variation with time (Vardavas, 1989). The objective function is thus taken to be (Sahoo et al., 2006):

f ¼

NP X ðho;i  hm;i Þ2

ð25Þ

i¼1

where ho and hm are the observed and the simulated groundwater levels, respectively, NP is the number of paired observed–simulated values. The optimised parameter values are given in Table 4. In the second calibration step, the value of the parameter, k, was determined for drains D2 and D3 using a manual calibration procedure, which aims to maintain the drainage discharge within observed ranges. The drain D1 may intercept the groundwater flow from the upstream palm plantation. Thus, it was not considered in the calibration and validation analyses. The accuracy of the model prediction is analysed using a correlation coefficient (R) and Nash–Sutcliffe efficiency coefficient (Neff). Results and discussion Model calibration The model was calibrated using groundwater levels and drain discharges obtained from 19th May to 24th July, 1995. Initial water

Table 3 Major soil parameters used in the model.

Saturated hydraulic conductivity (cm/h) Volumetric water content at saturation (cm3/cm3) Volumetric water content at field capacity (cm3/cm3) Volumetric water content at permanent wilting point (cm3/cm3)

Layer 0–20 cm

Layer 50–70 cm

4–19 0.26–0.29

5–39 0.26–0.28

0.07–0.12

0.07–0.09

0.03–0.04

0.02–0.04

contents of the irrigated and non-irrigated lands were assumed to be at field capacity (hFC) and permanent wilting point (hWP), respectively. Initial groundwater level was assigned for each CU based on the water table depths measured on 19th May, 1995. These depths varied between 0.72 and 1.46 m with an average of about 1.14 m. Groundwater level Fig. 7 shows the time series of simulated and observed groundwater levels in three locations during the calibration period. The observed hydrographs show that water table fluctuation has a large spatial variability, which can be explained by the irrigation water distribution to basins which is unequal. This fluctuation varies in magnitude and duration according to the location of piezometers in the parcel. Hydrographs observed in piezometers P2, P3 are much more fluctuating than that observed in P1. The first piezometer group is located close to the concrete canal. In this zone, basins are cultivated with two levels of crop palm trees and market garden crops. These basins being those that receive the highest irrigation water depths. Piezometers P1 is located at about 100 m from the concrete canal. Basins located in this zone are under-irrigated due to seepage water in the earthen canal systems. As expected, the simulated groundwater levels showed a relatively peaky response when irrigation occurred. Rapid rises of about 20 cm were observed, and then the water table fell back rapidly to around its pre-irrigation elevation. The observed water table peak magnitudes are under-predicted by the model. This can be explained by two main factors. The first is the uncertainties in the estimation of irrigation water depths which are assumed uniform. However, several crops as well as seepage water in the earthen canal systems induce a large spatial variability of this input data. Field observations show that several basins located in the

54

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

Table 4 Values of model parameters. Parameter

Unit

Measurements

Pmax Wsat WFC WWP Kc b n k

m/day m m m – mn+1/day – day1 day1 –

2 0.438 0.140 0.054

c l

Literature

Calibration

0.94 1.82  107 4 15  103 7  104 0.2

southern sector of the parcel were not irrigated during several irrigation events. More irrigation gauging stations are required to get an accurate assessment of irrigation water discharge and duration at the basin level. The second factor inducing the under-prediction of the water table peak magnitude is the uncertainty in estimating the specific yield of the shallow groundwater. Using small values of this parameter produces fast increase of the groundwater level during the irrigation events. However, high values produce lower and flatter peaks of water table hydrograph. During the recession phase succeeding the irrigation ceases, small values of the specific yield produce fast decreases of groundwater levels. It is clear that the observed water table hydrographs are much more fluctuating than the simulated hydrographs. There are two explanations for the disagreement between the predicted and observed groundwater levels. The first being water losses at the irrigation outlet. These losses are frequently observed when surrounding parcels are irrigated. They may induce a significant water table increase mainly in area close to the irrigation outlet. Our field observations show that the farmer opens the irrigation gate continually during the summer season, so that basins located at proximity of the irrigation outlet should receive the maximum quantity of irrigation water. These losses are unknown. Thus, they are not considered in the simulation. The water dynamics in the test parcel is the consequence of several hydrologic processes taking place in this ecosystem and in the upstream areas. Thus, the second explanation of for that disagreement may be the under-estimation of the groundwater inflow from the upstream areas. The intense irrigation of the upstream palm plantation may induce a significant increase in this flow and thus a high recharge of the shallow groundwater. As the groundwater inflow is assumed constant, this assumption may result in the underestimation of its real value. The daily measurement of the groundwater level in the upstream areas is required to get proper values of the groundwater inflow to the test parcel. The model simulated the water table recession observed in piezometer P7 with acceptable accuracy. The differences between the measured and the simulated groundwater levels are lower than 0.11 m. For the other piezometers, OASIS_MOD under-predicts the rate of the water table recession following the end of irrigation. This can be attributed to under-estimation of the value of the parameter (c), which determines the groundwater flow. Using higher values for this parameter induces a general under-estimation of the observed groundwater levels. Table 5 summarises the statistical parameters comparing the predicted and the observed groundwater levels. The maximum differences between these values range between 0.09 m in piezometer P7 and 0.20 m in P4. The first piezometer is located in an uncropped area. The second is located at about 60 m from the concrete irrigation canal. In this case, the difference between the predicted and observed values can be explained by the irrigation data uncertainties. The ranges of the correlation coefficient and of the

coefficient of efficiency are 0.31–0.72 and 0.09–0.68, respectively. The low values of these criteria obtained in some piezometers can be explained by several factors including the poor quality of the irrigation data, the under-estimation of the groundwater inflow from upstream areas and uncertainties when estimating some model parameters. Drain discharge A comparison of measured and predicted drain discharges during the calibration period is shown in Fig. 8. The observed hydrographs show that drain discharge occurred continually but at vastly different rates according to the position of the pipe drain in the parcel. The average drain discharges measured in drains D1, D2 and D3 are 0.2, 0.6 and 0.2 mm/day, respectively. The drain D2 is located at the central part of the parcel in area managed with several cultivated basins. This drain may intercept the seepage water that originates from an earthen irrigation canal located at its proximity. However, drains D1 and D3 are located at the border between the test parcel and an uncropped area for the first drain and an abandoned parcel for the second. These two drains were installed at a distance of about 20 m from the first basins series and at about of 30 m from the earthen canal system. It is clear that the model prediction averages the field observations of drain discharges. For drains D2 and D3, the mean absolute deviations between the simulated and observed drain discharges are 1.5 and 0.6 m3/day, respectively, and the Neff values are 0.17 and 0.29, respectively (Table 5). The values of Neff are less than 0.5 showing a disagreement between the predicted and observed drain discharges. The coefficient of correlation between these two variables is 0.58 for drain D2 and 0.70 for drain D3. The OASIS_MOD model under-predicts the values of the maximum drainage discharges when irrigation occurs. The discrepancies between the observed and simulated hydrographs can be explained by several factors including the uncertainties of drain discharge measurements and the model simplifications. As pipe drains are usually inaccessible in the junction boxes, the manual method used for the measurement of the drain discharges may under-estimates the real values. Drainage discharge is assumed to occur when the groundwater level is higher than the topographic level of the drainage system. This assumption may under-estimate the drain discharge, because the intercepted seepage water from earthen irrigation canal systems is not considered. Model validation The model was validated through independent measurements of groundwater levels and drain discharges from 25th July to 30th September, 1995. Initial water storages within the irrigated and uncropped lands were obtained from the calibration simulation. Initial groundwater levels were obtained on the base of values of the water table depths measured on 24th July, 1995. These depths varied between 0.77 and 1.55 m, with an average of about 1.20 m. Groundwater level Measured and simulated groundwater levels are shown in Fig. 7 for the validation analyses. In all piezometers, the measured hydrographs show two different phases. The first is a declining phase of the groundwater level which was observed between 25th July and 4th September, 1995. A similar trend was observed in several Tunisian oases during the summer season as the consequence of the increase in the water uptakes by soil evaporation and crop transpiration (Mtimet and Hachicha, 1995). The second phase began on 4th September, 1995 and continued until the 30th of September, 1995. It is characterised by a water table peak whose magnitude varies in the range 0.12–0.54 m. The model captured well

55

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

(a) (mm)

Rainfall/Irrigation

160

(b) 160

140

140

120

120

100

100

80

80

60

60

40

40

20

20

0 19/5

29/5

8/6

18/6

28/6

8/7

0 24/7

18/7

3/8

13/8

23/8

2/9

12/9

22/9

1995

1995

Groundwater levels (m)

(m)

51.7

(P1) 51.6 51.5 51.4 51.3 51.2 19/5

29/5

8/6

18/6

28/6

8/7

18/7

52 51.9 51.8 51.7 51.6 51.5 51.4 51.3 51.2 51.1 51 24/7

(P1)

1995

3/8

13/8

52.2

(P2)

52.1 52 51.9 51.8 51.7 51.6 29/5

8/6

18/6

28/6

8/7

18/7

52.5 52.4 52.3 52.2 52.1 52 51.9 51.8 51.7 51.6 51.5 24/7

12/9

22/9

(P2)

3/8

13/8

23/8

2/9

12/9

22/9

1995

1995

(m)

(m)

51.5

(P7) 51.4 51.3 51.2 51.1 51 19/5

2/9

(m)

(m)

51.5 19/5

23/8

1995

1995

29/5

8/6

18/6

28/6

8/7

18/7

1995

51.8 51.7 51.6 51.5 51.4 51.3 51.2 51.1 51 50.9 50.8 24/7

(P7)

3/8

13/8

23/8

2/9

12/9

22/9

1995

Fig. 7. Simulated (continuous line) and observed (discontinuous line) groundwater levels, expressed in meters, in three locations. (a) Calibration period. (b) Validation period.

the falling trend of the groundwater level observed between the 25th July and the 4th September, 1995. It reproduced the increasing trend in the groundwater level when irrigation occurs with acceptable accuracy. However, it is unable to predict accurately the water table peak registered the 4th September, 1995 perhaps resulting from intense irrigation in upstream palm plantation. The groundwater inflow from the upstream areas is assumed constant in time. This assumption is not realistic mainly when this ecosystem is intensively irrigated. Table 5 shows the statistical parameters obtained during validation. The disagreements between observed and simulated groundwater levels registered on 4th September, 1995 at some piezometers resulted in obtaining low values of R and Neff. The first parameter varies between 0.16 and 0.72. The second ranges from

0.08 to 0.34. These results show that the quality of the model prediction varies from one piezometer to another according to the precision in estimating the irrigation water depths and the groundwater inflow from upstream areas. Drain discharge A comparison of measured and predicted drain discharges is shown in Fig. 8 for validation analysis. The two hydrographs are quite similar for drain D2, showing a triple peak due to the influence irrigation events. The model accurately reproduced the time to start of drainage discharge increase occurring one day after irrigation starting. However, discrepancies were found in the timing of peak discharges. The predicted curve shows that the drain discharge attains maximum value two days after irrigation stops.

56

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

After calibration and validation, OASIS_MOD model is used in a hypothetical scenario to simulate the impact of this irrigation interval on the water and salt dynamics within the test parcel for the period from July 1989 to December 1996.

Table 5 Statistical analysis for calibration and validation. Calibration

lm

Validation

ls

R

Neff

lm

ls

R

Neff

Piezometer P1 51.43 P2 51.89 P3 52.46 P4 51.77 P5 51.57 P6 51.23 P7 51.20

51.25 51.89 52.41 51.76 51.62 51.37 51.20

0.67 0.58 0.60 0.63 0.36 0.60 0.72

0.68 0.39 0.25 0.49 0.03 0.14 0.09

51.30 51.82 52.29 51.60 51.48 51.16 51.23

51.25 51.82 52.35 51.62 51.56 51.16 51.22

0.60 0.56 0.41 0.72 0.54 0.61 0.16

0.1 0.15 0.03 0.35 0.23 0.31 0.15

Drain D2 10.4 D3 3.7

10.7 3.8

0.54 0.63

0.17 0.29

7.2 2.5

8.6 2.2

0.55 0.04

0.09 0.37

l, mean value. m, measured. s, simulated.

However, experiments realised in Tozeur oasis showed that the maximum subsurface drainage flow is attained six hours after the irrigation starting (Mtimet and Hachicha, 1995). Thus, the discrepancy in the timing of peak discharges can be explained by the choice of the computing time step. The model predicted that observed peak discharges are maintained for one day afterwards. For drain D3, the simulated hydrograph is a continue declining curve. In this case, the model under-predicts the maximum values of the drain discharges after irrigation starting. The total drainage discharges predicted during the validation period were 43 and 14 mm compared to the observed total of 33 and 11 mm for drains D2 and D3, respectively. The disagreements between observed and simulated drain discharge resulted in obtaining low values of R and Neff (Table 5). Numerical experiment In the Segdoud oasis, an irrigation schema with 5h_25 min duration is planned with regular irrigation interval (10–14 days).

(a)

The irrigated soil Fig. 9 shows the daily variability of the simulated soil salinity concentration and water stored in an irrigated land close the irrigation outlet from 29th May to 8th July, 1995. It is clear that the soil salinity increases as the soil water storage decreases. The soil salinity can vary from day-to-day, because it largely depends on meteorological conditions and irrigation practices. It changes from 13 g/ L before irrigation to about 9 g/L one day after irrigation stop. This decreasing trend is produced by irrigation with fresh water, due to leaching (Patel et al., 2000). The mixing of 4 g/L irrigation water with a soil solution of 13 g/L leads to its dilution. During the summer season (June–August), the duration of this decreasing trend is of about one day after irrigation event. With time, the soil water continues to be removed by evapotranspiration, the dilution effect diminishes, and the soil salinity beings to increase with time. The soil salinity reached seven days after irrigation is a similar value to the soil salinity before irrigation. During the summer season, the average increasing rate of the soil salinity concentration is about 0.5 g/L/day. The effect of irrigation is also clear in the soil water storage. Before irrigation, this storage is less than the storage at field capacity. When irrigation occurs, the soil water storage increases and attains maximum value of about 200 mm which represents about 50% of the storage at saturation. After irrigation stops, the soil water storage decreases due to the percolation and to

(b) (m3/day)

(m3/day) 14

(D2)

12 10 8 6 4 2 0 25/7

The input data In the scarce data, the initial soil water content is assumed at permanent wilting point. An arbitrary initial water table depth value of 2.1 m is used for each CU. This value was measured on May 1995 in an uncropped area located upstream of the oasis. The initial soil and groundwater salinity concentrations are assumed 3 g/ L. This value was obtained in a laboratory on soil samples taken in the field on April 1986 (Ben Marai, 1986).

4/8

14/8

24/8

3/9

13/9

20 18 16 14 12 10 8 6 4 2 0 19/5

23/9

(D2)

29/5

8/6

18/6

28/6

8/7

18/7

1995

1995

(m3/day)

(m3/day)

5

(D3)

7

(D3)

6

4

5 3

4 3

2

2 1 0 25/7

1 4/8

14/8

24/8

1995

3/9

13/9

23/9

0 19/5

29/5

8/6

18/6

28/6

8/7

18/7

1995

Fig. 8. Simulated (continuous line) and observed (discontinuous line) drain discharge, expressed in (m3/day), in two pipe drains (D2 and D3). (a) Calibration period. (b) Validation period.

57

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

evapotranspiration. Three days later, the soil water content attains the field capacity, and thus the percolation processes. If an irrigation interval of about ten days is maintained constant during the year, the monthly average of irrigated soil salinity shows cyclical seasonal fluctuations with higher values during the summer season and lower values after the autumn–winter rains (Fig. 10). These fluctuations are explained by the rising of salts normally contained in the soil and by the amount of salts

added through irrigation during the spring–summer period when evapotranspiration prevails, and by leaching of the salts during the autumn–winter period, when percolation prevails (Tedeschi and Menenti, 2002). During the last period, fresh water may be available in very limited quantities during the few rainfall events. In arid and semi-arid zones, this fresh water can be used to flush out salts accumulated in the root zone during the summer season (Mtimet and Hachicha, 1995). Using an irrigation interval of about

20 Soil-water storage

Soil water storage (mm)

400

Soil salinity concentration (g/L)

450

18

Soil salinity concentration

16

350

14

300

12

250

10 200

8

150

6

100

4

50

2

0 29/5

0 3/6

8/6

13/6

18/6

23/6

28/6

3/7

8/7

Date Irrigation

Irrigation

Irrigation

Irrigation

Fig. 9. Simulated daily water storage and salinity concentration of irrigated land from 29th May to 8th July 1995.

(a)

16

155

14 150

12

145

10

140

8 6

135

4

April-97

Feb-97

Aug -96

Feb-96

Aug-95

nov-96

0 May-96

2

125

nov-95

130 May-95

Soil water storage (mm)

(b)

18

Soil salinity concentration (/L)

Soil salinity concentration (g/L)

Soil-water storage (mm)

160

Water table level

7.4

groundwater salinity concentration

52.05

7.2

52

7

51.95

6.8

51.9 6.6

51.85

6.4

51.8

Groundwater salinity concentration (g/L)

Water table level (m)

52.1

6.2

51.75 51.7 nov-96

Aug-96

May-96

Feb-96

nov-95

Aug-95

May-95

6

Fig. 10. Simulated monthly averages of (a) water storage and salinity concentration of irrigated land from May 1995 to April 1997 and (b) groundwater level and groundwater salinity from May 1995 to November 1996.

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

ten days leads to soil salinity increases during the summer season. During monsoon months, using this interval maintains the soil salinity at about 5 g/L. This can be very useful to define irrigation strategies with saline water (e.g. amount of saline water applied in each irrigation, time between two irrigations) that can avoid risks of irreversible salination of these sandy soils (Moreno et al., 2004). The simulated results show significant inter-annual variability of the soil water storage and salinity concentration due to variability in irrigation water depth and reference evapotranspiration (Fig. 11). The annual averages of soil water storage decreased from 240 mm in 1990 to 146 mm in 1995 due to the reduction in the irrigation water allocation. During this period, the irrigated area inside the test parcel has increased from 600 to about 1200 m2 resulting in a decrease in the average irrigation water depth from 120 to 60 mm per irrigation event. The annual averages of soil salinity have increased from 3 g/L in 1989 to 9 g/L in 1995 due to the decrease in leaching. This increasing trend has already been observed in the experimental parcel of Ksar Gheriss, situated in the center of Tunisia (CRUESI, 1970). The experiment, which is one of the few longterm experiments being conducted in Tunisia, was monitored between 1964 and 1969 on a sandy soil and using saline irrigation water (4 g/L). The obtained results showed that using an irrigation interval of 10 days has led to soil salinity increases from 3mmhos/ cm in September 1965 to about 10 mmhos/cm in July 1968. In July 1989, when irrigation started in the Segdoud oasis, the irrigation discharge and duration were enough to ensure the irrigation of all basins. Over the years, the basins area has increased, so that the irrigation allocation was insufficient to bring water to all basins. This problem has led to the abandonment of several basins located in southern portion of the test parcel. Thirteen basins are assumed abandoned annually during the period between July

(a)

1990 and 1994. The simulated soil salinity in the test parcel for June 1990, 1992, 1993, 1994 and 1996 in the settings appear in Fig. 12. The results show remarkable differences in the salinity conditions across the test parcel. High saline soil dominates the southern portion of the parcel. This contrasts with the northern portion characterised by relatively low soil salinity. Several basins located in the first portion are assumed abandoned since July 1991. Thus, soils are less well leached in the absence of irrigation. The soil salinity within these abandoned areas rises constantly until equilibrium is reached between input of salt by capillary rise and leaching after heavy rainfall. Fig. 12 shows two isolated areas of relatively high soil salinity that develop since June 1994 and 1996. These areas represent two basins that were abandoned in July 1993. The same figure shows an isolate area of relatively low soil salinity within the high salinity southern portion of the test parcel in June 1996. This area represents a cultivated basin which was maintained irrigated until June 1996. The simulation showed that the soil salinity increases yearly in the whole parcel area. The highest increasing rate occurred in the southern portion. The area of soils affected by salination increases yearly as observed by the increase in the number of abandoned basins. In June 1992, only about 2% of the whole parcel area had the soil salinity higher than 15 g/L. The area affected by that level of soil salinity has increased to about 50% of the whole parcel area in June 1996. The shallow groundwater Fig. 10 shows the predicted monthly average groundwater level and groundwater salinity from June 1995 to December 1996. It is clear that these two characteristics have a cyclical seasonal variation. The groundwater level has an increasing trend from August to February as the result of the evapotranspiration decrease. During this wet period, the soil water content is above the field capacity,

10

260 Soil water storage Soil salinity concentration

Soil water storage (mm)

240

9

220

8

200 7 180 6 160 5

140 120 1990

(b)

1991

1992

1993

1994

1995

Soil salinity concentration (g/L)

58

4 1996

water table level (m) Groundwater salinity concentration

52

7 6

51.8 51.7

5

51.6

4

51.5 3

Groundwater salinity concentration (g/L)

Water table level (m)

51.9

51.4 51.3 1990

1991

1992

1993

1994

1995

2 1996

Fig. 11. Simulated annual averages of (a) water storage and salinity concentration of irrigated land and (b) groundwater level and groundwater salinity over the 7 (1990– 1996) years.

59

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

and thus the percolation rate attained maximum values. The lower groundwater level occurs during the summer season due to increased evapotranspiration. Groundwater salinity has a different trend and increases from May to September. During this hot and dry period, the accumulated salts in the soil profile are frequently carried back to the shallow groundwater after the addition of irrigation water to the land surface. The shallow groundwater salinity has a decreasing trend from October to April. During this rainy period, the soil salinity attains minimal value. Thus, the percolation process may induce the dilution of the groundwater solution. The annual average groundwater level simulated by OASIS_MOD model shows an increasing trend from the year 1990 to

m

(a) 0

20

40

N

1995 (Fig. 11). This trend is explained by the increase in the cultivated area inside the test parcel. As the irrigation water depth is assumed constant through years, the annual increase in cultivated area inside the test parcel may increase the amount of water that drains to the shallow groundwater. The annual average groundwater salinity has increased from 3 g/L in 1989 to about 7 g/L in 1995 as a consequence of the soil salinity increase. The increasing trend of the annual averages groundwater level and groundwater salinity has been noticed by Ben Marai and Belarbi (1995). The observed groundwater depth has decreased in the test parcel from 1.07 m in September 1994 to 0.91 m in 1995. During this period, the observed groundwater salinity has increased from 6.8 to 7.4 g/L.

(b)

g/L

0

40 m

0

40 m

4 8 12 16 20 24

g/L

0

40 m

g/L

(d)

(d)

5

(c)

5

5

5

10 10

15 15

25

20 20

25 25

45

30 30

(e)

0

g/L 5

40 m

(f) ☼

g/L

0

40 m

5 10

15

15

2

25

35

35

45

45 55

55

Irrigation outlet Piezometer Fig. 12. Spatial distribution of measured groundwater level in 19th May 1995 (a) and simulated soil salinity concentration (g/L) with OASIS_MOD in the test parcel, on five dates (b): 30th June, 1990, (c) 30th June, 1992, (d) 30th June, 1993, (e) 30th June, 1994 and (f) 30th June, 1996.

60

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61

Conclusion OASIS_MOD is a conceptual, daily, semi-distributed hydrologic model which has been developed to predict the daily soil salinity, groundwater level and drain discharge within modern Tunisian oases given the daily irrigation water depth and salinity concentration, rainfall and reference evapotranspiration. The spatially semi-distributed approach adopted in this model allowed a reliable simulation of spatial and temporal variation of irrigation management within an oasis ecosystem. The application of the model to the oasis scale provides estimates, with good accuracy, of the irrigation water depth and frequency at each parcel separately. The lack of available irrigation data at 10,000 basins that constitute the cultivated area inside the oasis is problematic for the model application at the oasis scale. For instance, the model will help to anticipate the effects of irrigation interval and climatic conditions on soil salinity and groundwater level at the parcel scale. It is also useful in the planning of flood irrigation with saline water and the management of the drainage systems. The model tests are conducted with field data collected at a 1.5 ha parcel of farmland situated in the Segdoud oasis. It involves calibration and validation analyses on independent data sets. The calibration results show that the model under-predicts some water table peaks when irrigation occurs. It is concluded that the bias is due to inconsistencies in the irrigation data at the basin level. The model prediction averages the observed drainage discharges. The discrepancy between the simulated and observed hydrographs is explained by the model simplification and uncertainties of drain discharge measurements. The validation results show that deviations between observed and simulated groundwater level have increased to about 0.5 m due to under-estimation of the groundwater inflow from the upstream palm plantation. The daily measurement of the water table depth in the upstream area outside of the model domain is required to get proper values of the groundwater inflow to the test parcel. OASIS_MOD is calibrated and validated using data monitored during the summer season 1995. The application model should be calibrated and validated to a longer period of record in order to capture additional climatic variability. The simulation results show that the soil salinity in modern Tunisian oases depends on the climate condition and the availability and the management of the irrigation water. The irrigation interval is found to have an important effect on the components of the water and salt balances in the soil zone and on the groundwater level. Using an irrigation interval of approximately 10 days causes decrease soil water storage and increase soil salinity during the spring and summer seasons. Using this frequency leads to groundwater level increase during the autumn and winter seasons. Thus, a higher irrigation frequency during the summer season should result in higher water storage and less salt accumulation in the soil zone. A lower irrigation interval during the autumn and winter seasons should cause a water table decrease. References Allen, R.G., Pereira, L.S., Raes, D., Smith, M., 1998. Crop evapotranspiration guidelines for computing crop requirements. FAO Irrig. Drain. Report modeling and application. J. Hydrol. 285, 19–40. Ashraf, M., 1997. Dynamics of soil water under non-isothermal conditions. PhD Thesis, University of Newcastle upon Tyne, UK. Ashraf, M., 2000. Water movement through soil in response to water-content and temperature gradients: evaluation of the theory of Philip and de Vries (1957). J. Eng. Appl. Sci. Pakistan 19, 37–51. Askri, B., 2002. Modélisation des processus de salinisation des sols irrigués en zones arides. Cas de l’oasis de Segdoud. PhD Thesis, University Tunis II, Tunisia, 170 p. Bahceci, I., Dinc, N., Tari, A.F., Agar, A.I., Sonmez, B., 2006. Water and salt balance studies, using SaltMod, to improve subsurface drainage design in the KonyaCumra Plain. Agric. Water Manage. 14, 261–271.

Belloumi, M., 2006. A stochastic frontier approach for measuring technical efficiencies of date farms in southern Tunisia. Agric. Res. Economies Review, 1–12. Ben Marai, M., 1986. Etude pédologique de périmètre du Sagdoud. Etude no. 1 2051/ E, Direction des sols, Ministère de l’Agriculture, Tunisie, 6 pp. Ben Marai, M., Belarbi, A., 1995. Surveillance de la salinisation et contrôle de la nappe du périmètre Segdoud. Direction des sols, Ministère de l’Agriculture, Tunisie, 10 pp. Ben Mohamed, M., 2002. Geothermal utilization in agriculture in Kebeli region, southern Tunisia. GHC Bulletin, 27–32. Boughton, W.C., 1968. A mathematical catchment model for estimating runoff. J. Hydrol. N.Z. 7, 75–100. Bouraoui, F., Vachaud, G., Normand, B., 1997. A distributed physical approach for surface–subsurface water transport modeling in agricultural watersheds. J. Hydrol. 203, 79–92. Brunel, J.P., Ihab, J., Droubi, A., Samaan, S., 2006. Energy budget and actual evapotranspiration of an arid oasis ecosystem: Palmyra (Syria). Agric. Water Manage. 213, 213–220. Connell, L.D., Gilfedder, M., Jayatilaka, C., Bailey, M., Vandervaere, J.P., 1999. Optimal management of water movement in irrigation bays. Environ. Modell. Softw. 14, 171–179. CRUESI, 1970. Research and training on irrigation with saline water. UNESCO/UNDP, 208 p. Doorenbos, J., Pruitt, W.O., 1977. Guidelines for Predicting Crop Water Requirements. FAO, Rome, p. 144 (FAO Irrig. Drain. Paper 24, Revised). Eilers, V.H.M., Carter, R.C., Rushton, K.R., 2007. A single layer soil water balance model for estimating deep drainage (potential recharge): an application to cropped land in semi-arid north-east Nigeria. Geoderma 140, 119–131. Gardner, W.R., 1958. Some steady-state solution for the unsaturated moisture flow equation with application to evaporation from a water table. Soil Sci. 85, 228– 232. Horchani, A., 2007. Water in Tunisia: a notional perspective. In: Proceedings of workshop in Tunisia. Agric. Water Manage., 88–96. Jorenush, M.H., Sepaskhah, A.R., 2003. Modelling capillary rise and soil salinity for shallow saline water table under irrigated and non-irrigated conditions. Agric. Water Manage. 61, 125–141. Kadri, A., Van Ranst, E., 2002. Oasis crop production constrains and sustainable development strategies. Science et changement planétaires. Sécheresse 13 (1), 5–12. Kassah, A., 1996. Les oasis tunisiennes: Aménagement hydro-agricole et développement en zone aride. Série géographique (No. 13), Centre d’études et de recherches économique et sociales, Tunis. Konukcu, F., Gowing, J.W., Rose, D.A., 2006. Dry drainage: a sustainable solution to waterlogging and salinity problems in irrigation areas? Agric. Water Manage. 83, 1–12. Lamsal, K., Guna, N.P., Saeed, M., 1999. Model for assessing impact of salinity on soil water availability and crop yield. Agric. Water Manage. 41, 57–70. Lasram, M., 1990. Les systèmes agricoles oasiens dans le Sud de la Tunisie. Options Méditerranéennes 11, 21–27. Liu, Y., Pereira, L.S., Fernando, R.M., 2006. Fluxes through the bottom boundary of the root zone in silty soils: parametric approaches to estimate groundwater contribution and percolation. Agric. Water Manage. 84, 27–40. Mermoud, A., Seytoux, H.J., 1989. Modélisation et observation du flux hydrique vers la surface du sol depuis une nappe peu profonde. Hydrol. Continent 4 (1), 11– 23. Moran, M.S., Peters, C.D., Watts, J.M., Mc Elroy, S., 2004. Estimating soil moisture and the watershed scale with satellite based radar and land surface models. Can. J. Remote Sen. 30 (5), 805–826. Moreno, F., Vaz, R., Fernandez-Boy, E., Cabrera, F., 2004. Simulating the composition of the in situ soil solution by the model EXPRESO: application to a reclaimed marsh soil of SW Spain irrigated with saline water. Agric. Water Manage. 66, 113–124. Mtimet, A., 1983. Soils of Tunisia. Options Méditerranéennes 34, 243–258. Mtimet, A., Hachicha, M., 1995. Salinisation et hydromorphie dans les oasis tunisiennes. Science at changement planétaires. Sécheresse 4 (6), 319–324. Ninghu, Su, Matthew, B., Louise, M., Alfred, H., 2005. Simulating water and salt movement in tile-drained fields irrigated with saline water under a serial biological concentration management scenario. Agric. Water Manage. 78, 165– 180. Patel, P.M., Prasher, S.O., Bonnel, R.B., 2000. Effects of water table depth, irrigation water salinity, and fertilizer application on root zone salt building. Can. Agric. Eng. 42 (3), 111–115. Raes, D., Deproost, P., 2003. Model to assess water movement from a shallow water table to the root zone. Agric. Water Manage. 62, 71–91. Rassam, D., Cook, J.F., 2002. Numerical simulations of water flow and solute transport applied to acid sulphate soils. J. Irrig. Drain. Engineering (March/ April), 107–115. Renschler, C.S., Cachron, T.H., Diekkringer, B., 2001. Regionalization methods for watershed management – hydrology and soil erosion fro, point to regional scale. In: Stoo et al. (Eds.), Sustaining the Global Farm. Selected Paper from 10th International Soil Conservation Organization Meeting Held May 24–29, 1999 at Purdue University and USDA-PRS National Soil Erosion Research Laboratory. Richter, M., 1986. The climate of oases and its ecological significance. Appl. Geogr. Dev. 29, 44–72.

B. Askri et al. / Journal of Hydrology 380 (2010) 45–61 Rodriguez, F., Andrieu, H., Morena, F., 2008. A distributed hydrological model for urbanized areas – model development and application to case studies. J. Hydrol. 351, 268–287. Roost, N., Cai, X.L., Turral, H., Molden, D., Cui, Y.L., 2008. Adapting to intersectoral transfers in the Zhanghe irrigation system, China. Part II: impacts of in-system storage on water balance and productivity. Agric. Water Manage. 95, 685–697. Sahoo, G.B., Ray, C., De Carlo, E.H., 2006. Calibration and validation of physically distributed hydrological model, MIKE SHE, to predict streamflow at high frequency in a flashy mountains Hawaii stream. J. Hydrol. 327, 94–109. SANYU Consultants Inc., 1996. Projet d’amélioration des périmètres irrigués dans les oasis du sud en république de Tunisie. Agence japonaise de coopération internationale. Rapport Final, 300 pp. Sellami, M.H., Sifaoui, M.S., 2008. Modelling of heat and mass transfer inside a traditional oasis: experimental validation. Ecol. Modell. 210, 144–154. Sghaier, M., 1984. Identification et analyse des systèmes de production agricole dans les oasis de Nefzaoua. INAT, Tunisia. Singh, R., 2004. Simulation on direct and cyclic of saline waters for sustaining cotton– wheat in a semi-arid area of north-west India. Agric. Water Manage. 66, 153–162.

61

Singh, R., Subramanian, K., Refsgaard, J.C., 1999. Hydrological modelling of a small watershed using MIKE SHE for irrigation planning. Agric. Water Manage. 41, 149–166. Srinivasulu, A., Sujani, Rao Ch., Lakshmi, G.V., Satyanarayana, T.V., Boonstra, J., 2004. Model studies on salt and water balances at Konanki pilot area Andhra Pradesh, India. Irr. Drain. Systems 18, 1–17. Tedeschi, A., Menenti, M., 2002. Simulation studies of long-term saline water use: model validation and evaluation of schedules. Agric. Water Manage. 54, 123– 157. Vardavas, I.M., 1988. A simple water balance daily rainfall-runoff model with application to the tropical Magela Creek catchment. Ecol. Modell. 42, 245–264. Vardavas, I.M., 1989. A fibonacci search technique of model parameter section. Ecol. Modell. 68, 65–81. Wang, Y.R., Kang, S.Z., Li, F.S., Zhang, L., Zhang, J.H., 2007. Saline water irrigation scheduling through a crop-water-salinity production function and a soil-watersalinity dynamic model. Pedosphere 17 (3), 303–317. Xu, P., Shao, Y., 2002. A salt-transport model within a land-surface scheme for studies of salinization in irrigated areas. Environ. Modell. Softw. 17, 39–49.