Intensively exploited Mediterranean aquifers : resilience and proximity to critical points of seawater intrusion

Introduction Conclusions References


Introduction
The generally slow rate of water movement and storage are innate advantages of groundwater resources in terms of resilience to climate variations and relative to surface water resources.Groundwater is often of potable quality, and thus does not require expensive treatment, permitting scaled development upon demand and use of infrastructure that is normally of lower cost than that of surface water supplies (Taylor et al., 2009).However, aquifers in general and especially coastal aquifers are often exploited so intensively that their natural hydrological regime is strongly disturbed and may be thrown out of balance (Custodio, 2010).In particular, sea intrusion changes are highly intrusion into a coastal aquifer may occur in response to even small sea level and/or groundwater management changes (Mazi et al., 2013).The possibility of such acceleration into full sea intrusion is a major concern for the resilience and sustainability of coastal populations that depend on groundwater for their water supply.Several aquifers along the densely populated Mediterranean coasts are already suffering sea intrusion.Acceleration of this phenomenon could also be particularly large here, because the Mediterranean region, and especially its semi-arid areas, is likely to be seriously affected by a decline in water resources (Kundzewicz and Döll, 2009).Human vulnerability to decrease of renewable groundwater resources induced by climate change, as projected by climate models up to 2055, is expected to be higher for the North African rim of the Mediterranean Sea (Döll, 2009).
In this paper, we study three coastal aquifers in the Southeastern Mediterranean region that are already suffering seawater intrusion, with the aim to determine how close these aquifers are to reaching critical seawater intrusion extents in relation to different possible groundwater pumping locations (Koussis et al., 2010a, b) and/or to tipping points for complete seawater intrusion in the aquifer (Mazi et al., 2013).The studied aquifers are (Fig. 1): (1) that underlying the Nile Delta (the Nile Delta Aquifer, referred to as NDA), which is of profound importance to Egypt, (2) the Israel Coastal Aquifer (ICA), which is considered to be the most important freshwater source in Israel, and (3) the Cyprus Akrotiri Aquifer (CAA), which is a major source of potable and irrigation water in Cyprus.These aquifers have been subject to extensive exploitation and modified hydrologic regime by human activities and constructions (e.g.dams), and are located in a sensitive region of the world, both in climatic and political respect.
We here review the salinization history and current status of these aquifers, and quantify their resilience/vulnerability to current and future exploitation schemes and/or hydroclimatically changed recharge conditions, in view of their respective tipping points and other critical extents of seawater intrusion.Figures

Back Close
Full 2 Materials and methods

Site descriptions
In this section we summarize the physical characteristics of the three regional coastal aquifers (Fig. 1) and present the manner of their exploitation and the resulting current and possible future development conditions.

The Nile Delta Aquifer (NDA)
The Nile Delta area (22 000 km 2 ) is among the largest deltas in the world.This area is important as the outflow area of one of the largest rivers of the world, and because it concentrates 45 % of the total population, 45 % of the arable land, 50 % of the industrial production, 40 % of the agricultural production and 60 % of the fish catch in Egypt (Mikhailova, 2001).
The Nile Delta has its apex near Cairo, approximately 200 km from the Mediterranean Sea, and its base at the Mediterranean coastline extending from Alexandria to Port Said, a distance of about 240 km (Fig. 2a).The area of the Nile Delta is bounded to the east by the Ismailiya and the Suez Canals, to the west by the El Nubariya Canal, and to the north by the Mediterranean Sea.The Damietta branch of the Nile, flowing to the NE, and the Rosetta branch, flowing to the NW, divide the NDA in a middle, an eastern, and a western part.
The average annual rainfall in the Nile Delta ranges from 180 mm at the coast to ca. 26 mm in the Cairo area, occurring mostly in the winter (Sherif and Singh, 1999).The alluvial NDA underlies the Delta plain and accounts for freshwater resources exceeding 400 km 3 (Mikhailova, 2001).The aquifer is recharged by infiltration from a network of irrigation conduits, from excess irrigation water and limited precipitation, and by flows from the Nile Valley aquifer.The main aquifer consists of unconsolidated coarse sands and gravels (with occasional clay lenses), having an estimated porosity varying between 25 and 40 % and thinning towards the east and west sides of the Delta.The Introduction

Conclusions References
Tables Figures

Back Close
Full near-surface deposits change in the seaward direction into impervious clays and silt.
The aquifer is mostly unconfined, but it is considered leaky and partially semi-confined close to the sea.The aquifer base is a clay aquiclude that is inclined towards the north with a slope of 0.3-0.4%.The depth of the aquifer at the Mediterranean coast can reach 1000 m, but near Cairo it is only about 200 m (Sherif, 1999).
Groundwater flows in a radial pattern from the apex of the Delta to the sea and its quality deteriorates towards the north mainly as a result of seawater intrusion due to over-pumping, but also due to leakage of domestic wastewater (septic tanks), discharges of untreated industrial wastewater (from industries located in the greater Cairo area) and fertilizers used in agriculture.The piezometric head field in the aquifer fluctuates mildly (±1 m) (Kashef, 1983) in an annual cycle.A total volume of about 2.4 km 3 (Mikhailova, 2001) is extracted annually from the aquifer for all uses.But an accurate balance for the Nile Delta cannot be calculated due to unreliable or missing data (Mikhailova, 2001), and only a few components can be estimated with some certainty.Sherif et al. (2012) report that, today, seawater intrusion has reached 100 km from the coast.

