Conceptual and Mathematical Modeling of a Coastal Aquifer in Eastern Delta of R. Nestos (N. Greece)

: In this paper, the development of the conceptual and groundwater ﬂow model for the coastal aquifer system of the alluvial plain of River Nestos (N. Greece), that suffers from seawater intrusion due to over-pumping for irrigation, is analyzed. The study area is a typical semi-arid hydrogeologic environment, composed of a multi-layer granular aquifers that covers the eastern coastal delta system of R. Nestos. This study demonstrates the results of a series of ﬁeld measurements (such as geophysical surveys, hydrochemical and isotopical measurements, hydro-meteorological data, land use, irrigation schemes) that were conducted during the period 2009 to 2014. The synthesis of the above resulted in the development of the conceptual model for this aquifer system, that formed the basis for the application of the mathematical model for simulating groundwater ﬂow. The mathematical modeling was achieved using the ﬁnite difference method after the application of the USGS code MODFLOW-2005.


Introduction
Coastal aquifer systems are considered a water resource reservoir of ample importance, since almost 50% of global population lives within a zone of several kilometers away from the coast and therefore depend largely on them for their water supply. However, these environmental systems are considered hydrologically sensitive as they are vulnerable to salinization due to seawater intrusion, as a result of overexploitation. The aforementioned argument is proved by numerous case studies which are reported in scientific literature [1][2][3][4][5][6][7][8], proving that ecosystem services of coastal aquifers have become limited. The above problem becomes even more pronounced in semi-arid and arid climatic zones, such as the conditions within the Circum-Mediterranean , where limited groundwater resources may pose serious socio-economic impacts.
Based on the above, it can be realized that sustainable use of groundwater is considered a key issue for the optimized management of water resources of these systems, which in turn relies solely on extensive monitoring of hydrologic processes. The latter, is usually conducted within a model-driven approach, as it usually provides a mathematical model with input data in order to predict the future response of the aquifer under different optimizations scenarios.
The investigation of such complex physical systems usually demands a series of field measurements for the development of a reliable conceptual model, which in turn will result in a reliable mathematical model. It is also argued that the more accurate the conceptual model, the less efforts are needed for the application of the mathematical model. Accurate conceptual model is directly related to a thorough investigation of the natural system involving subsurface investigations (groundwater well logging, geophysical investigations, pumping tests), hydrological (groundwater level monitoring, surface and unsaturated zone studies) and hydrochemical measurements (groundwater sampling and analyses). This paper describes the use of a complete set of field hydrogeological investigations for the analysis of a coastal aquifer system in semi-arid conditions, in order to form the basis for the development of a conceptual aquifer model. Then, based on the above, the research involved the application of a mathematical groundwater flow simulation model.

Location of the Study Area
The research area is surrounded by Nestos river on the west and by the Aegean Sea on the south (Figure 1). The broader region morphology, on the south and south-east, is identified by low relief and some shallow pits, thus creating minor ponds, that are fragment of the eastern Nestos River delta wetland. The research area is mainly coved by clay, sandy clay, sand and in some cases pebbles ( Figure 2). Nestos River stands dominant as the most significant surface water body of this region. It is situated on the west of the research area, and at the same time, it is considered to be one of the major recharge source of freshwater for the inshore unconfined aquifer of the research area [46][47][48]. Laspias stream is situated along the eastern boundaries of the research area where the degraded industrial and sewage treatment effluents are discharged. The composite image of Figure 1, shows: (i) the monitoring network of groundwater wells of the study area, (ii) the location of geoelectrical sounding measurements sites along the boundary conditions in conjunction with the resulted geological sections, (iii) axes of geological material graduation, (iv) locations of the upper layer of sand-gravel material (thickness: 5-20 m) and locations of upper layer of clay-sand material (thickness: 15-20 m).
The research area is part of the National Park of Eastern Macedonia-Thrace, that was founded in 2008 and included the already protected areas of Nestos Delta, Vistonidas Lagoon, Ismarida Lagoon and their wider area. The study area ( Figure 2) faces serious environmental degradation due to: − The annual decrease of groundwater level of phreatic aquifers due to reduction of natural recharge with fresh water. − The seawater intrusion, occurring on surface water bodies as well as on groundwater aquifers. − The quality degradation of fresh water from different sources such as: applied fertilizers and pesticides, irrigation return flows and disposal of municipal and industrial waste. − The use of wetlands for grazing and fishing. − The rural development of the coastal zone.

