Arab J Geosci (2011) 4:645–653 DOI 10.1007/s12517-010-0261-8

ORIGINAL PAPER

Seawater intrusion in Syrian coastal aquifers, past, present and future, case study Khomine Abedelrahem Allow

Received: 26 May 2010 / Accepted: 6 December 2010 / Published online: 4 January 2011 # Saudi Society for Geosciences 2010

Abstract Overexploitation of shallow aquifers on the Syrian coast, north of Latakia (Damsarkho) for irrigation and tourism has caused an intrusion of seawater. The seawater intrusion into this aquifer has been presented by a three-dimensional finite element model using the FEFLOW numerical code. This conceptual model is based on field and laboratory data collected during the period 1966–2003. Meteoric infiltration and flows from the adjoining carbonate aquifer recharge the aquifer; natural outflow occurs through a diffuse flow into the sea; and artificial outflow occurs through intensive extraction of groundwater from wells. Water exchanges in the aquifer occur naturally (leakage) and artificially (multi-screened wells). The model was calibrated for transient conditions. The model helped establish that seawater intrusion is essentially due to withdrawals near the coast during the irrigation season and that it mainly occurs in the Damsarkho plain. The effects of hypothetical aquifer exploitation were assessed in terms of salt budget. Keywords Seawater intrusion . FEFLOW . Groundwater . Coastal aquifer . Syria

Introduction In many coastal areas, the growth of human settlements, together with the development of agricultural, industrial and touristic activities, has led to the overexploitation of the aquifers. Such overexploitation induces a rise in the K. A. Allow (*) Mineralogy, Petrology and Geochemistry Department, University of Szeged, 6722 Egyetem utca 1, Szeged, Hungary e-mail: [email protected]

freshwater–saltwater interface (seawater intrusion) and thus, a degradation of the chemical quality of the groundwater. Under natural conditions, the geometry of the saltwater wedge depends on the hydraulic properties of the aquifer, on the physical properties of the two fluids (Henry 1964a, b) and on the geometry of the aquifer (Abarca et al. 2007). This phenomenon is observed at the Syrian coast, north of Latakia city (Damsarkho plain). The coastal plain of the Damsarkho is currently experiencing an intrusion of seawater owing to an irrational exploitation of the aquifer through a series of hundreds of wells of different types and depths with different pumping rates. The first real study of this problem was carried out in 1997 by Abed Rabo who specified the location and extension of this intrusion by chemical analyses and by measuring the groundwater levels; the next study was done in 2003 by Abo Zakhem who used electrical conductivity and isotopic elements to investigate this phenomenon, and these two studies have explained and confirmed the existence of the intrusion. The deterioration of the groundwater quality is currently a limiting factor for local economic growth. Agriculture has either been completely abandoned or has been directed towards crops which can tolerate saltwater but are of inferior quality. The intrusion of seawater in coastal aquifers was first conceptualised independently by Badon-Ghijben (1889) and Herzberg (1901), assuming hydrostatic equilibrium, immiscible fluids and the existence of a sharp interface between freshwater and saltwater in a homogeneous unconfined aquifer. They found that the depth (Zs) of the freshwater–saltwater interface below sea level is given by: Zs ¼

rf hf rs  rf

where ρf is the density of the freshwater, ρs, the density of the saltwater and hf, the elevation of the groundwater table

646