The Israel Coastal Aquifer (ICA)
The ICA (Fig. 2b) is considered the most important freshwater source among Israel's major water reservoirs; the next two are the Lake of Galilee and the Mountain aquifer.The ICA covers an area of 1900 km 2 , with a length of 126 km and an average width of 15 km (Assouline and Shavit, 2004).It is a seaward sloping aquifer of inclination 0.01 (1 % slope) (Yechieli et al., 2010) and its maximum thickness at the coast reaches 200 m, thinning to a few meters towards the eastern (land) boundary.
The aquifer belongs to the Kurkar Group, consisting of inter-layered sandstone, calcareous sandstone, siltstone and red loam, alternating with continental and marine clays that overlie impervious marine clays of the Saqiye Group of Pliostocene age (Yechieli et al., 2010).East of the shoreline, to a distance of 5-8 km, clay interlayers subdivide the aquifer into four sub-aquifers.Certain of the lower sub-aquifers 13821 Introduction

Conclusions References
Tables Figures

Back Close
Full are confined while the upper ones are phreatic (Yechieli et al., 2010).The general flow direction in the aquifer is from the east towards the Mediterranean Sea to the west.The average annual precipitation ranges from 550 mm yr −1 in the north, to 300 mm yr −1 in the south, while the infiltrated natural recharge is about one third to one half of the precipitation (Aberbach and Sellinger, 1967).
The chlorinity of groundwater in the aquifer did not exceed 100 mg L −1 Cl − during the 1930's (Kass et al., 2005).Since then, however, the aquifer has been under extensive exploitation, and during the 1960's groundwater withdrawal was up to 400×10 6 m 3 yr −1 (or twice the natural recharge of the aquifer); in 1962, about 480×10 6 m 3 yr −1 of groundwater were extracted from the aquifer (Aberbach and Sellinger, 1967).In some locations, over-exploitation of groundwater has caused the water level to drop by 15-20 m and seawater to intrude 2.5 km from the coast (Aberbach and Sellinger, 1967).Assouline and Shavit ( 2004) report monitoring data for groundwater quality that indicate a steadily increasing mean salinity; at 2004, groundwater salinity was around 200 ppm Cl − and according to some predictions may reach 300 ppm Cl − by 2020.
To enhance the water balance of the aquifer Israel practices artificial recharge, mostly in the winter, through wells and spreading grounds, with water conveyed annually by the National Water Carrier (NWC) from the Lake of Galilee (earlier also floodwater).Water brought by the NWC in 1973 was reported to be 80 × 10 6 m 3 yr −1 (Sellinger and Aberbach, 1973), corresponding to 40 mm yr −1 on average over the whole aquifer area; other reports estimate this contribution to be 100-150 × 10 6 m 3 yr −1 (Melloul and Zeitoun, 1999), corresponding to 50-80 mm yr −1 .Furthermore, 300 × 10 6 m 3 yr −1 of Introduction

Conclusions References
Tables Figures

The Cyprus Akrotiri Aquifer (CAA)
The CAA (Fig. 2c) extends over 40 km 2 in the southern edge of the Akrotiri basin (78 km 2 ), which forms the southernmost peninsula of Cyprus.The Akrotiri basin has a smooth and hilly relief, with average altitude over 200 m a.m.s.l. in the northern part of the basin.The Salt Lake, located in the middle of the peninsula, is a topographic low, with mean water level below sea level, receiving drainage from the Akrotiri basin (Mazi, 2000).
The Akrotiri basin and aquifer geology, hydrological budget, and water management scenarios and conditions have been studied extensively by Mazi (2000), Koussis (2001), Mazi et al. (2004a, b), and Koussis et al. (2010bKoussis et al. ( , 2012)).The alluvial aquifer consists of river deposits, gravels with layers of marls, and sands with boulders interbedded with lenses of silts and clays.The impermeable base of the aquifer (alternating marls, chalks, chalky marls and marly chalks) slopes towards the south at about 1.7 % and the saturated aquifer thickness ranges from 10 m (water table at 50 m a.m.s.l.) at its northern edge to more than 100 m near the Salt Lake in the south.The mean annual precipitation  over the basin ranges from over 500 mm at its northern end to 420 mm near the Salt Lake.
Originally, the CAA was replenished by rainfall, leakage from Kouris River in the adjacent basin, inflows from its northern boundary and return irrigation flows.The closure and operation of Kouris Dam in 1988 has caused a 66 %-decrease in the replenishment of the aquifer by Kouris flows (Balasha and Phedonos, 1992).In response to declining groundwater levels and attendant seawater intrusion in the aquifer, the Water Development Department of the Republic of Cyprus initiated artificial recharge.Artificial recharge was applied to spreading grounds and ponds in the main aquifer Introduction

Conclusions References
Tables Figures