Geological Setting
The research area is situated in Holocene delta environment deposits with a thickness of some tens of meters built by River Nestos and sub-streams. These deposits consist of alternate clay, sand and silt layering, have resulted in a broad range of formational and depositional tasks, that produced a very diverse geological domain. In addition, due to delta marshes the existence of organic clay at some spot is identified. The advancement of the eastern Nestos delta in conjunction with flooding conditions has been contributory in shaping low potential aquifers in the research area.
Recent as well as previous [46] field investigations were conducted (i.e., installation of monitoring wells, geophysical measurements, borehole loggings) in order to analyze the geological setting of the study area. The aforementioned studies showed that within the northern part of the research area, alternate clay-sand layering expands downwards down to a depth of approximately 30 m. Next, a marly layer of about 50 m thickness exists, while at a depth lower than 80 m, similar clay-sand layers are also present [46].
The covered area of east delta plain is approximately 176.4 km 2 , and only 60% or about 106.63 km 2 are cultivated. In addition, the coastal saline zone includes uncultivated regions with a total coverage of about 45 km 2 . Taking into account the respective hydrological figures of Nestos Delta flatland (infiltration approximately 15% and mean yearly rainfall from 1965 to 1996: 546.9 mm), the annual infiltration was evaluated at about 10.8 × 10 6 m 3 for the potential cultivated area of 131.4 km 2 (176.4 − 45.0 = 131.4). The total irrigated area expands up to 89.89 km 2 , and the 35.0 km 2 of which cover the watering needs from the nearby River Nestos. The remaining areas cover watering necessities by pumping groundwater. Bearing in mind the variety and the proportions of the cultivations in the study area, the yearly water utilization was reckoned at about 27 × 10 6 m 3 . Considering the values of the above fundamental factors of groundwater hydrologic equilibrium (infiltration and water utilization), a shortage of approximate 16 × 10 6 m 3 water was estimated [46]. A minor part of this deficiency is restored by Nestos riverbed infiltration, and the rates rely partly on the length of the old smaller and buried riverbeds and partly on the distance from the actual active riverbed.
In the research area ten (10) geoelectrical soundings were conducted, with the main target being to identify the stratigraphy of the region ( Figure 1). This task was executed as stated by the Schlumberger method for measuring soil resistivity. The maximum electrodes spread was 800 m (AB/2 = 400 m) for G1 and G5 cases and 1000 m (AB/2 = 500 m) for all the others. The survey spots were determined so that they can cover parts of the research area where the stratigraphy was not well identified. Thus, the survey spots are situated across River Nestos (3 surveys), across stream Laspias (3 surveys) and finally across the northern border of the research area, about 9.5 km from the seashore (4 surveys).
The results of the geoelectrical soundings correlate directly with stratigraphy data. The G1 survey was the only case in which a clay-sand layer remained undetected by the sounding. It is worth mentioning that any interpretation of a geoelectrical sounding must take into account surface geology data as well as borehole intersections.
High conductivity values during the geophysical investigations within that part of the study area, are attributed to percolation of saline surface water that is originated from Laspias stream (as described in more details §6). The results of the G9 and G10 surveys are in cohesion with all the other researches and surveys in this region. The G9 and G10 spots are situated next to the seashore, where the respective indicators and quality specifications of the groundwater have the worse values due to the intrusion of seawater. The G8 survey (located 5 km away from the seashore), came up with saline groundwater state, a somehow unexpected outcome. This finding should be examined further because there might be a potent affiliation with the geothermal field of the broader region that exists deeper below.
Finally, it must be mentioned that the presented resistivity values can only be validated in this research region, and so the span of hereby values is smaller, in contrast with the value range depicted in relevant international bibliography, for the geological formations of the research area.

