Groundwater Modeling of Multi-Aquifer Systems Using GMS

The Nubian Sandstone complex in the western desert is a part of the major regional Nubian aquifer system comprising West Africa. The different rock strata deposited in localities forming the units of the aquifer system. Farafra Oasis lies in the northern part of the Dakhla basin. In its southern region, the Nubian Sandstone (deep aquifer) is overlained by Dakhla shale but in the central and northern regions by fissured chalky limestone (shallow aquifer). The two overlaying aquifers in Farafra Oasis represent a typical hydrogeological model of a huge multilayered artesian basin extending over the territory of Egypt. The Post Nubian aquifer played an essential role in the development of Farafra Oasis for a long time through the many springs issued from this aquifer. The rapid drilling process of deep wells started in 1960s led to stop flowing of many springs and wells plus the depletion in discharges and pressure of many others. Therefore, there is a real danger of either dewatering or increasing the water depths to uneconomic lifting depths for both the shallow and deep aquifers. A two-dimension flow model GMS (Groundwater Modeling System) was used to investigate this problem. Application of the present conditions indicated that drawdowns in the Post Nubian aquifer range from 5 m to about 9 m. The second scenario tries to sustain the groundwater utilities in the Post Nubian aquifer through a group of procedures. Accordingly, drawdowns are expected to range from 5 m to 8.6 m in the Post Nubian aquifer. According to this scenario, 3 m decline in the Nubian Sandstone aquifer followed by declining in the Post Nubian by about 1 m.


Introduction
Many aquifers in nature loss or gain water through adjacent confining beds of relatively low permeability. However, such leaky aquifers are often only part of multiple-aquifer systems. When water is withdrawn or recharged in one particular aquifer the head distribution in the entire system will be influenced. In such a system several aquifers can be distinguished, each separated from the others by aquitards (semi pervious layers). These aquifers and aquitards together form a single leaky system when recharge or discharge at one place influences the head distribution in all other parts, but is not felt at its boundaries. In practice, however, the system considered may be restricted by boundaries with a negligible influence during the period of interest. When water is withdrawn from one or more layers of a multiple-aquifer system, the induced drawdown will be dependent on the hydraulic properties of all aquifers and aquitards. Analytical solutions have recently been presented for several types of steady-state flow in systems comprising any number of layers [1]. In an attempt to combine the advantages of an analytical approach with the capability of numerical models to include heterogeneity, Hemker [2] developed a hybrid analytical-numerical solution for transient well flow in vertically heterogeneous aquifers. The radial components of flow are treated analytically, while the finitedifference technique is used to compute the vertical flow components in the horizontally layered aquifer. The resulting drawdown equations in the Laplace domain also account for the effects of a finite diameter pumped well and wellbore storage. To determine transmissivity values for aquifers and hydraulic resistances for aquitards, field investigations have to be conducted which usually involve one or more aquifer tests. As the presented multiple-aquifer well flow solution can only be used for steady state drawdowns, extrapolation of the observed time--drawdown relationship for each piezometer can provide the required information. This implies that for such an analysis a relatively long pumping period should be chosen, especially when aquitard resistances have to be obtained with some accuracy. The availability of analytical solutions for multiple-aquifer transient well flow will provide improved methods for aquifer test evaluation, probably substantially reducing the required pumping period. These practical implications and the improved possibilities to predict the response of a multiple aquifer system to abstractions have provided a stimulus to further extend the eigenvalue approach in well flow hydraulics. Saafan et al., [3] applied multi-objective genetic algorithm (MOGA) model in deep Nubian Sandstone aquifer in El-Farafra oasis, Egypt to develop the maximum pumping rate and minimum operation cost as well as the prediction of the future changes in both pumping rate and pumping operation cost. They concluded that the optimal pumping rate and the corresponding drawdown range from 190699 to 179423 m 3 /day and 6.13 to 8.34 m respectively. Moharram et al. presented optimization model based on the combination of the MODFLOW with GA. This model was used in deep Nubian Sandstone aquifer in El-Farafra oasis, Egypt El-Farafra oasis-Egypt to obtain optimal pumping rate during different simulation period.