Back Close
Full area, when water from external sources was available.Nowadays, the CAA is replenished by various supplementary sources: (a) releases from the Kouris Dam, (b) Kouris reservoir losses, and (c) artificial recharge (effluents from Limassol's wastewater treatment plant and water from the Garyllis and Yermasogia reservoirs, located outside the Akrotiri basin).Almost no artificial recharge took place during the drought years of 1998 and 1999, however, and groundwater pumping was then also reduced.Mean annual groundwater extraction from the CAA was 14.5 × 10 6 m 3 yr −1 during 1967-1977(Jakovides et al., 1982)), but fell to only 10 × 10 6 m 3 yr −1 in 2000 (Mazi et al., 2004).
Our study focuses on the area of Zakaki, west of Limassol (the second largest city of Cyprus).During past dry years, groundwater from this eastern part of the CAA was used as additional potable water source, with groundwater salinity at 500 m from the coast reaching several thousand ppm of total dissolved solids (TDS) in the early 1990's.As a consequence, after 1997, groundwater withdrawals were reduced by 25 %, to counteract deterioration of groundwater quality, with salinity then stabilizing at 1000 ppm TDS in boreholes at about 1000 m from the coast.

Aquifer similarities and differences
The aquifers studied here are all important freshwater sources for the local populations and for the sustainability and prosperity of their economies.The natural regimes of the aquifers, replenishment from rainfall, and inflows from neighbouring freshwater sources, have all been altered by human activities during the last 50 yr, in the NDA through the Assuan Dams, in the CAA through the Kouris Dam, and in the ICA and CAA also through artificial recharge.All three are further unconfined sloping aquifers, from which large quantities of fresh groundwater are extracted annualy.In the NDA, groundwater pumping is dispersed over the area of the whole Delta; for the ICA and the CAA, more concentrated pumping zones can be identified.
Among the three aquifers, the NDA area is the largest (22 000 km 2 ), its depth exposure to the sea is greatest (1000 m), the inclination of its base is lowest (0.3-0.4 %) and it is much more permeable (mean hydraulic conductivity K = 100 m d −1 ) than the Introduction

Conclusions References
Tables Figures

Back Close
Full other two.The intruding seawater front has here advanced many tens of kilometers inland from the coast in the last 50 yr, and this movement has accelerated in the last 20 yr, reaching today a distance of ∼ 100 km from the coast.Another important difference from the other two aquifers is further the absence of a serious attempt to manage the NDA, which is in contrast to the importance of this aquifer for the regional water supply; the lack of management is also manifested in the absence of credible data for calculation of an overall water balance for NDA.
For direct aquifer comparison, the ICA area is about 1900 km 2 , its depth exposure to the sea is about 200 m, the inclination of its base is 1 % and its mean hydraulic conductivity K = 30 m d −1 .The aquifer is recharged by precipitation and artificial recharge applied through wells and spreading grounds.No boundary inflows replenish the aquifer.
An extensive monitoring network is in place and measurements indicate that the seawater front reaches today to 2.5 km from the coast.Several studies have been performed to estimate the permissible withdrawal from the ICA.
The CAA is further relatively small (40 km 2 ), its exposure to the sea is 50 m, the inclination of its base is 1.7 % and its mean hydraulic conductivity is K = 28 m d −1 .The aquifer is recharged by precipitation and by freshwater inflows from its northern boundary, and also by artificial recharge.In CAA there is an established network for monitoring the groundwater level and quality.

Modeling approach
This study focuses on a first-order vulnerability assessment of seawater intrusion in the three studied aquifers.These aquifers have also been studied previously by several scientists using variable-density models such as 2-D-FED: NDA Sherif and Singh (1999), FEFLOW: NDA Sherif et al. (2012), ICA Yechieli et al. (2010), SUTRA: CAA (Prieto, 2005;Prieto et al., 2006).However, the present study is the first that addresses their Introduction

Conclusions References
Tables Figures

Back Close
Full vulnerability in direct relation to the high non-linearity of the seawater intrusion process with associated tipping points for complete aquifer intrusion (Mazi et al., 2013), in conjunction with consideration of critical seawater intrusion to different possible groundwater pumping locations (Koussis et al., 2010a, b), and comparatively so across the different aquifers.
Present-day personal computers are capable of solving the coupled nonlinear equations governing variable-density flow and transport.Indeed, this is often done for specific studies at particular sites, like those of the present aquifer cases.However, the high computational demands of detailed simulations are an obstacle, and the usefulness of such simulations has also been questioned on the grounds that their results often lack reliability due to the common paucity of data for adequately characterising the detailed hydrogeology (hydraulic parameters) and concentration field (Sanford and Pope, 2010).In addition, the boundary conditions and hydrologic forcing (recharge and pumping) are often not known over time sufficiently.Thus, in the absence of credible data, the benefits of using detailed and complex numerical variable-density models are undermined by the guesswork (or non-unique fitting) and simplifications needed for the estimation of multiple aquifer parameters and their variability in space and time.
In comparison, analytical solutions of steady-state sharp-interface flow of oneflowing fluid are evaluated with efficiency and are suitable tools for first-order, regionalscale assessments of coastal aquifer vulnerability to seawater intrusion.They are also useful for screening the impact of various groundwater pumping alternatives in management scenarios (Koussis et al., 2012;Mantoglou, 2003).After developing a sound conceptual model, properly schematising aquifer geometry, and carefully estimating dominant hydraulic parameters and forcing quantities, analytical solutions can be used for systematic analyses, elucidating the way in which various physical and management parameters impact the seawater intrusion behaviour of coastal aquifers, and determining likely intrusion responses to various forcing conditions (Werner and Simmons, 2009;Ferguson and Gleason, 2012;Mazi et al., 2013).Introduction

Conclusions References
Tables Figures

Back Close
Full From this perspective, we use here the generalized analytical model of Koussis et al. (2012) (see Appendix A for the mathematical details) for steady interface flow in sloping phreatic aquifers, taking also into account a recent correction that alleviates the hydrodynamic defect of Dupuit-type sharp-interface solutions in phreatic coastal aquifers related to ignoring the finite outflow gap (rather than idealized singular point on a sharp-interface cross-section) through which submarine groundwater discharge occurs (Koussis et al., 2013).The solution of Koussis et al. (2012) accounts for the generally present and usually hydraulically significant aquifer slope that has previously been ignored in analytical sharp-interface solutions of seawater intrusion (e.g.Strack, 1976).This discharge-potential model approximates the gravity-driven flow component and represents the aquifer geometry in a schematised yet realistic manner.To compare our results to those obtained with variable-density models, we assigned to the sharp interface the 50 %-salinity line (17 500 ppm) of the transition zone, as suggested by Reilly and Goodman (1987) and adopted generally.With the use of this generalised model and approach, Koussis et al. (2012) have previously identified important criteria useful in the management of an aquifer, and Mazi et al. (2013) have shown the high non-linearity of seawater intrusion and associated tipping points that should not be crossed in order to avoid rapid loss of control and complete sea intrusion in an aquifer.
The model considers an inclined (sin ϕ) unconfined homogeneous aquifer with depth H sea at the coast, and with length L and hydraulic conductivity K .The aquifer is further recharged uniformly at the rate r , in addition to receiving an inflow, q b (for prescribedflux boundary condition, FB), or having a constant-head boundary, h L , that determines inflow through its land boundary (for prescribed-head boundary condition, HB).Furthermore, in the case of inland FB condition, the model solution can account for groundwater pumping in a collector trough (line sink, or well gallery) at a distinct distance w from the coastline, that penetrates the aquifer completely and draws groundwater at the rate q w .In the case of inland HB condition, dispersed groundwater pumping can be accounted for by reducing the recharge rate r with the pumped groundwater amount.

Conclusions References
Tables Figures

Back Close
Full the location of the sea intrusion interface toe, T , and the important (Destouni et al., 2008) but typically uncertain (Prieto and Destouni, 2011) submarine discharge, q SD .See further Fig.A1 in Appendix A for more details on all key parameters, boundary conditions and solutions discussed here.

Analysis approach
For the present analysis of aquifer cases with inland FB condition, and with a distinct groundwater pumping location, w , and pumping rate q w , we introduce the term q norm = −(q o + q w )/(KL), which is the normalised remaining groundwater flow from upgradient of w after the pumping of q w , where q o = −r (L−l w )+q b is the total groundwater flow from up-gradient of the pumping location w (and with q b , q o and all other groundwater flow terms defined as negative in the seaward direction).In order to make calculation results for different such case studies directly comparable with each other, they are in the following section shown and discussed in terms of the normalised position of the sharp interface toe, T /L, as function of q norm .For the analysis of aquifer cases with inland HB condition, and dispersed groundwater pumping, resulting T /L is instead shown as function of the consistently normalised submarine discharge q SD /(KL), with q SD = q b −rL+q w = q b −r n L, where r n = (r +q w /L) is apparent net aquifer recharge after dispersed pumping of groundwater at rate q w .The resulting curves of T /L vs. q SD /(KL) or q norm are then referred to as the performance curves of seawater intrusion (Koussis et al., 2012).
For each case study, we have further for the analysis selected typical cross-sections that are parallel to the groundwater flow in the different aquifers.For the NDA, where the flow is radial, we have selected two such cross-sections, one along the Middle Nile Delta Aquifer (mNDA), in the north-south direction, and one in the East Nile Delta Aquifer (eNDA), following a northeast direction (Fig. 2a).For ICA, the cross-section has an east-west orientation (Fig. 2b), while in CAA the orientation is northwest to southeast (Fig. 2c).Parameters used in the cross-section modeling of each case study Introduction

Conclusions References
Tables Figures

Back Close
Full are summarized in Table 1, and the cross-section conceptualizations are shown in Fig. 3 (a-b for NDA, c for ICA, d for CAA).
For the NDA, the lengths of the two cross-sections are as follows: the mNDA section has a length of 175 km and the eNDA section has a length of 180 km (Fig. 3a and b).Using the surface of the Mediterranean Sea as datum, the piezometric head at the inland boundary is at 14 m for both sections.The hydraulic conductivities of the two sections were found by interpolating the curves in Fig. 7a  The slope of the aquifer base and the aquifer depth at the coast were determined for the mNDA profile-section as 0.3 % and 740 m, respectively, and for the eNDA profilesection as 0.4 % and 910 m, respectively (Fig. 3a and b).
In their recent paper, Sherif et al. (2012) performed simulations with FEFLOW for five horizontal layers (areal sections) of the NDA, reaching to the depth of 400 m; without explicitly reporting the aquifer recharge, they report the seawater intrusion to exceed 100 km (1000 ppm) from the coast.As no information is openly available for the NDA recharge, and in order to use the latest published data on the Nile Delta sea-intrusion, we have used the present analytical model to calculate the likely current location of the interface toe, T , by calibrating the net aquifer recharge, r n , in the definition of submarine groundwater discharge q SD .Through this calibration, the sharp interface toe for the mNDA was positioned at T = 100 km, according to the intersection of the constant slope aquifer base with the 50 %-salinity line of profile 1 in Fig. 10 of Sherif et al. (2012) (the 50 %-salinity line was interpolated from the 1000 ppm and 35 000 ppm lines in the salinity transition zone).The thus calibrated net recharge value r n was found to be about 10 mm yr −1 , incorporating the aquifer recharge from all possible (natural and artificial) sources minus the known magnitude of groundwater extractions from the mNDA.
Because the extraction locations are dispersed in the NDA, there is no single representative pumping location, w , that we could use in the modeling and the illustration of resulting intrusion performance curves for this aquifer case; performance curves for the NDA sections are therefore shown in terms of T /L as function of q SD /(KL).Full ther summarizes the modeling parameters for the NDA sections, with r n for the eNDA being assumed approximately equal to that calibrated for the mNDA based on results reported by Sherif et al. (2012).
For the dispersed present pumping in the NDA cross-sections, we calculated then first the interface toe position T for various exploitation conditions of the aquifer (at various times in the past), assuming that the head at the inland boundary remained unaltered (at today's levels) at 14 m.We further also performed a sensitivity analysis of the combined effects of calculated net recharge changes from former to current pumping rates, in conjunction with possible groundwater head changes at the inland boundary.The latter could be brought about by changing management conditions in the Nile River Valley aquifer, or, more generally, by a change in the aportionment of the Nile waters among neighbouring countries, and/or by climate change.
For the ICA, our cross-section conceptualization is based approximately on the work of Yechieli et al. (2010), but considers a 20 km long cross-section ending in a perfect cliff, with an aquifer depth at the coastline of 200 m, slope of the impermeable aquifer base of 1 %, hydraulic conductivity K = 30 m d −1 and no boundary inflows (Fig. 3c).The total aquifer recharge considered here is 240 mm yr −1 , with 200 mm yr −1 derived from natural recharge and 40 mm yr −1 from artificial recharge.We further consider a pumping rate of 285 × 10 6 m 3 yr −1 (3000 m 3 m −1 yr −1 ) through a fully penetrating trough at 3000 m from the coast.Our analytical modeling applied to current conditions yields then the sharp-interface toe position at 2.6 km, in agreement with the intrusion stated by Aberbach and Sellinger (1967).The complete set of the parameters used for the ICA modeling is further summarized in Table 1.
For the CAA, the one-layer cross-sectional conceptualization is shown in recharge, r , and the inflows through the aquifer's landside boundary, q b , were estimated by Mazi et al. (2004) from the hydrological balance of the Akrotiri basin for the period 1972-1992 at the following mean values: r = 92 mm yr −1 and q b = −549 m 2 yr −1 (m 3 per metre width).For the same period, Koussis et al. (2012) translated the annual groundwater withdrawals in the Zakaki area of 2×10 6 m 3 to the 2-D-equivalent pumping rate at a representative well gallery of q w = 500 m 3 m −1 yr −1 for use in the profile model of the aquifer.The full parameter set used in the CAA modeling is further summarized in Table 1.
For both the ICA and CAA cross-sections, with their more well-defined average present pumping locations, w , we have calculated and illustrate in the following their respective intrusion performance curves in terms of T /L as function of q norm = −(q o +q w )/(KL), i.e., the normalised remaining groundwater flow from up-gradient of w after the pumping of q w .We further analyzed the effects of different possible pumping locations on the intrusion performance curves.

Model result interpretation of critical intrusion points
In the ICA and CAA cases with inland FB condition (prescribed boundary inflow q b ) and distinct pumping location w , the associated analytical solution allows for evaluating whether a hydraulic divide (with maximum hydraulic head) is formed at a location div down-gradient of w .The ratio T /l div , which expresses the relation between the interface toe location and the divide location, depends then on and is inversely proportional to the submarine discharge, q SD = q b − rL + q w .As pumping increases, q SD decreases and div moves seawards, i.e., decreases.If a hydraulic divide is formed down-gradient of w , i.e., if div /l w ≤ 1, then the pumping will draw groundwater also from the downgradient region ( w − l div ) between the coast and w .A critical tipping point of intruding seawater that is just on the verge of invading the whole coastal aquifer will here be Introduction

Conclusions References
Tables Figures

Back Close
Full reached if the interface toe reaches the groundwater divide, i.e., if T /l div = 1, as identified by Mazi et al. (2013).
If no groundwater divide is formed down-gradient of the pumping location, i.e., if div /l w ≥ 1, then the condition T /l w = 1 constitutes another point of critical seawater intrusion.This condition means that the invading seawater reaches the distinct well location, w , so that the well is lost before the intruding seawater reaches the tipping point of complete aquifer intrusion.The critical intrusion point T /l w = 1 is reached if the total groundwater flow q o = −r (L − l w ) + q b from up-gradient of the pumping location w is too small to counteract seawater intrusion to w .The pumping of q w may further draw groundwater only from up-gradient, or from both up-gradient and down-gradient of w .When the value of the normalised remaining groundwater flow after pumping, q norm = −(q o + q w )/(KL), is positive, the pumped water is only from up-gradient; when q norm is negative, some of the pumped water is from down-gradient of w .The flow q norm cannot become negative if a divide is not formed down-gradient of the pumping location.If the pumping location is moved further up-gradient, the pumping of q w may draw water from both up-and down-gradient of w , such that a divide can form between the sea and w if sufficient submarine groundwater discharge q SD = q b − rL + q w remains after that pumping.However, just moving a pumping well farther inland does not by itself ensure that the well will not be invaded by intruding seawater; such protection can only be achieved if a hydraulic divide indeed is formed down-gradient of w , and if the tipping point T /l div = 1 of complete aquifer invasion is not reached by the toe extending to that divide.
In the case of the inland HB condition (prescribed head at the boundary), the solution allows evaluating whether a hydraulic divide is formed at a location div down-gradient of the inland boundary.If a hydraulic divide exists in the aquifer and the net recharge r n decreases (e.g.due to increased dispersed groundwater pumping), the submarine groundwater discharge q SD will decrease and the groundwater divide location div will move landward, such that it may then reach the end of the aquifer, so that div /L = 1.This means that a tipping point of the type T /l div = 1 for complete aquifer intrusion Introduction

Conclusions References
Tables Figures

Back Close
Full will not be reached in this case.However, if the hydraulic head at the inland boundary decreases (e.g., due to groundwater pumping, or climate-driven recharge decrease in the hydrological catchment up-gradient of the inland aquifer boundary), the divide location div will move seaward, allowing for the interface toe of sea intrusion to possibly meet it, and the critical tipping point of T /l div = 1 to be reached such that complete sea intrusion occurs.
In both cases of possible inland boundary conditions (FB and HB) there exists thus the possibility of reaching critical tipping points of complete sea intrusion in the coastal aquifer, of the type T /l div = 1.In the flux boundary (FB) case with distinct pumping location, w , the other type of critical sea intrusion point, in which the intruded seawater toe reaches the well location, i.e., T /l w = 1, may also occur.Which of the two types of critical sea intrusion points that may be reached first, and thus be limiting for aquifer resilience to seawater intrusion, depends on actual aquifer conditions and their changes, as quantified further in the following section.

Aquifer-specific calculation results
According to our calculations, prior to the 1950's, when pumping from the NDA was negligible, the interface toe location was at 14 km in the mNDA and at 23 km in the eNDA.In the 1990's (Sherif and Al-Rashed, 2001) the pumping from the aquifer was 1.92 billion m 3 yr −1 and, according to our calculations, the toe moved then to 26 km in the mNDA and to 40 km in the eNDA.Increased pumping after the year 2000 to 2.4 billion m 3 yr −1 (Mikhailova, 2001) caused the interface to move to its current position, i.e., 100 km from the coast in mNDA and 120 km from the coast in eNDA.The strong non-linearity of the intrusion performance curves implies that the increased pumping has brought the NDA to a marginal situation (Fig. 4).More specifically, the submarine discharge q SD is currently 4.2 m 3 m −1 d −1 for the eNDA, and 3.9 m 3 m −1 d −1 for the mNDA, with the associated groundwater divide and interface toe locations then being at ∼ 154 and 116 km, respectively, in the eNDA, and at ∼ 143 and 98 km, respectively, in the mNDA.Thus, the ratio of the toe to the divide 13833 Introduction

Conclusions References
Tables Figures

Back Close
Full location, T /l div , which is a key parameter for aquifer resilience/vulnerability to sea intrusion, is 0.755 for the eNDA and 0.684 for the mNDA, which explains the higher vulnerability of eNDA.As r n decreases, the hydraulic divide will reach the end of the aquifer, div = L, for r n = 6.5 mm yr −1 in the eNDA and for r n = 5.7 mm yr −1 in the mNDA.
The ratio T /l div will then be 0.69 and 0.64 in the eNDA and mNDA, respectively, showing that div moves faster than T under decreasing r n .Because of the landside HB condition, with piezometric head at 14 m, the aquifer will not be completely invaded by seawater even if net recharge r n diminishes to zero; under those conditions T ≈ 140 km in the eNDA and T ≈ 131 km in the mNDA (Eq.2, Appendix A).In other words, even for zero r n , the interface toe T cannot here intrude to the inland boundary at L = 180 km in the eNDA and L = 175 km in the mNDA, if the prescribed hydraulic head and the implied inflow there are sustained at current levels (Fig. 4).However, if water management practices and/or climate change in the hydrological catchment upstream of inland aquifer boundary would decrease the hydraulic head at that boundary, the tipping point T /l div = 1 of full sea intrusion could be reached, even if the net recharge would be maintained at current levels (inset, Fig. 4).
The sensitivity analysis of the combined effect of recharge and boundary head change for the eNDA and mNDA shows that, if net aquifer recharge and associated submarine groundwater discharge q SD had remained at their former greater values, the head at the inland boundary would not play such a big role for the vulnerability of the aquifer to seawater intrusion.Specifically, a change of ±10 % in the boundary head would only cause the interface toe to fluctuate by ±350 m at the most in the mNDA section and by −860 m to +420 m in the eNDA section (Fig. 5a and Table 2).However, under current conditions of net recharge limited to about 10 mm yr −1 , a +10 % head change causes the interface to retreat by 6.5 km in both sections, while a −10 % head change causes the interface to advance by 5.3 km in the mNDA and by 7.7 km in the eNDA (Fig. 5b and Table 2).In particular, a decrease in the inland boundary head of 22 % in the mNDA and of 17 % in the eNDA can bring the NDA to its tipping point of complete aquifer intrusion (Fig. 4, inset).Therefore, under current recharge conditions, Introduction

Conclusions References
Tables Figures

Back Close
Full it is critical for the resilience and sustainability of the Nile Delta aquifer to control the upgradient catchment hydrological conditions so that sufficient piezometric groundwater head is maintained at the inland aquifer boundary.
For the ICA and CAA cross-sections, resulting sea intrusion performance curves are shown in terms of normalized intrusion toe position T /L as function of the normalized remaining groundwater flow q norm = −(q o + q w )/(KL) (Fig. 6).In general, depending on which condition is the limiting one for each case, the resulting performance curves terminate either at the tipping point of complete aquifer intrusion T /l div = 1, or at the limit of well intrusion T /l w = 1.If the tipping point condition T /l div = 1 is the limiting one, Fig. 6 shows also the corresponding value of T /l w < 1.If the well intrusion condition T /l w = 1 is instead limiting, the performance curves are continued (simulating the condition of pumping seawater) as dotted lines from that limiting point until the tipping point of complete seawater intrusion in the aquifer.
For current pumping rates (corresponding to the q norm = −(q o + q w )/(KL) value of the red dashed lines, Fig. 6; with q w = 3000 m 3 m −1 yr −1 in ICA and q w = 500 m 3 m −1 yr −1 in CAA) and pumping locations (black curves; w /L = 0.15 ( w = 3 km) in ICA, Fig. 6a; w /L = 0.33 ( w = 1 km) in CAA, Fig. 6b), the well intrusion condition T /l w = 1 (black filled circle) is the limiting one for both aquifers.The ICA (Fig. 6a) is then closer to that limit than the CAA (Fig. 6b).Increasing the pumping rate q w at the same well location will decrease q norm = −(q o + q w )/(KL) and move the aquifer closer toward the well intrusion limit T /l w = 1 (along the black curves, Fig. 6).The risk of crossing that limit can be reduced by moving the pumping location further inland (shifting from the black to, e.g., the green or blue curves, Fig. 6) but the tipping point of complete aquifer intrusion T /l div = 1 (open square symbol, Fig. 6) may then become the limiting condition instead of the well intrusion limit.Table 3 summarizes for the different well location choices illustrated in Fig. 6 the highest possible groundwater pumping rate, (q w ) max (and other associated variable values), until a critical sea intrusion limit is crossed.
In the ICA simulations for pumping locations w = 6000 m and w = 8000 m, the maximum pumping rate (q w ) max = 3603 m 3 m −1 yr −1 is limited by the condition T = l div = Introduction

Conclusions References
Tables Figures

Back Close
Full 4987 m.In the CAA simulations for pumping locations w = 1400 m and w = 2000 m, the maximum pumping rate (q w ) max = 697.6 m 3 m −1 yr −1 is limited by the condition T = l div = 1385 m.We note then that moving pumping locations more "up-gradient" does not, without further specification, automically imply well protection from seawater intrusion and increase of exploitable water volume.What matters is the specific location of the hydraulic divide to be down-gradient of the pumping location (critical point of sea intrusion to well) and the interface toe not to reach the divide (tipping point of complete intrusion).In the specific ICA and CAA cases, there is then no gain at all, in either exploitable water volume or aquifer protection, of moving pumping locations further up-gradient than w = 6000 m in ICA and w = 1400 m in CAA.

Conclusions
Directly comparable performance curves have here been developed for analyzing the resilience/vulnerability of different coastal aquifers to seawater intrusion.Two different types of critical sea intrusion limits are then identified: the well intrusion limit T /l w = 1, and the complete intrusion limit T /l div = 1.The risk of crossing the first limit can be reduced by moving pumping locations further inland, but then the tipping point of complete aquifer intrusion may be crossed instead.The measure of moving pumping locations more inland in coastal aquifers will offer no gain in more exploitable water volume or improved aquifer protection from sea intrusion beyond the limiting location of the prevailing groundwater divide in a coastal aquifer.
The site-specific vulnerability assessments show that the advance of seawater currently seriously threatens the Nile Delta Aquifer (NDA).The inland boundary head must here be sustained at any cost, as even a 10 % decline will cause sea intrusion advancement by around 5-8 km and a 20 % decline may bring the NDA to its tipping point of complete aquifer intrusion.These indications call for careful further study for establishing reliable hydrologic conditions for the NDA.Introduction

Conclusions References
Tables Figures

Back Close
Full In the Israel Coastal Aquifer (ICA), the current pumping locations do not allow for increased groundwater abstractions.The maximum pumping could here be increased by about 20 % if the pumping location were moved more inland, at least at w = 5 km from the coast, thereby allowing a groundwater divide to form down-gradient.In the Cyprus Akrotiri Aquifer (CAA), current forcing conditions do not directly threaten groundwater sustainability by increased seawater intrusion.Maximum pumping could here be raised by up to around 28 % above the present rate.If pumping locations were here moved more inland (to about w = 1.4 km from the coast), the maximum pumping rate would increase by up to 40 % over the present rate.

One-dimensional model of interface flow in phreatic coastal aquifers on sloping base
Here, we summarise the analytical model of Koussis et al. (2012), referring to that paper for the theory leading to the equations used here.Consider, then, steady flow with a fresh-sea water interface in a shallow phreatic aquifer of length L and hydraulic conductivity K , resting on an impervious inclined base (angle ϕ against the horizontal) and recharged at the uniform rate r (see schematic Fig. A1).The axis follows the aquifer base, from the vertical projection on the aquifer base of the intersection of the sea level and the coast, at the depth H sea below the sea surface, where the datum is placed.The depth normal to the aquifer base is h, the z axis vertical, and the hydraulic head is φ = h cos ϕ + z b , with z b the bed elevation; dz b /dl = sin ϕ.The sea constitutes a constant-head boundary, where φ(l = 0) = H sea .Koussis et al. (2012) use a discharge potential for the gravity-driven flow in which the depth is approximated as a constant h o (i.e., the gravity-driven flow −Kh sin ϕ ≈ −K h o sin ϕ) and assume that the freshwater flows over a stagnant seawater wedge.Introduction

Conclusions References
Tables Figures

Back Close
Full  2013) have solved the governing equation subject to flux-control and head-control conditions (Werner and Simmons, 2009) referring to management schemes applied in the aquifers; here we adopt the terms flux-boundary condition and prescribed-head boundary condition to denote the landside boundary conditions in the aquifers.The related solutions developed by Koussis et al. (2012) and Mazi et al. (2013) apply here.
The equation governing the flow in the entire aquifer is d 2 Φ/dl 2 = −r cos ϕ±Σσ i (σ i = point sources/sinks), with the discharge potential Φ defined differently in the freshwater zone 1 and in the interface zone 2 (see Mazi et al., 2013).Integrating, and implementing the boundary conditions for prescribed flux or prescribed head at the inland boundary, yields the discharge potential solution for each boundary condition.These solutions, combined with the discharge potential at the sea intrusion toe, lead in the flux-control case (prescribed flux at the inland boundary) to the quadratic Eq. (A1) and in the headcontrol case (prescribed head at the inland boundary) to the cubic Eq. (A2) for the location of the sea intrusion toe T : where δ = (ρ s − ρ f )/ρ f , with ρ s and ρ f the respective salt-and freshwater densities, increases as H sea increases and is the smaller of the two positive roots (a negative root is rejected).Koussis et al. (2012) have solved the flux-control case (flux-boundary condition) in which the aquifer receives a constant boundary inflow q b = q(l = L), plus recharge r , and includes a pumped collector trough of infinitesimal width (line sink or well gallery), located at = l w , which penetrates the aquifer completely and draws groundwater at the rate q w .The flow that reaches the pumped collector trough from up-gradient is then Our interest is in the profile of the discharge potential between the sea and the well gallery, where interface flow occurs, and particularly in the location of the sea intrusion interface toe, T .This is given by the quadratic Eq. (A1) in which with the physically meaningful solution being the root with negative square root.Koussis et al. (2013) corrected the 1-D (Dupuit) interface models to also include a vertical gap ζ o through which submarine discharge occurs (Destouni and Prieto, 2003;Prieto andDestouni, 2005, 2011).The size of ζ o was determined by adapting analytical estimates of the outflow gap of 2-D potential-flow interface solutions for infinitely thick aquifers.This correction allows predicting a more realistic location of the interface toe on the aquifer base (shorter penetration of a deeper interface), after reducing the depth at the aquifer outflow section H sea by ζ o to [H sea ] = H sea − ζ o .The validity of this approach was verified through comparisons with results of variable-density models.We adopt this methodology here as well, to better reflect the true near-shore groundwater dynamics in the Mediterranean aquifers that we are examining.For notational simplicity, henceforth H sea stands for the corrected value, by ζ o , of the initial sea depth [H sea ] .
In aquifers with flux-boundary and pumping from a well-defined location w , if a hydraulic divide is formed between the interface toe and the pumping location, such that 13839 Introduction

Conclusions References
Tables Figures

Back Close
Full T ≤ l div ≤ l w , the divide location can be calculated from l div = l w − (q o + q w )/r , (A5) while the head at the divide is calculated from with the physically meaningful solution being the root with the positive square root.
In aquifers with prescribed-head boundary and without a distinct pumping location (but still possible dispersed pumping that may be represented by reduced recharge r ), if a hydraulic divide is formed inside the aquifer, its location can be calculated (Mazi et al., 2013) from where z b denotes bed elevation and the hydraulic head is given by  Full  Full  3. Variable quantification related to the results of Fig. 6 for the Cyprus Akrotiri Aquifer and the Israel Coastal Aquifer.Quantified variables include groundwater exploitation (current pumping rate q w , possible maximum pumping rate (q w ) max , pumping location w in relation to aquifer length L, intrusion toe position T /L and (l T /L) max for q w and (q w ) max pumping, respectively, and relative location div /L of a possible groundwater divide formed down-gradient of w .
Discussion Paper | Discussion Paper | Discussion Paper | Fig. 3d.The length of the aquifer is L = 3000 m, its hydraulic conductivity K = 28 m d −1 , the slope of the aquifer base is sin ϕ = 0.017, and the depth of the sea to the aquifer base at the coastline is H sea = 50 m; a representative pumped collector trough is located at w = 1000 m.Conditions are typical of the eastern area of the unconfined aquifer in the Akrotiri peninsula, the Zakaki area bordering the city of Lemessos.The natural Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Mazi et al. ( and h L is the depth of flow at the inland boundary.Only the root with the negative square root has physical meaning for Eq.(A1).For Eq. (A2), the value of the depth h o must be estimated through iteration, e.g. by iterating on h o = [h(l T ) + h L ]/2.The single meaningful root of Eq. (A2) lies on the branch of the function T (H sea ) on which T Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Discussion Paper | Table

Fig. 4 .Figure 5 .
Fig. 4. Sea intrusion performance curves for the Middle and East Nile Delta Aquifer crosssections.Symbols mark the calculated normalized intrusion toe position, T /L, as function of the normalized submarine groundwater discharge q SD /(KL) under constant boundary head conditions for present and past groundwater pumping and associated net recharge conditions (continuous lines in main graph and inset), for which the tipping point of full intrusion (open squares) cannot be reached.The inset shows comparative results (dashed lines) for variable head conditions at the boundary under constant net recharge conditions at current levels, for which the tipping point of full intrusion (open squares) can be reached.

Fig. 5 .Figure 6 .
Fig. 5. Sea intrusion performance curves for the Middle and East Nile Delta Aquifer crosssections for different groundwater head conditions at the inland boundary.Normalized intrusion toe position, T /L, is shown as function of the normalized submarine groundwater dischargeq SD /(KL) for: (a) former values of net recharge and associated q SD /(KL) (range of higher q SD /(KL) values in the x-axis scale), and (b) current values of net recharge and associated q SD /(KL) (range of lower q SD /(KL) values in the x-axis scale).

Fig. 6 .
Fig. 6.Sea intrusion performance curves for different normalized well locations in: (a) the Israel Coastal Aquifer (ICA) and (b) the Cyprus Akrotiri Aquifer (CAA).Normalized intrusion toe position, T /L, is shown as function of the normalized remaining groundwater flow q norm from the normalized pumping location w /L.Absolute w values for the curves shown here are then 3 km (black line; current condition), 6 km (green line) and 8 km (blue line) for ICA, and 1 km (black line; current condition), 1.4 km (green line) and 2 km (blue line) for CAA.Black dotted line: continued pumping of seawater.

Table 1 .
Characteristics and parameters of the conceptual model sections (Fig.3) of the three Mediterranean aquifers.

Table 2 .
Interface toe location, T , for different conditions of net recharge rate and groundwater head at the inland boundary for the Nile Delta Aquifer sections.