Hydrometeorological Conditions
In Figures 3-5, hydrometeorological data are presented, for the period 1966 to 2013. Figure 6 shows the average monthly rainfall for the periods 1966-2012, 1985-2012, 1995-2012, 2000-2012. After analysis and processing of all the above diagrams of the preceding figures, the following outcomes obtained [52]: A peak value of precipitation is observed in November and December 1996 (the big rainfall event and the devastating floods of 30 November and 1 December, 1996). Characteristic is, for the season, the high precipitation value for the summer months     These outcomes are valuable considerations for the evaluation of natural groundwater recharge fluctuations of the research area, where alluvial aquifers are being formed.

Hydrogeological Setting
The area is located outside the zone with a large proportion of permeable material (Nestos paleogeographic axis) and therefore, the materials hosted expected to have small particle diameter (sands and clays).
Two hydrogeological systems are identified in the research area, which are situated within the quaternary coastal and alluvial deposits of the broader region [46]. The first hydrogeological system is a shallow one, and it consisted by of phreatic to semi-confined aquifers with a thickness of approximately 30 m. The main recharge of this hydrogeological system is mainly from precipitation and partly from nearby small streams and canals. The second deeper hydrogeological system of the study area, incorporates alternate confined aquifers with a thickness of approximately 200 m. Recharge of the second system derives partly from river Nestos infiltration through old buried river beds, and partly from lateral groundwater inflows originate from the neighboring hydrogeological basin of Vistonida lagoon.
Estimation of the groundwater hydraulic parameters of the research area aquifers after analyzing the pumping test data from 11 selected wells in the broader research area derived in values for [46]: (1) transmissivity (T), ranging from 4.0 × 10 −4 to 1.1 × 10 −2 m 2 /s, (2) storage coefficient (S), ranging from very low values to 10 −3 , rendering the aquifers of the research area as mainly confined westward, and in several sites as semi-confined. The major groundwater flow direction (Figures 7 and 8) is from northwest towards south with a little flow from northeast and central parts of the research area towards south. The first aquifer system is recharged mainly by rainfall and slightly by percolation from the upstream hills aquifers. The shallow aquifer system in the region over the past decade, was exploited by a large number of shallow wells (up to 1000) with small diameter, up to a depth of 15 m. Nowadays, few of them are remaining operational and they have been replaced by deeper wells up to 50 m depth. At the eastern part of the study area, the limited groundwater potential, the drainage constructions, as well as the particular geological structure have conducted to the groundwater and soil salinization in the region. When the area receives large amounts of precipitation, the effect of increasing of the shallow aquifer systems discharge is occurred, also taken into account that the infiltration rate for the region is over 15%. The second aquifer system (deep aquifer) is recharged to a great extent by groundwater percolation inflows through the Nestos paleogeographic axes (buried old stream beds), and partly by lateral groundwater inflows coming from the adjoining Lake Vistonida hydrogeological basin. This deeper aquifer system is being exploited by deep wells (about 50) the depth of which reaches 200 m.
In this study, the groundwater level fluctuations were monitored in the study area during the period 2007-2009 and in 2013 and piezometric maps of groundwater system of the study area were compiled for the upper phreatic aquifer system and the underlying deeper pressurized aquifer system of the study area, such as charts of groundwater level fluctuations. Based on the interpretation of these diagrams and maps (Figures 7 and 8 show the piezometric maps of the phreatic aquifer for the wet and dry periods of 2014) the following general conclusions were drawn: The main groundwater recharge source of the upper unconfined aquifer system occurs mainly from N-NW part of the study area from the River Nestos and old riverbeds, as well as from the local irrigation network and the northeastern part of the Laspias stream (dashed polylines in Figures 7 and 8).