above sea level. Due to molecular diffusion and hydrodynamic dispersion, fresh and salt water are actually miscible liquids; the contact between the two fluids is, therefore, a transition zone rather than a sharp interface. The situation is further complicated by the fact that the saltwater intrusion itself changes the fluid density, so this parameter varies in space and time as a function of changes in concentration, temperature and pressure in the fluid. Furthermore, the porous medium itself is usually stochastically heterogeneous. In order to properly reproduce the mechanism of saltwater encroachment, a variable density flow and transport modelling approach is therefore currently adopted. This work presents a conceptual and numerical model for simulating the hydrodynamics of the multi-aquifer system of the plain in Damsarkho. The model will be used for hydrochemical simulations and the correct management of local water resources. This aquifer system provides a good example of the situation in present-day coastal plains. The study area The coastal plain of Damsarkho is located in the north of Latakia, NW Syria, and consists of marine and alluvial sediments covering an area of about 40 km2. The study area is characterised by a Mediterranean climate with a wet winter (November–April) and a dry summer. The average rainfall varies from 800 to 1,000 mm/year (Selkhozpromexport 1979). The Damsarkho plain is characterised by the presence of a good aquifer, which consists of loose sands, gravelly sands, sandstone, limestone and sandy gravelly clay. These deposits are in direct connection with the seawater. The geologic map and the studied area are illustrated in Fig. 1. The average thickness of these aquifers is about 25 m. In the past, before the salinisation of the groundwater, the area was affected by a natural hydraulic gradient where the groundwater moved seawards, while the seawater only penetrated a little landward due to the hydrodynamic balance between the fresh and saline waters, and all that was before the intensive exploitation of groundwater (Abed Rabo 2000). The threat from seawater intrusion started at the beginning of the 1970s and at the beginning of the agricultural growth, when the interest in this problem did not exceed reports and general recommendations. With time, the demand for water has increased, and the number of unlicensed wells has also increased, from which they were irresponsibly extracting groundwater, leading to an exacerbation of the problem. All of that led to a decrease in the level of the groundwater and therefore a decrease in the abundance of fresh water that is flowing seaward. Measurements of the groundwater table indicated that the groundwater level in some wells is lower than the sea level and, of course, this has caused the seawater to move into the land through permeable formations. This gradual increase in salinity has

Arab J Geosci (2011) 4:645–653

been noted at the centre and at the limits of the plain, especially in the summer (Abed Rabo 2000). The phenomenon of salinisation has been neglected, which has led to a deterioration of 35% of the citrus trees and a destruction of large areas in which there used to be cultivated fruit trees and crops, and also in the deterioration of the physical and chemical properties of the soil because of the high concentrations of sodium ions in the wells (Abed Rabo 2000). The main source of water for irrigation is rainwater and water from wells. In the winter, the dependence is usually on rainwater for irrigation, but, in summer, it relies only on the groundwater. There are more than 700 wells in the Damsarkho plain; the depths of those wells are between 2–40 m, and the average pumping rate is about 15 L/s. The methods of irrigation in the region are the traditional and the submerged crater pit. The water that restores the Damsarkho aquifer comes primarily from the direct infiltration of rainfall and the deep percolation of excess irrigation and surface water, and partially from underground recharge from the eastern part of the aquifer. The aquifers in the Damsarkho plain are shallow groundwater reservoirs in which the movement of water is almost entirely seaward.

Conceptualisation Geometry and structure The main geometric–structural and hydrogeological characteristics of the Damsarkho multi-aquifer system of the coastal plain were reconstructed on the basis of general geologic data and 20 wells/boreholes. The system is comprised of five layers (sandy clay, limestone, sandy clay with gravels, marl and dolomite). The grid of the modelled area is illustrated in Fig. 2. These layers sometimes combine to form a single-layer aquifer. The geologic sections of these layers are illustrated in Fig. 3. The aquifer system is mostly recharged by the direct infiltration of precipitation falling in the modelled area. It is also recharged by the flow of water through contact with the alluvial deposits and cavernous limestone (lateral and vertical flows in the eastern sector). The natural outflow discharges into the sea. The intensive withdrawal from the aquifer through wells is the artificial outflow. Water exchange within the aquifer occurs due to both natural (leakage) and artificial causes (multi-screened wells). Piezometric surface and hydraulic properties The piezometric levels in 29 wells were measured monthly from February 1997 to February 1998 (Abed Rabo 2000). Note that, owing to continual groundwater extraction from the various wells in the area, it was only possible to measure

Arab J Geosci (2011) 4:645–653

647

Fig. 1 Geological map and study area (modified from Ponikarov 1966)