The Study Area
Farafra Oasis is one of many morphotectonic depressions present in the Western Desert. It lies in the heart of these depressions between latitudes 26° 30 \ , 27° 30\ and longitudes 27° 30 \ and 29°00 \ inclosing an area of about 86200 Km 2 ( Figure 1). It forms an irregular triangular shape with an apex to Bahariya Oasis and base towards Dakhla Oasis. Groundwater is the only water source in this remote area. It is supplied from two main aquifers; the shallow Post Nubian aquifer and the deep Nubian Sandstone aquifer.
The Post Nubian aquifer (chalky limestone aquifer) played a vital role in the development of Farafra Oasis for a long time through the many springs issued in this aquifer. These springs formed the only source of water that was used by the inhabitants for domestic and irrigation purposes. Moreover, this water forms the easiest and cheapest water source in the Oasis. Since 1960, rapid drilling process of deep-water wells started to invest the high potentiality of the well-known Nubian Sandstone aquifer. Thus, huge quantities of groundwater are discharging from this aquifer (about 150 million m 3 / year, Ali, [4] to irrigate the continuously increasing reclaimed land (49,000 feddan) until now. The extensive exploitation of groundwater led to the depletion of the aquifer storage. Therefore, there is a real danger of either dewatering or increasing the water depths of wells to uneconomic lifting depths for both the shallow and deep aquifers. Nowadays, many springs dried out and the remnants are expected to face the same problem in the near future. Accordingly, many social and economic problems appeared ( Figure 1).
The present study aims at assessing this problem from the hydrogeological point of view. In this respect, it tries to decrease the rapid head decline and sustains the groundwater supplies of the Post Nubian aquifer through the application of some management plans using mathematical modeling techniques.

Hydrological settings
According to the previous, hydrological studies i.e. [5][6][7] and the lithology of the old and new drilled wells; Farafra Oasis is distinguished into two distinctive aquifers; the Post Nubian aquifer and the Nubian aquifer ( Figure 2).

Post Nubian aquifer
This aquifer is composed of limestone, chalky limestone with shale interbeds. A thin dolomitic bed characterizes the base of this aquifer. The aquifer thickness shows relatively small variations. It ranges from 90 m in the south at Abu Minqar area to about 170 m in north and west with an average value 130 m. Different oriented and mutually intersecting joints and fissures cut this aquifer. These structures provide the avenue for the upward movement of groundwater from deep horizons under artesian pressure [8]. Groundwater of this aquifer is supplied through number of natural springs and very few wells. Water, of these springs, is utilized by inhabitants for domestic and agricultural purposes.
The number and discharges of these springs were inventoried and compared with the previous data. It was found that the numbers and discharges of these springs had intensively decreased (Figures 3  and 4). Spring numbers were decreased from 67 springs at year 1962 to 11 springs now while the discharge of Ain Tenin for example was decreased from 156 m 3 /day at year 1962 to about 13 m 3 /day at present inventory. Although Ain El Balad showed a big anomaly, its discharge was estimated as 2050 m 3 /day indicating a deep fault structure origin. The high temperature of the spring water (36°C) enhances this result. On the other hand, Most of the dried springs are located at south where the flowing ones are located at extreme north ( Figure 5). This situation may be attributed to the topography of the depression surface, which slopes from south to northeast i.e. from +135 m to +50 m.
Based on the groundwater heads and pressure measurements of some non-flowing and flowing springs, groundwater level contour map of the Post Nubian aquifer was available ( Figure 6), groundwater levels of the karst aquifer are ranging from 105 m at El Sheikh Marzouq area to 60 m in the extreme north. Groundwater generally flows from south to north and from southeast to northwest. The smooth, regular and spaced contour lines reflect the homogeneity of the chalky limestone aquifer (Figures 3-6).

Nubian Sandstone aquifer
This aquifer consists of thick alternating sequence of coarse classics sediments of sandstone and clay beds. Accordingly, three successive water-bearing formations are well defined; Taref, Sabaya and Six Hills Formations, named Zone A, Zone B and Zone C ( Figure 2   beds that form aquicludes slow down the vertical movement of the groundwater through the aquifer layers. As a result, the sandstone beds form single aquifers in the regional aquifer complex. The bottom level of the Nubian aquifer in the study area varies from -1700 m to about -2300 m [9]. The number of wells tapping this aquifer increased from 18 wells in 1960s to about 140 wells in the present time. Hence, the pumping from the aquifer was increased in the last decade to reach about 145 million m 3 /year (Ali, 2003). This figure is expected to increase in the near future due to continuous increment in drilling of wells.
The aquifer is generally a high potential aquifer due to the large saturated and prevalence of sand facies in the aquifer sediments. According to Hamad, [10], transmissivity of Zone A recorded 148.6 m 2 /day with an average permeability 1.3 m/day while zone B recorded 1613 m 2 /day and 5.7 m/day. On the other hand, the third productive zone (zone C) recorded the highest transmissivity and permeability values of the aquifer (1642 m 2 /day and 7 m/day) reflecting the thick sandy aquifer section.
Based on the pressure measurements of some flowing wells and the depth to water in some non-flowing wells, a groundwater level contour map was constructed ( Figure 7). From this map, groundwater generally flows from south to north with the direction of head decline i.e. from +125 m at Abu Minqar area to about + 88 m at north Qasr El Farafra area. Two cones of depression were detected El Sheikh Marzouq and north Qasr El Farafra area indicated by the depression cones formed in such areas. This situation reflects the extensive exploitation of groundwater in those areas. These depression cones form local groundwater flow directions to the center of these cones. Sahl El Baraka and Sahl Qarawen areas showed small drawdown, no depression cones were formed indicated by regular contour lines reflecting the small well number and recent pumping in that areas (Figure 7).

Modeling of Groundwater
The model describes groundwater flow of constant density under non-equilibrium conditions in a heterogeneous and anisotropic medium according to the following equation (Bear [11]). K XX , K yy and Kzz are values of hydraulic conductivity (L T -l ); along the x, y, and z coordinate axes; h is the potentiometric head (L); W is the volumetric flux per unit volume and represents sources and/or sinks of water (T -l ); Ss is the specific storage of the porous material (L -1 ); and t is time (T). The model was used the finite difference approach to solve the groundwater flow equation.
Three dimensional, numerical time dependent flow model of finite difference GMS program was used to investigate the hydrodynamic impacts of the present discharging rates of the deep Nubian and Post Nubian aquifers on the groundwater levels in Farafra Oasis. Furthermore, it was used to study the proposed investigated plans to sustain the groundwater supplies from the Post Nubian aquifer, which faces acute depletion problem. Our approach involves the development of a two-layer groundwater flow model through the establishing of the boundary conditions, initial conditions and calibration of every layer in the model.

Boundary and initial conditions
According to the hydrogeological setting of the study area, two aquifer systems are present. Each one has its individual characteristics. The first (Post Nubian aquifer) is an unconfined aquifer of one homogenous layer. Its top is the ground surface and its bottom is the level of the marked dolomitic limestone bed that defines the border between the Nubian and Post Nubian aquifer. On the other hand, the Nubian Sandstone aquifer (the second layer) is a confined aquifer type. It was considered as one complex layer of vertical and horizontal variations in sediments and hydraulic parameters. The aquifer top is represented by the previously mentioned dolomitic limestone bed while the aquifer bottom is defined by the bottom of the third productive zone (Zone C) where the deepest wells were immersed (-1200 m). A grid of 3600 cell (60 columns and 60 row) was constructed to cover the model domain (Figures 8 and 9). This grid is refined in the area of springs for the first layer and wells for the second layer where detailed information on the aquifer properties is available. The groundwater system was built for the two layers by assigning the hydraulic conductivity, transmissivity, porosity, leakance and storativity to the grid cells using all information of the previous and present work.
Referring to the boundary conditions (Figures 8 and 9), no-flow boundary prescribed the eastern and western sides of the modeled area since it matches the flow directions and traces the major faults between the scarps and the depression area. General Head boundaries (GHB) are applied to the southern and northern boundaries of the model domain. The hydraulic heads and hydraulic conductance of both layers were assigned to boundary cells. For the upper layer (Post Nubian aquifer), the southern boundary was applied as 120 m head in the boundary and 520 m 2 /day as hydraulic conductance while the northern boundary was simulated as 50 m as head in the boundary and 520 m 2 /day as hydraulic conductance. On the other hand the southern heads of the lower layer (Nubian aquifer) was applied as 130 m as head on the boundary and 1250 m 2 /day as hydraulic conductance while they were assigned as 60 m and 750 m 2 /day at the northern boundary. These heads varies during the simulation process according to different stresses applied on the modeled area (Figures 8 and 9).

Model calibration
The model was firstly run under steady state conditions till the year 2004 then under transient conditions from 2004 to 2008. The calculated heads for both the first and second layer were compared by the measured ones at 2008. Big differences between the calculated and the measured heads are present, so the calibration of the model is critical and important. Calibration was carried out by a trial-and-error method through modifying the hydraulic parameters of the aquifer until a satisfied difference between the calculated and measured heads reached ( Figure 10). The difference was lowered to 0.05 m for the first layer and to 0.65 m for the second layer. As a result, head distribution maps resulted from the model are closely related to the actually measured maps (Figures 6 and 7).
The calibrated hydraulic conductivity values of the first layer range from 9 m/day to about 25 m/day. These values reflect the intensive secondary porosity characterize the Post Nubian aquifer and clarify Successful transient calibration depends mainly on the good estimation of hydraulic conductivities obtained from the steady state condition. Generally, specific storage for confined aquifer is the main parameter that is changed during the transient calibration. In the process of calibration of transient state, specific storage values were modified on a trial and errors basis, until a good match between the observed heads of years 2004 and 2008 and the calculated heads were achieved. The range of the resulted specific yield for post aquifer after the final calibration of the transient state was found to be varying from 0.18 to 0.26. The range of the resulted specific storage for Nubian sandstone after the final calibration of the transient state was found to be varying from 0.00001 to 0.0001. It can be seen that in general, there is good agreement between the observed and simulated head.

First Scenario
This scenario investigates the impact of the present discharging rates (490,000 m 3 /day) from about 140 water well, drilled in Nubian Sandstone aquifer and about 3000 m 3 /day from about 25 spring (Flowing, pumped springs and wells) engrossed in the Post Nubian aquifer. The model was run under these stresses and heads were estimated for every time step (next 20 years) to the upper and lower layer by the result extractor tool of the GMS program.
The head distribution maps of the Post Nubian aquifer (First layer) showed gradual decline in the groundwater levels indicated by the gradual migration of contour lines southward ( Figure 11). Contour line 85 passes through Ain Qarawen and immediately south Ain Shagy at the first year and it slowly moves southward to about 12 Km south Ain Shagy at the end of the simulation time (20 years). The smooth and displaced contour lines indicate aquifer homogenity. No changes in the groundwater flow directions were recorded due to the regional head decline in the whole aquifer.
Drawdown calculations in the Post Nubian aquifer (Table 1) indicate that the drawdowns range from 5 m to about 9 m. The maximum drawdowns are expected in Ain Tenin and Ain Shagy reflecting the direct impact of the intensive pumping of water from the Nubian Sandstone aquifer in these areas (Farafra and Elsheikh Marzouq). On the other hand, the minimum drawdowns are expected in the extreme north in the area of Ain El Makfi where no pumping occurs. The drawdowns also indicated that the springs that are still flowing in the northern part are expected to sand up within the next 3 years where the water level of these springs are about 2 m above the ground surface ( Figure 11).
On the other hand, head distribution maps ( Figure 12    Therefore, many wells are recommended to be drilled in the aquifer that will be discussed in the second scenario ( Figure 13).

Second Scenario
This scenario tries to sustain the groundwater utilities in the Post Nubian aquifer by substituting the drying springs by shallow wells (19 well with depths 150 -200 m and discharge rate 1000 m 3 /day/well) and cleaning up the flowing ones to enhance their productivities. Furthermore, keeping the present discharging rates at 490,000 m 3 / day with rearranging the wells immersed in the Nubian aquifer by substituting the stop flowing wells by new deep wells (depths >1000 m , Zone C) in the eastern and southeastern parts (Baraka and Qarawen areas) where the aquifer thickness and characteristics are reliable. No more wells immersed in the first productive zone since it forms the main recharging source of the Post Nubian aquifer. In addition, reducing the discharging rates of some high productive wells to about 5000 m 3 /day/well especially in Abu Nuss area.
The model was run under these new stresses and heads were estimated for the next 20 years for the upper and lower layer. The head distribution maps of the Post Nubian aquifer (upper layer) showed relatively gradual decrease in the groundwater levels indicated by the slow migration of contour lines southward (Figure 14), i.e. contour line 85 slowly moves southward to about 10 Km south Ain Shagy at the end of the simulation time (20 years). This behavior may be attributed to the new wells imposed in the aquifer that maintains the drawdown rates closely related to that of the first scenario. Drawdowns are expected to range from 5 m to 8.6 m ( Table 2). The maximum drawdowns were recorded in Ain shaggy and Tenin reflecting the new pumping activities in such areas. Groundwater flow directions were not changed through this scenario ( Figure 14).
On the other hand, head distribution maps of the Nubian Sandstone aquifer in this scenario (Layer 2) showed relatively small drawdowns compared with the declines occurred in the first scenario ( Figure 15). The maximum drawdowns were recorded at Abu Nuss area as 35.8 m while the minimum drawdowns were recorded at Abu Minqar area as 5.1 m ( Table 2). Depression cones formed at Abu Nuss and Marzouq areas enlarged with time progress but with smaller rate than the first scenario. These drawdowns are direct result to the number and rate of well discharges at these areas. The drawdowns in Baraka and Qarawen areas showed small difference compared with the first scenario in spite of the increase of the well numbers. This may be attributed to the control of the well discharges and good well distribution. Drawdowns recorded about 12 m at the end of the simulation time.
The drawdowns in Post Nubian and Nubian Sandstone aquifers showed continuous relationship ( Figure 16). Three meters head decline in the Nubian Sandstone aquifer is followed by decline in the Post Nubian by about 1 m. The drawdown rate in the Post Nubian aquifer through the simulation time (20 years) not exceed 7 m and this is an acceptable value regarding the proposed extraction quantity (19000 m 3 /day) which represent about 4% of the total water exploited from the Oasis (Figures 15 and 16).

Conclusion and recommendations
According to the study, a real danger of either dewatering or increasing the water depths to uneconomic lifting depths faces both the shallow and deep aquifers in Farafra Oasis. Many springs dried out and the discharging rate of the flowing ones intensively decreased, i.e Ain Tenin discharging rate decreased from 156 m 3    • Cleaning up the flowing springs to enhance their productivities.
• Keeping the present discharging rates at 490,000 m 3 /day with rearranging the wells drilled in the Nubian aquifer by substituting the stop flowing wells by new deep ones (depths >1000 m, Zone C) in the eastern and southeastern parts (Baraka and Qarawen areas) where the aquifer thickness and characteristics are reliable.
• No more drilling of wells in the first productive zone since it forms the main recharging source of the Post Nubian aquifer.
• Reducing the discharging rates of some high productive wells to about 5000 m 3 /day/well especially in Abu Nuss area.