Hydrochemical Setting
The groundwater samples from the research area were collected in polyethylene bottles following the procedure described by Rainwater and Thatcher [53]. The chemical analyzes were performed within 24 h, from the moment of sampling, in the Laboratory of Technical Geology of the Department of Civil Engineering of the Democritus University of Thrace. The measurements of cations and anions (with the exception of that of chlorine ions) were carried out using a spectrophotometer. The determination of chlorine ions (Cl − ) was performed by the titration method, as described in detail in the relevant manual of the American Public Health Association [54], with the gradual addition of 0.1 N AgNO 3 . All determinations concerning measurements of specific electrical conductivity values (µS/cm), temperature ( • C) and active acidity (pH), were performed in the field using portable instruments. Sodium (Na + ) and potassium (K + ) ions were determined using a flame photometer. In addition, during the pumping period July-August 2009, 30 samples were taken for the determination of deuterium (δD ‰) and oxygen-18 (δO-18 ‰) isotopes. The samples were taken in 100 mL polyethylene vials and packed in a suitable device-transfer refrigerator. These determinations were carried out by the Institute of Applied Geosciences of the Technische Universität Darmstadt, Germany.
To monitor the groundwater quality of the survey area, a network of appropriate 30 sampling wells was selected, in such a way in order to correspond to the underinvestigation aquifer and to cover the area as well as possible. The research periods basically include the pumping periods (July-August) of the years 2006, 2007, 2008, 2009, 2013 and 2014 where chemical analyzes were carried out.
The values of specific electrical conductivity EC of the samples ranging from 417 to 5100 µS/cm for 2009 and from 311 to 3750 µS/cm for 2013 ( Table 1). The higher values (>4000 µS/cm) are observed in SE part of the study area, while the majority of the wells presents specific electrical conductivity of around 1000 µS/cm. For the wells where low values of specific electrical conductivity are observed, are these which are parallel to the irrigation canal and recharged directly with good water quality from the river Nestos. The concentration values of Cl − ranging from 15.6 mg/L to 1320.2 mg/L for 2009 and from 1.8 mg/L to 819.8 mg/L in 2013 with around 30% of the values of water samples clustered close to the mean value (216.64 mg/L) or above it. Values less than 25 mg/L are apparent in the western part of the study area. Large concentrations of Cl − (>200 mg/L, higher than acceptable limit) make the local water undrinkable.
The concentrations of nitrate (NO 3 ), are ranging from 0 mg/L to 60 mg/L for 2009 and from 0 mg/L to 120 mg/L for 2013, with three (3) of the thirty (30) boreholes to be out of the acceptable limits (allowable: 50 mg/L), while the concentrations of nitrite (NO 2 ), ranging from 0.013 mg/L to 0.358 mg/L, with five (5) of the thirty (30) wells out of the acceptable limits (permissible limit: 0.1 mg/L).
Regarding the suitability of the water samples analyzed for irrigation, based on TDS (total dissolved solids), SAR (sodium adsorption ratio), the concentrations of sodium (Na), chlorine (Cl − ), bicarbonate ion (HCO 3 − ), the value criteria %E.sp (%Na, alkalinizing degree) and conductivity EC for 2009 and 2013 and evaluating the results according to ratings-rankings by Ayers [55], Richards [56] and Wilcox [57], generally all samples identified as problematic.
It is recognized that: − The majority of water samples classified into classification classes in which there is a high or very high salinity. This shows the seriousness of the problem of seawater intrusion into the groundwater aquifers of the area and thus makes it necessary to take measures to control salinity and crop only salinity resistant plants. − The largest proportion of the samples, the effect of salinity on fruiting classed with growing problems class, and there are also samples taken from specific wells in which the effect of salinity on fruiting classified as serious problems. It is obvious that the effect of salinity on groundwater, and consequently on irrigation is higher in areas closer to the sea. This fact is in direct consequence with the emerge of various problems in agriculture, land degradation and thus gradually reduction of agricultural yield. − According to the criterion of Wilcox water quality, 14 wells classified as 'good', 5 wells as 'accepted', 3 wells as 'doubtful' and 7 wells as 'excellent'. In part of the study area recharged by the river Nestos and the irrigation canal, waters are characterized from 'excellent' to 'good', while in the rest of the region are characterized from 'acceptable' to 'doubtful'.
In research area, the TDS (total dissolved solids) ranging from 267.3 mg/L to 3.269.1 mg/L in 2009 and from 199.3 mg/L to 2403.7 mg/L for 2013. Rates TDS > 3000 mg/L were measured in the water of just a well. Regarding the suitability of the analyzed groundwater sample for drinking from domestic animals the overall amount of the samples are characterized as appropriate for 'watering' all domesticated animals.
To avoid erroneous diagnosis of seawater intrusion due to a temporary increase in TDS, Revelle [58] proposed as a criterion of seawater penetration, the use of the ratio Cl/(CO 3 + HCO 3 ) in meq/L, which is also known as Revelle coefficient (R). Values of Revelle (R), which fluctuate from 1 to 10, indicates slight to moderate pollution from seawater intrusion, while values greater than 10 may be regarded as evidence of a serious pollution of seawater intrusion [51]. Based on Revelle values (R), it is concluded that the water of the studied area is characterized from 'good' to 'slightly contaminated' by seawater intrusion and 'serious polluted' in some places of the SW part.
From the analysis of Piper and Durov charts (Figures 10 and 11) as well as the correlation plot of SO 4 + HCO 3 vs. Cl (Figure 12), obtained from 30 sample wells, it is concluded that the majority of the samples are in the range of conservative mixture between freshwater and saltwater, and some of the samples show significant salinity problems due to possible seawater intrusion. Regarding salinity (affecting water for plant growth), as small to medium (SM) in the majority of samples, and as severe (S) in a few of them (4), 2.
Regarding the permeability (affecting the rate of water infiltration in soil), as minimal (N), small and average (SM), 3.
Regarding the specific ion toxicity (affecting sensitive plants), as minimal (N) to severe (S), with typical the surface irrigation where several samples of groundwater assessed as small to average degree (SM) in restriction on use (DRU).