the dynamic water level, thus the hydraulic head data do not reflect the natural equilibrium conditions, especially in the summer months when the demand for water increases due to irrigation and tourism. In our study, we have used the measurements of water levels for May 1997 as the initial hydraulic head in the area because we have no previous data. For the recharge of the aquifer by infiltration, we have used the data of rainfall from 1966 to 2003 (unpublished data, Syrian Irrigation Ministry, 2003). To forecast the intrusion until 2020, the least values of rainfall were used to calculate the recharge, and as for the lateral flow through limestone, we have used the notes and data from the study of Abed Rabo (2000).

Numerical modelling Software The groundwater numerical flow model was developed using FEFLOW (finite element subsurface flow system)

working under both steady-state and transient conditions. For theoretical and practical information concerning use of the software and the solution of equations, readers can refer to the respective manuals (Diersch 1998). The finite element method was adopted for its flexibility and capacity to simulate complex geometric forms and to refine the nodal grid around points and/or single lines (observation points, coastline, etc.) Discretisation and boundary conditions and model parameters The aquifer system was discretised using a grid of triangular elements made up of nodes and elements covering an area of about 20 km². A higher degree of refinement was adopted for the control points (observation wells which were placed in the grid as fixed nodes representing their true position) and along the coast (important boundaries for the correct simulation of flow). The three-dimensional grid consists of five layers

648 Fig. 2 A three-dimensional grid of the area modelled

Fig. 3 Schematic cross-section of the modelled area through first and third axes of wells, showing the geological structure of the modelled area depending on the data of Abed Rabo (2000)

Arab J Geosci (2011) 4:645–653

Arab J Geosci (2011) 4:645–653 Table 1 The hydraulic properties of the five layers in the modelled area

649 Layer

Porosity (1)

Hydraulic conductivity (m/s)

Specific yield (1)

0.25 0.18 0.30 0.25 0.18

5.00e-4 1.50e-8 9.00e-4 4.05e-7 7.80e-8

0.07 0.05 0.25 0.25 0.05

Sandy clay Limestone Sand, clay and gravel Marl Dolomite

corresponding to the previously mentioned hydrogeological layers. The boundary conditions assigned to the numerical model were directly derived from the conceptual reconstruction of the aquifer system. A constant head value was specified to nodes along the coastline, where groundwater is in contact with the free surface of the sea. A hydraulic head at sea side was identified to grid nodes coinciding with the topographic surface by using the facility seawater surface. The flux was assigned at the eastern boundary to the limestone (Abed Rabo 2000); to determine this flux, we have used the power function parameter (i.e., the flux changes in time). Flow conditions were also used to determine inflows or outflows, which varied in time and space, and were thus used to simulate meteoric infiltration in the model (areal recharge). Effective precipitation was initially estimated according to Abed Rabo (2000) on a daily basis, depending on the simulation conditions (transient). The discharge condition was determined for the wells that exploited the aquifer, simulating a withdrawal for water extraction at individual nodes, which was estimated using statistical data and distributed accordingly to the density of distribution of the wells. The ability of FEFLOW to simulate the hydraulic connection between aquifer layers through multi-screened wells was used to determine the discharge rates. Due to the Fig. 4 Schematic cross-section of the modelled area showing the most important boundary used in the model

lack of reliable data on surface water and groundwater withdrawals for agricultural purposes, these quantities were estimated on the basis of the information provided from survey of the area (Abed Rabo 2000). We have also used the time power function facility to determine the discharge of wells and for the recharge rate from precipitation. The hydraulic properties (hydraulic conductivity, porosity and specific yield of the five layers of the model) were estimated using the ENVBRINS (groundwater data management programme); these parameters are shown in Table 1. The discharge of wells ranged between 3 and 20 L/s; the exploitation was during the dry months. The concentrations of solids have been assigned a value of 0 g/L at the eastern boundary, and 35 g/L at the sea side; the most important boundaries are shown in Fig. 4. Model input parameters In the model, the following parameters were used; the concentration and density of freshwater are 0 g/L and 1,000 kg/m3, respectively; the concentration and density of seawater are 35 g/L and 1,025 kg/m3, respectively; the density ratio is 245*10e-4; and the longitudinal dispersivity and transverse dispersivity are 5 and 0.5 m, respectively.