Stable Isotopic Signatures of Groundwaters
The results for the isotope analysis reveal that the groundwater from of the study area is isotopically heavy. This can be seen from the plot of δD as a function of δ 18 O (Figure 13), where most points fall close to the global meteoric water line (GMWL: δ 2 H = 8δ 18 O + 10) with heavier isotopic values. The plot shows a pattern of two clusters of points; one very close to GMWL and the other away from it. The former represents isotope results for the wells located to the east of the study area close to the Aegean Sea. This implies that the aquifer in the part of the study area is subjected to recharge mainly from the sea by seawater intrusion. Ideally, oceanic water or seawater would fall below the meteoric water line but the observed result is due to dilution during mixing with freshwater in the aquifer. Deviation from the meteoric water line further left indicates the part of the recharge that originates from the freshwater of the Nestos River.
Considering the local meteoric water line for the Mediterranean region which the Mediterranean meteoric water line (MMWL: δ 2 H = 8δ 18 O + 22), all the groundwater samples fall below the MMWL meaning that the groundwater aquifer contains isotopically heavy water which can be attributed to seawater intrusion resulting from continued over pumping.

Description of the Numerical Code (Finite Difference Method)
The groundwater flow of the studied aquifer system was simulated by using the finite difference method, after the application of the USGS MODFLOW-2005 code. MODFLOW-2005 [62] simulates steady and non-steady flow in an irregularly shaped flow system in which aquifer layers can be confined, unconfined, or a combination of confined and unconfined ( Figure 14).
The three-dimensional movement of ground water of constant density through porous media is expressed by the partial-differential equation [63] as follows: → K xx , K yy , and K zz are values of hydraulic conductivity along the x, y, and z coordinate axes, which are assumed to be parallel to the major axes of hydraulic conductivity (L/T); → h is the potentiometric head (L); → W is a volumetric flux per unit volume representing sources and/or sinks of water, with W < 0.0 for flow out of the ground-water system, and W > 0.0 for flow into the system (T −1 ); → SS is the specific storage of the porous material (L −1 ); → t is time (T).