650

Arab J Geosci (2011) 4:645–653

Fig. 5 Simulated water concentration in the Damsarkho plain (a 1976, b 1986, c 1998, d 2010, e 2020). Note that the intrusion happens through three axes, and the largest extent of the toe is in the southwestern and northwestern parts

Results and discussion Based on the available data (Abed Rabo 2000) and the first seawater intrusion present in Damsarkho plain (unpublished data, Syrian Irrigation Ministry, 2003), the presented time of the model was from 1966 to 2020. Seawater intrusion until 1998 After the running the model, we have found that the intrusion occurs through three axes; the axes of seawater intrusion are illustrated in Fig. 5. The first one is in the

southwestern part of the area; the second and third are in the northwestern part. The location of seawater intrusion is dependent on the place of toe, the geological structure through the axes of intrusion and the amount of extracted water and these conditions are available in the two previously mentioned locations. The model helped establish that water withdrawal near the coast during the irrigation season was the main cause of seawater intrusion, and the most important factor that determined the place and extent of the saline edge was the hydraulic conductivity of the formations. This model is still not completely reliable because of incomplete knowledge of the system and few

Arab J Geosci (2011) 4:645–653

opportunities to fully check the results. However, it can already provide important information on the general evolution of the system under different stress conditions. In particular, the model can be quantitatively used to assess the impact of variations in the amount of water withdrawn and the position of injection wells that could help remediate saltwater intrusion. For the validation of the model, we have compared it with a previous study, and the results were similar to the published data in 2003(Abou Zakhem and Hafez 2003); the published report has shown that the seawater intrusion in the area presented was through three axes, as our model has also shown. Seawater intrusion until 2020 The aim of the developed model was to predict the impacts of future withdrawals on water conditions until 2020, including saltwater intrusion in the Damsarkho plain. In order to run the model, the lowest average rainfall was used, assuming the effect of climate change; the other properties were kept constant. The effect of pumping beside the little recharge gained by the aquifers lead to further intrusion of saline water into the

Fig. 6 The depression in the Damsarkho plain in 2010 according to the model (Image from Google Earth)

651

subsurface aquifers. The predicted location of the saline toe in the area in 2020 is similar to those in the previous years, but it has a greater extent landwards. The 35 g/L concentration of the intruding water in 2010 has reached the first, 24th and 25th wells, while, in 2020, the 35 g/L concentration may be also noticed in the second and 26th wells. Depression in the area After extracting the data and information from the model for 2010, the depression of the water level has been obtained in the area using Surfer 8. The values of depression range between 0 and 1.8 m; the depression has happened due to extraction of water from the wells; here, it is necessary to repeat that the type of wells are shallow wells, where the depth ranges between 2 and 40 m, although the amount of pumping from these wells does not exceed 20 L/s. The contours of the depression are illustrated in Fig. 6. As already mentioned, the water from the wells has been used for irrigation and drinking; it is worth noting that the highest values of the depression were in the southwestern and northwestern parts of the Damsarkho plain, which were fully compatible with the location of the seawater intrusion into the area.

652

Arab J Geosci (2011) 4:645–653

Table 2 Computed salt quantity intrusion through sea side, the quantum of extracted salt from wells and the imbalance over the five periods of time Years

Inflow through seaside (ton)

1976 1986 1998 2010 2020

36,452 30,134 29,371 26,959 26,455

Outflow through seaside (ton)

Outflow through wells (ton)

Imbalance (ton)

−20499.6 −18024.1 −15872.8 −13791.5 −13360.2

−6.53173 −1584.69 −2068.25 −2469.38 −2745.61

15,299.43 10,525.71 11,429.95 10,698.12 10,349.19

Note that the imbalance has a positive value

Salt budget To have a view about the salinity and its evolution in this aquifer, we have extracted the amount of salt that entered through the sea side during the mentioned time period. The salinity imbalance was also extracted, and it had a positive value which means the salinity increased in the area. The salinity is shown in Table 2, which demonstrates that the salinity in the Damsarkho plain is gradually increasing. The increase in salinity through the axes of intrusion was different, as shown in detail in Fig. 7. The highest values of salinity were through the first and fourth axes of the wells, while the medium increase was through the third axis, and the low one was through the second axis of the wells. This

Fig. 7 The increasing of salinity through the axes of wells (I, II, III, IV). The greatest increase was through the first and fourth axes, which is due to the influence of the geological structure

difference is due to many factors such as depth of well, amount of extracted water, thickness and the degree of homogeneity of the aquifer. The amount of extracted water at the first and fourth axes was the largest among all other wells. In addition, the limestone which is the essential aquifer has the largest thickness through first axis.

Conclusion The deterioration of groundwater quality in the Damsarkho coastal aquifer is mainly the result of seawater intrusion and the local upcoming into the aquifer. The hydrogeological nature of the formation was helped by the presence of the

Arab J Geosci (2011) 4:645–653

good permeability of the formation that was in contact with the sea. This intrusion is due to the lowering of the freshwater level and in relation to excessive groundwater abstraction. FEFLOW is a useful tool for predicting the transition zone of the the Damsarkho coastal aquifer. The result of the model will be used for identifying the location and number of monitoring sites of seawater intrusion, for predicting the future expansion of the intrusion if excessive abstraction continues, and for determining the best solutions to prevent the intrusion. By using FEFLOW, we could determine the value of the salt budget and, most importantly, we could specify the imbalance of the aquifer, which has a positive imbalance favouring salt. In addition, the seawater intrusion according to the model occurred through three axes where the geological structure and hydrogeologic properties were the factors that favoured the intrusion in addition to the extreme abstraction. The coastal aquifer should be continuously monitored in order to determine the progress of the saline intrusion and carefully check the model's predictions. This can be accomplished by establishing a detection monitoring network using conventional observation wells for groundwater sampling, water level measuring and to specify the salinity of groundwater. The prevention of further intrusion of seawater must be started immediately before the hazards of saline water expand to most of the area, and that can either be by injection wells or by building subsurface barriers along the beach, especially where the toe increases.

653

References Abarca E, Carrera J, Sanchez-vila X, Voss CI (2007) Quasi horizontal circulation cells in 3D seawater intrusion. J Hydrology 339:118– 129 Abed Rabo R (2000) Study of seawater intrusion into fresh groundwater in the coastal area. In: MS Thesis, civil engineering. Damascus University, Damascus Syria Abou Zakhem B, Hafez R (2003) Environmental isotopes study of seawater intrusion in the coastal aquifer (Syria), AECS-G\FRSR 285 (internal report). Damascus, Syria Badon-Ghijben W (1889) Note in Verband met de Voorgenomen Put boring Nabij Amsterdam {Notes on the probable results of well drilling near Amsterdam}.Tijdschrift van het Koninklijk Instituut van Ingenieurs, The Haag, 8–22 Diersch H-JG (1998) A shock-capturing finite-element technique for unsaturated-saturated flow and transport problems. Computational Methods in Water Resources XII, Vol. 1 Computational Methods in Contamination and Remediation of Water Resources 207–214 Henry HR (1964a) Effects of dispersion on salt encroachment in coastal aquifers. Water-Supply Paper 1613-C, Sea Water in Coastal Aquifers: C71-84, U.S. Geological Survey Henry HR (1964b) Interfaces between salt water and fresh water in coastal aquifers. Water-Supply Paper 1613-C, Sea Water in Coastal Aquifers: C35-70, U.S. Geological Survey Herzberg A (1901) Die Wasserversorgung einiger Nordseebaeder {The water supply of selected North Sea towns}. ZF Gasbeleucht Wasserversorg 44:815–819 Ponikarov VP (ed) (1966) Explanatory notes of the geological map of Syria, scale 1200000, sheets 1-37 XIX,1-36 XXIV Selkhozpromexport (1979) Hydrogeological and hydrological surveys and investigation in four areas of Syrian Arab republic, vols. 1 and 2. Georgian State Institute for Design of Water Resources Development Projects, Tbilisi, USSR (1979) p. 173