Spatial Discretization
The developed groundwater flow model expands in an area of 140 km 2 with the aquifer being represented by a 30 m thick layer (Figure 15

Boundary Conditions
The boundary conditions are used to establish the hydrological processes taking place in the investigated aquifer ( Figure 16). In this study, the MODFLOW-2005 packages used to represent the boundary conditions are:

•
General head boundary (GHB): This package is used to introduce the general groundwater flow conditions in the area. The hydraulic head used is 16 m.a.s.l. and is introduced in the northern boundary of the model area. The selected head is based on the piezometric maps mentioned before (Figures 7 and 8). The hydraulic conductance used in the simulation is 5 m 2 /d. • River (RIV): In the western part of the study area the Nestos river is simulated using the river package. The river depth for the in the reference period is 1.57 m., while the hydraulic conductance in the beginning of the simulation is 30 m 2 /d. The package is used for introducing the pumping for irrigation in the study area. The number of pumping wells in the area is very large (>2000), so for optimal data management the pumping was distributed across the model area ( Figure 17). The irrigation period is from June to August.

Hydraulic Parameters of the Model
The hydraulic parameters values are based on the pumping test analyses from various locations within the study area. The evaluation gave a range between 0.78 and 15.22 m/day, with an average value of 7.92 m/day. The same set of pumping tests is used for estimating the storativity coefficient, with extracted values varying between 1 × 10 −3 and 2.8 × 10 −2 with a mean value of 6.79 × 10 −3 .

Time Discretization
The time discretization of the model is based on the hydrogeological conditions that were defined during the field survey. The model has 13 stress periods (each corresponding to a calendar month) with daily time steps in each period, with the first one having steady state flow conditions and the rest transient (Table 2). Pumping was implemented in stress periods 6, 7 and 8, which correspond to June, July and August respectively. The observations used in the model are assigned in the end of the 4th and 10th stress periods.

Mathematical Model Results
After the groundwater flow model was calibrated the water budget showed that Nestos River is the main contributor to the model, with a total of 27.61 MCM/a, which accounts 50.61% of the total inflows in the model (Figure 18). Inflows from the northern (GHB) boundary are calculated 3.26 MCM/a (5.98%), while water coming from precipitation recharge is 16.11 MCM/a and irrigation return flows water is 7.57 MCM/a, which correspond to 29.53% and 13.88% respectively. A small amount of water (0.1 MCM/a) also comes from the boundary representing the sea (CHD). The vast majority of outflows from the system come from the water that is pumped for irrigation, which is 50.46 MCM/a, representing 90.07% of the total outflows. This leaves only 5.02 MCM/a flowing towards the sea (8.95%), while smaller amount of water outflow towards the Nestos River boundary (0.29 MCM/a, 0.52% of outflows) and the Laspias drain boundary (0.26 MCM/a, 0.46% of outflows). The model results were evaluated using the statistical method of root mean square error, RMSE, following the mathematical expression below: where z hi are measurements of heads h i at n locations and times.

Conclusions
The case study of the eastern delta plain of River Nestos (N. Greece) has played an important role to the economic development of the region. However, this was the reason for the quantitative and qualitative degradation of its phreatic and confined groundwater aquifer, after a long period of intense pumping conditions. Land reclamation through the construction of a network of drainage canals, aggravated the local piezometric conditions through lowering of the water table within the coastal part.
The research area is surrounded by Nestos river on the west and by the Aegean Sea on the south. The broader region morphology, on the south and south-east, is identified by low relief and some shallow pits, thus creating minor ponds, that are fragment of the eastern Nestos River delta wetland. The research area is mainly coved by clay, sandy clay, sand and in some cases pebbles. Nestos River stands dominant as the most significant surface water body of this region.
A large part of the wider area of the eastern Nestos delta is estimated as problematic because of its hydraulic linkage to the sea and also to the limited recharge of the groundwater and the overexploitation of the aquifers. These factors have contributed to the continuous salinization of groundwater and soils.
The research area is situated in holocene delta environment deposits with a thickness of some tens of meters built by River Nestos and sub-streams. These deposits consists of alternate clay, sand and silt layering, have resulted in a broad range of formational and depositional tasks, that produced a very diverse geological domain. In addition, due to delta marshes the existence of organic clay at some spot is identified. The advancement of the eastern Nestos delta in conjunction with flooding conditions has been contributory in shaping low potential aquifers in the research area.
In the research area ten (10) geoelectrical soundings were conducted, with the main target being to identify the stratigraphy of the region (Figure 1). This task was executed as stated by the Schlumberger method for measuring soil resistivity. The survey spots were determined so that they can cover parts of the research area where the stratigraphy was not well identified. The results of the geoelectrical soundings correlate directly with stratigraphy data.
Rainfall data from the area is analyzed in this manuscript, for the period 1966 to 2013, and the average monthly rainfall for the periods 1966-2012, 1985-2012, 1995-2012, 2000-2012 is presented. After processing of all the relevant diagrams of the preceding figures, useful outcomes were obtained.
Two hydrogeological systems are identified in the research area, which are situated within the quaternary coastal and alluvial deposits of the broader region [46]. The first hydrogeological system is consisted by of phreatic to semi-confined aquifers with a thickness of approximately 30 m. The main recharge of this hydrogeological system is mainly from precipitation and partly from nearby small streams and canals. The second deeper hydrogeological system of the study area, incorporates alternate confined aquifers with a thickness of approximately 200 m. Recharge of the second system derives partly from river Nestos infiltration through old buried river beds, and partly from lateral groundwater inflows originate from the neighboring hydrogeological basin of Vistonida lagoon.
Estimation of the groundwater hydraulic parameters of the research area aquifers after analyzing the pumping test data from 11 selected wells in the broader research area derived in values for [46]: (1) transmissivity (T), ranging from 4.0 × 10 −4 to 1.1 × 10 −2 m 2 /s, (2) storage coefficient (S), ranging from very low values to 10 −3 , rendering the aquifers of the research area as mainly confined westward, and in several sites as semi-confined.
Chemical analyzes and in-situ measurements were performed for the needs of this research. The examined parameters were: Ca 2+ , Mg 2+ , Na + , K + , HCO 3 − , Cl − , SO 4 2− , NO 3 − , NO 2 − , NH 4 + , PO 4 3− , SiO 2 , Fe 2+ , Mn 2+ , alkalinity P, alkalinity M, hardness (permanent, non-permanent, total), pH, electrical conductivity (EC). From the analysis of Piper and Durov charts as well as the correlation plot of SO 4 + HCO 3 vs. Cl, obtained from 30 sample wells, it is concluded that the majority of the samples are in the range of conservative mixture between freshwater and saltwater, and some of the samples show significant salinity problems due to possible seawater intrusion.
Based on SAR, TDS, EC values and concentrations of Na + , Cl − and HCO 3 − of water samples in 2009 and 2013, the restriction level on use was estimated (DRU, degree of restriction on use) [59,60], and characterized [52,61]: (i) regarding salinity (affecting water for plant growth), as small to medium (SM) in the majority of samples, and as severe (S) in a few of them (4), (ii) regarding the permeability (affecting the rate of water infiltration in soil), as minimal (N), small and average (SM), (iii) regarding the specific ion toxicity (affecting sensitive plants), as minimal (N) to severe (S), with typical the surface irrigation where several samples of groundwater assessed as small to average degree (SM) in restriction on use (DRU).
The results for the isotope analysis reveal that the groundwater from of the study area is isotopically heavy. This can be seen from the plot of δD as a function of δ 18 O (Figure 12), where most points fall close to the global meteoric water line (GMWL: δ 2 H = 8δ 18 O + 10) with heavier isotopic values. The plot shows a pattern of two clusters of points; one very close to GMWL and the other away from it. The former represents isotope results for the wells located to the east of the study area close to the Aegean Sea. This implies that the aquifer in the part of the study area is subjected to recharge mainly from the sea by seawater intrusion.
The synthesis of the above, in combination with existing data (geological, hydrological, meteorological etc.), provided the basis for the development of the conceptual model for the study area, prior the application of the mathematical model.
The study of the coastal aquifer of Nestos phreatic coastal aquifer involved also the development of a mathematical model for the simulation of groundwater flow. The applied numerical code is MODFLOW-2005, using the finite difference method and taking into account all the available hydrogeological information for the setting of the boundary conditions of the model. The model involved several packages (GHB, CHB, RIV, DRN, WEL, RCH, HOB) for the simulation of the relevant hydrological processes and stresses on the investigated aquifer, resulted in satisfactory results.
The application of the above, proved that the development of an accurate mathematical model strongly relies on a reliable conceptual model which in turn depends largely on a series of field and lab measurements usually involving time and economic costs. However, it was realized that the modeling efforts can seriously lowered when the conceptualization process is of high quality.