Abstract

Groundwater depletion poses a major threat to global groundwater resources with increasing trends due to natural and anthropogenic activities. This study presents a surface-groundwater framework for water resources modelling of ill-posed problems in hydrogeologically data-scarce areas. The proposed framework is based on the application of a conceptual water balance model and composed of surface hydrological (UTHBAL) and groundwater flow simulation with the integration of a Newton formulation of the MODFLOW-2005 code (MODFLOW-NWT) and PEST suite modules. The groundwater simulation includes a preprocessor tool for automated calibration and a post-processor tool for automated validation. The methodology was applied to a rural region of Central Greece, Lake Karla Basin, which is degraded due to groundwater resources overexploitation to cover irrigation water demands. The aquifer is modelled focusing on a precise simulation–validation procedure of the conceptual model. The groundwater model was calibrated with the calibration preprocessor tool for spatially distributed hydraulic conductivity with the pilot points method. The calibration process achieved satisfactory results as validated by the post-process analysis of observed and simulated water levels. The findings for the groundwater budget indicate that the groundwater system is still under intense pressure even though farming activity in recent years has turned to less water-intensive crops.

HIGHLIGHT

  • My research deals with the pilot points method on groundwater modelling in an area with scarce hydrogeologic data.

Graphical Abstract

Graphical Abstract
Graphical Abstract

INTRODUCTION

Groundwater resources are vital in providing water for domestic, agricultural, and industrial use globally. Over two billion people rely on groundwater as their primary water source, while half or more of the irrigation water for agriculture is supplied from underground sources (Pisinaras et al. 2013; Famiglietti 2014). However, despite its critical importance to global water supplies, inadequate management practices are applied for monitoring groundwater resources. Furthermore, it is crucial to estimate the pressure on surface and groundwater systems due to groundwater depletion. Climate change and related changes on the components of the hydrological cycle complicate the task of sustaining groundwater supplies for the future. Hence, about 20% of the world's aquifers systems are estimated to be under over exploitation status and with increasing trends due to the climate change phenomenon and population growth (Upton et al. 2020).

Understanding why groundwater depletion has occurred may be the first step towards sustainable groundwater resources management. Numerical simulation models are effective and valuable tools for computing water budgets of regional groundwater systems and the planning and policy making for groundwater resources management (Konikow & Kendy 2005). The analysis of groundwater system through simulation is the main process before a researcher/modeller deals with the management of groundwater resources (Doherty & Hunt 2010). The simulation problems in groundwater hydrology are part of the complex environmental issues that are difficult to solve. These problems in groundwater modelling become extremely difficult when in addition to the lack of field data problems also meet a lack of quality data (Singh 2014).

One of the most important issues for water resources management is developing strategies for groundwater modelling that are adaptable to data scarcity and field measurements and so begins the difficulty of the inverse problem. According to Zhou and her associates (2014) many approaches have been proposed to resolve the inverse problem in hydrogeological modelling. Inverse methods are categorized as direct and indirect ones. The inverse problem in hydrogeological fields is difficult to solve because it is ill posed and the calibration parameters are more than the observations (Anderson et al. 2015). Common indirect methods used to solve the inverse problem are the geostatistical approach, the maximum likelihood method where zonation is employed, and pilot points. The most advanced methods where the stochastic approach is introduced are the self-calibrated method, the Markov chain Monte Carlo method, the gradual deformation method, the probability perturbation method, and the ensemble Kalman filter. The geostatistical approach has several advantages, however it is not appropriate for use in the inverse model where the hydraulic conductivity variable poses spatial distributions with significant heterogeneity, forcing the transition from linearized to non-linearized approaches (Zhou et al. 2014).

Zonation assumes an abrupt spatial change in parameter values, whereas the pilot points method produces smoothly distributed parameters and flexibility compared to the zonation approach (Anderson et al. 2015; Rajabi et al. 2018). These strategies are particularly important in arid and semi-arid areas where access to data is poor and data collection is difficult, such as the Lake Karla aquifer, Greece. This study analyzes the groundwater resources status of the Lake Karla aquifer focusing on a precise simulation–validation procedure of the groundwater model. A surface hydrological simulation tool was applied to calculate the groundwater recharge while the MODFLOW with a Newton formulation was employed to simulate groundwater flow. The groundwater model was calibrated for spatial distributed hydraulic conductivity based on the pilot points method, since there are limited data for this parameter.

This paper presents a new framework for calibrating and validating the groundwater flow simulation applying two automated tools with the use of pilot points. The first relates to calibration post-process analysis and is called PLPROC and the second one deals with the post-process analysis of the results and is called mod2obs (Doherty 2016, 2018). Both tools are part of the PEST: Model-Independent Parameter Estimation and Uncertainty Analysis, an open-source, public-domain software suite (Doherty 2020a, 2020b). The PLPROC tool has been applied in the Precipice Sandstone aquifer system of the Surat Basin, Great Artesian Basin, in Australia, by Hayes et al. (2020). The overall purpose of this paper is to provide a useful framework in surface and groundwater modelling dealing with scarce hydrogeological information in rural areas.

Study area

The study area of about 500 km2 is located at Eastern Thessaly, Greece and occupies the plain zone of Lake Karla watershed (Figure 1). The mean annual values of precipitation, temperature and potential evapotranspiration of the region are 560 mm, 14 °C and 775 mm respectively (Mylopoulos & Sidiropoulos 2016). The grainy, phreatic aquifer consists of alluvial deposits (Quaternary formations). The thickness varies in places reaching its maximum value of 200 m at the centre of study area. The presence of Mavrovounio Mountain at the east, where low permeable bedrocks prevail at the borderline with the lowland creates a no flow boundary. At the west, the permeable formations of Chalkodonio Mountain (south-west) and of the plain (west north) allow a hydraulic connection of adjacent aquifers with Lake Karla aquifer. The main crops are wheat, cotton, alfalfa, maize, vegetables, and trees. Lake Karla's aquifer is a characteristic example of an agricultural area in Greek territory, where surface water systems are limited and have relatively small capacity, in contrast to the great irrigation water demands imposed by the intensive agricultural activity of the area. As a result, the largest part of the region's water needs for irrigation is covered by groundwater with the use of pumping wells, most of them being private wells (Mylopoulos & Sidiropoulos 2016).

Figure 1

The study area map indicating the geological formation, the irrigation zones, the hydraulic heads observation points and the hydraulic conductivity sample values locations.

Figure 1

The study area map indicating the geological formation, the irrigation zones, the hydraulic heads observation points and the hydraulic conductivity sample values locations.

These agricultural practices and the absence of collective irrigation networks result in pumping non-renewable groundwater resources leading to groundwater depletion. Agriculture is the primary income source for most of the local population and for decades agricultural production in the region was based on subsidies from the Common Agriculture Policy (CAP) of the European Union. A typical example is the cultivation of cotton, a high-water-demanding crop which, for years, is based on the provided subsidies without considering sustainable management of groundwater resources (Alamanos et al. 2020). Previous studies have shown that the aquifer is in over exploitation condition and effective restoration practices are needed for sustainable water resources management (Sidiropoulos et al. 2015, 2019).

METHODS AND APPLICATION

The methodology is described in Figure 2 with the use of a flowchart. The first step is running the UTHBAL surface hydrology tool (Loukas et al. 2007). One of the outcomes of UTHBAL is groundwater recharge that is used as input parameter in the MODFLOW_NWT code (Niswonger et al. 2011). The preprocessor PLPROC was used during the parameterization while the optimization process follows the PEST methodology. Because the inverse problem is an ill-posed problem the regularization mode is used to prevent overfitting. Finally, the objective function results are examined. If the objective function values show a clear decrease, then post-process analysis follows (mod2obs). When there is no observed decrease, we can conclude that the groundwater hydrology model is unstable and needs to be re-evaluated.

Figure 2

Flow chart of the methodology.

Figure 2

Flow chart of the methodology.

The most common statistical criteria for the performance of a simulation model/tool are the Mean Absolute Error (MAE), Root Mean Square Error (RMSE) and for transient conditions is the Nash–Sutcliffe Efficiency (Ef). An acceptable criterion for the successful completion of calibration is a RMSE deviance of 10% of the water levels' total range. MAE is more robust in estimating the effects of outlier values of the spatial distribution than the RMSE, and therefore it is usually smaller than RMSE. However, if the calibration procedure should be continued, depends on score of the Nash–Sutcliffe Efficiency (Ef) (Anderson et al. 2015).

Estimation of groundwater recharge and irrigation demands covered by the groundwater system

The surface hydrological model (UTHBAL), is a surface hydrological simulation tool that gives a representation of the rainfall–runoff interactions of a catchment area. UTHBAL uses the components of the hydrological cycle, monthly time series of precipitation (P), mean temperature (T) and potential evapotranspiration (Ep) as inputs. The surface hydrological simulation tool disaggregates the total precipitation into rainfall and snowfall. UTHBAL distinguishes the total watershed runoff (Qsim), between the surface runoff (SR), the interflow (MR), and the baseflow (Qg) utilizing a soil moisture mechanism (Nmoist) and finally estimates the groundwater recharge to the subsurface systems (Loukas et al. 2007). In this study, UTHBAL has been applied in semi-distributed mode because the aquifer system covers the lower zone of the Karla watershed (e.g. Sidiropoulos et al. 2019). The groundwater recharge rate varied for each stress period.

The annual crop data from 2001 to 2016 come from the Integrated Administration and Control System (IACS), which is maintained by the Greek Payment Authority of Common Agricultural Policy (CAP) Aid Schemes (OPEKEPE, online). The crop pattern changes annually and the irrigation demands were estimated for every month from 2001 to 2016 with Near Irrigation Requirement (NIR). Hence, the pumping rates varied for the dry/pumping period.

Table 1 presents the estimation of irrigation water demands covered by groundwater extraction for the period 2001–2016. Table 2 shows the pumping rates for each zone and the respective irrigation area. It is worth noting that although the reservoir exists at present, it does not participate in the current study. The former natural Lake Karla had been drained in the 1960s, causing serious environmental problems in the area. A reservoir for the restoration of Lake Karla has been recently constructed to reverse the bad environmental status. The Lake Karla reservoir project construction has been completed in 2018 (Sidiropoulos et al. 2021). Thus, its contribution to the groundwater system is not considered during the simulation.

Table 1

Irrigation demands of each zone covered by the groundwater system

YearsZones (hm3)
1234567
2001 25.23 10.19 16.36 2.40 19.41 173.17 120.79 
2002 18.81 7.46 11.63 2.50 14.16 99.42 81.44 
2003 19.31 7.94 12.48 3.55 14.46 105.36 81.79 
2004 19.63 8.02 13.09 6.55 14.64 99.14 79.55 
2005 18.67 7.68 14.36 7.12 16.04 109.54 87.32 
2006 20.59 9.12 12.11 5.14 13.75 123.95 84.92 
2007 21.54 9.53 14.26 6.76 16.41 72.40 103.28 
2008 13.40 10.03 13.96 5.01 16.99 149.72 105.03 
2009 18.58 8.86 13.11 4.73 16.67 128.99 89.61 
2010 21.99 10.26 15.28 6.09 18.92 150.59 105.66 
2011 19.68 9.35 12.58 4.76 18.63 127.68 90.31 
2012 18.23 9.44 14.07 5.04 6.17 117.62 97.35 
2013 17.93 8.49 13.57 4.07 13.51 151.78 77.50 
2014 20.15 7.12 13.65 4.56 15.89 134.00 93.16 
2015 17.10 6.62 14.32 3.81 13.59 113.81 78.34 
2016 19.93 8.49 17.05 4.55 16.31 136.55 92.87 
YearsZones (hm3)
1234567
2001 25.23 10.19 16.36 2.40 19.41 173.17 120.79 
2002 18.81 7.46 11.63 2.50 14.16 99.42 81.44 
2003 19.31 7.94 12.48 3.55 14.46 105.36 81.79 
2004 19.63 8.02 13.09 6.55 14.64 99.14 79.55 
2005 18.67 7.68 14.36 7.12 16.04 109.54 87.32 
2006 20.59 9.12 12.11 5.14 13.75 123.95 84.92 
2007 21.54 9.53 14.26 6.76 16.41 72.40 103.28 
2008 13.40 10.03 13.96 5.01 16.99 149.72 105.03 
2009 18.58 8.86 13.11 4.73 16.67 128.99 89.61 
2010 21.99 10.26 15.28 6.09 18.92 150.59 105.66 
2011 19.68 9.35 12.58 4.76 18.63 127.68 90.31 
2012 18.23 9.44 14.07 5.04 6.17 117.62 97.35 
2013 17.93 8.49 13.57 4.07 13.51 151.78 77.50 
2014 20.15 7.12 13.65 4.56 15.89 134.00 93.16 
2015 17.10 6.62 14.32 3.81 13.59 113.81 78.34 
2016 19.93 8.49 17.05 4.55 16.31 136.55 92.87 
Table 2

Pumping rates of each zone and the respective irrigation area in hectares

ZonesArea (ha)Q (m3/h)
3,600 60–120 
1,490 35–100 
4,900 70–207 
1,400 105–250 
3,000 55–160 
21,000 40–120 
14,400 35–190 
ZonesArea (ha)Q (m3/h)
3,600 60–120 
1,490 35–100 
4,900 70–207 
1,400 105–250 
3,000 55–160 
21,000 40–120 
14,400 35–190 

Groundwater simulation

Groundwater flow simulation is achieved with the use of the MODFLOW_NWT code through the commercial GMS graphic user interface (Aquaveo 2019). MODFLOW_NWT code is a Newton formulation for MODFLOW-2005 (Niswonger et al. 2011). A no flow boundary was used in the east part to simulate the low permeable conditions of the aquifer. While in the west, general-head boundary was used to simulate the hydraulic connection of the aquifer with the adjacent ones. Conductance was estimated by the comparison of the observation wells on either side of the border ranging from 15 to 40 m2/d. The aquifer is represented as a single layer and the resulting grid consists of 12,500 active cells with 200 m side each. Grid refinement is observed at the wells' locations. This study examines the groundwater resources status extending the time series to November 2018 and applying the MODFLOW_NWT code instead of MODFLOW 2000. A transient groundwater flow model was developed and simulated from 1987 to 2018 with monthly changing irrigation needs from 2001 to 2016 also taking under consideration the annually different crop pattern application. Additionally, the groundwater flow was simulated with varying timesteps in accordance with the exact dates of water level observations. MODFLOW_NWT was chosen because the implemented method does not deactivate dry cells, but instead uses an upstream weighting approach to limit flow from each dry cell to an adjacent, partially saturated cell. It also uses additional smoothing for conductance and storage functions which allows the Newton numerical method application.

Investigation of the hydrogeology environment and application of a subsurface flow model

The development and application of a simulation model relies heavily on knowledge of the dominant hydrogeology conditions in the study area, the initial piezometry, the boundary conditions, the grid geometry, the trends in the groundwater levels and the assessment of hydraulic properties (Sundell et al. 2019). Once a groundwater model is developed and well calibrated using a technically sound hydrogeological structure for pre-process and post-process, its outcome estimates the groundwater budget on a regional scale and can be used to analyze future water management scenarios (Konikow & Kendy 2005).

Calibration and validation process

The calibration process attempts to bring the model response as close to the observed ones as possible. In the first stage, a trial-and-error method is employed whereas, in the next stage, the pilot points approach was chosen. The pilot point approach is followed due to the lack of field measurements of hydraulic conductivity. Figure 1 indicates that there are only 32 measurements of hydraulic conductivity in the study area which is not enough to correctly express the heterogeneity of such a complex aquifer system. Pilot points are used to express the heterogeneity of hydraulic variables in groundwater systems naturally without the problem of numerical instability in fields with limited measurements data (Doherty & Hunt 2010; Sundell et al. 2019).

As a general approach in groundwater modelling pilot points are located between the observation boreholes and the pumping wells in the direction of groundwater flow (Doherty & Hunt 2010; Baalousha et al. 2019). In this paper, spatial variations of hydraulic conductivity parameter were estimated using pilot points. The pilot points placement can be seen in detail in Figure 3(a) whereas Figure 3(b) is an illustration of groundwater flow. In our simulations 254 pilot points were not randomly chosen but manually set between the observation boreholes and pumping wells. However, it is worth noting that increasing the number of pilot points in the model lengthens the simulation duration. Specific yield value was set as 0.001 corresponding to the average value between fine sand and clay, according to (Johnson 1967). Hydraulic head values were calculated by following the inverse calibration process based on hydraulic conductivity variations.

Figure 3

Conceptual model and the location of pilot points (a) and the groundwater flow (b).

Figure 3

Conceptual model and the location of pilot points (a) and the groundwater flow (b).

The parameter estimation software suite (PEST) (Doherty 2020a, 2020b) is well known and used extensively in calibration of complex environmental systems analysis such as groundwater systems. Simulation tools in PEST have been developed to estimate appropriate values for pilot points for the model calibration using interpolation to the model grid. One such tool is the PLPROC pre-processor that interpolates pilot points to a MODFLOW grid and exports them to MODFLOW arrays. PLPROC runs as part of the ‘model domain’ and alters pilot point values during the calibration process using the Kriging interpolation method (Doherty 2016). In this study, the ordinary Kriging option was followed with no anisotropy, no nugget, and a sill of 1 log cycle. The main advantage of the PLPROC pre-processor choice is that it requires no additional implementation. However, the tool requires input files to be in a predefined format (Doherty 2016). Upon completion of the calibration preprocess, the template files of hydraulic conductivity are created automatically and then used during the PEST run. PEST tool is used in the model optimization process in the regularization mode. The general form of the measurement objective function (Equation (1)) is shown in Figure 2 
formula
(1)
where (hobs indicates the observation hydraulic head, hsim is the simulated hydraulic head, and w the weight that has been applied to the measurement. In this study, we followed the same process as Baalousha et al. (2019), where an equal weight of 1 have been assigned for all the observation points as all have the same importance. As the authors mention regularization mode should be used when the inverse problem is ill posed (Doherty 2020a, 2020b).

In the presented experiments the ill-posed inverse problem was solved with the application of Tikhonov regularization method. The study uses 32 observation wells monitoring the groundwater status during 1987 (Initial Conditions), while the calibration and the validation process was conducted using datasets from 10 observation wells that monitored the groundwater status during 2015 and 2018, respectively.

Post-processing analysis

Post-processing data were analyzed using mod2obs (Doherty 2018), a PEST individual module that compares the simulated hydraulic head with the calibration target. Mod2obs converts the raw model outputs to the model generated observations. Mod2obs works like most programs provided by the PEST suite. The user needs to have prepared the dimension files of the model grid and the coordinates file of the observation points. Then mod2obs outputs a bore sample file of the model generated heads over time, interpolated to the sites and times in the borehole sample file (Doherty 2018).

RESULTS

The model calibration was carried out based on hydraulic head data for the period between 2/11/1987 and 12/11/2015. Data were validated for the period between 12/11/2015 and 11/22/2018. The calibration target was 12/11/2015 and the validation target was 11/22/2018. Scatter plot of observed versus simulated groundwater levels in the study area for the calibration process and the validation process are presented in Figure 4(a) and 4(b), respectively. The parameters which show the goodness-of-fit parameters are the Nash–Sutcliffe Model Efficiency (Ef) (Equation (2)), the RMSE (Equation (3)), and the Coefficient of Determination (R2). The statistical results are presented in Table 3:
formula
(2)
where, is the mean observed value, Ym is the simulated value and Yo is the observed one at time t:
formula
(3)
where, Ym is the simulated value, Yo is the observed one, i is the code number of observation points, and n is the number of observation points.
Table 3

Statistical results for the calibration and the validation period

Statistical parametersCalibration periodValidation period
MAE [m] 1.14 1.11 
RMSE [m] 1.21 1.22 
Ef [1] 0.99 0.99 
Statistical parametersCalibration periodValidation period
MAE [m] 1.14 1.11 
RMSE [m] 1.21 1.22 
Ef [1] 0.99 0.99 
Figure 4

Scatter plot of observed versus simulated groundwater levels in the study area for the calibration (a) and the validation process (b).

Figure 4

Scatter plot of observed versus simulated groundwater levels in the study area for the calibration (a) and the validation process (b).

Figures 57 demonstrate two types of results. The calibration and the validation using the pilot points methodology that show satisfactorily results whereas a high drawdown in hydraulic heads is observed. The above findings strongly suggest and verify the presence of groundwater depletion. The phenomenon can also be observed when examining individually boreholes AD11, SR30 and SR31 which are in the central part of the study area individually including SR63a which is located in the southwest part of the study area (Figure 1). The exact drop figures for each borehole can be seen in the Table 4.

Table 4

Level drop for each borehole

BoreholesLevel drop
AD11 40 m 
SR30 70 m 
SR31 50 m 
SR63a 20 m 
BoreholesLevel drop
AD11 40 m 
SR30 70 m 
SR31 50 m 
SR63a 20 m 
Figure 5

The groundwater trends during the simulation period.

Figure 5

The groundwater trends during the simulation period.

Figure 6

Groundwater trends during the end of wet and the end of the dry period.

Figure 6

Groundwater trends during the end of wet and the end of the dry period.

Figure 7

The overall groundwater balance for the simulation period.

Figure 7

The overall groundwater balance for the simulation period.

Figure 6 shows the results achieved through calibration and validation process where the observed heads converge with the simulated ones before the pumping period and at the end of the pumping period. Large groundwater fluctuations in hydraulic heads are observed between the wet and dry periods. Statistical values for each borehole are shown in the Table 5.

Table 5

Statistical results for each borehole

BoreholesRMSEMAE
PZ7 1.4 m 1.22 m 
PZ67 1.24 m 1.07 m 
SR30 1.24 m 1.07 m 
SR31 1.39 m 1.12 m 
SR63 1.21 m 0.99 m 
BoreholesRMSEMAE
PZ7 1.4 m 1.22 m 
PZ67 1.24 m 1.07 m 
SR30 1.24 m 1.07 m 
SR31 1.39 m 1.12 m 
SR63 1.21 m 0.99 m 

Figure 7 shows the groundwater balance of the area for the period 1987–2018. The groundwater recharge variable is the sum of deep infiltration variable simulated by the UTHBAL model plus the returned irrigation flow was approximately accounted as 10% of the irrigation requirements. The groundwater balance is negative throughout the years, indicating the intensive groundwater abstraction over the region's capacity to replenish the outflows with adequate recharge to minimize the deficit water balance.

DISCUSSION

The simulation of the subsurface system of Lake Karla watershed has been presented in the past by Sidiropoulos et al. (2015, 2019). Sidiropoulos et al. (2015, 2019) used the MODFLOW 2000 code and 15 observation points for hydraulic conductivity parameter applying the geostatistical approach and the model ran from 1987 to 2009. In the current study the model was expanded forward to simulate the groundwater flow from 1987 to 2018 using an updated code, the MODFLOW_NWT mode, which resolves the problem of dry cells, uses 254 pilot points, and updates the groundwater pumping data and observation during the simulation, and utilizes the automated calibration–validation tools PLPROC and mod2obs. The methodology of pilot points deals with the data scarcity, the MODFLOW_NWT does not deactivate the model cells with low water level resulting in the simulation of groundwater flow, without the numerical error that is inserted with the assumption of dry cells by the MODFLOW 2000 code. In addition to the simulation of the groundwater flow and the evolution of all simulation and analysis tools that occur over the years, it is important to mention the difference in the groundwater levels. The simulated water levels of the historical period of 1987–2009 agree well with the simulation results of the previous studies. Furthermore, the groundwater balance indicated an improving trend during the period 2009–2018 but it remains under pressure and intensive and extensive groundwater abstractions. The varying timesteps of the simulation according to the varying dates of groundwater measurements depict clearly the actual and real conditions of the groundwater with increased timely and spatial accuracy.

CONCLUSIONS

This study analyzes the groundwater resources status of the Lake Karla aquifer focusing on a precise simulation–validation procedure of groundwater model based upon new and innovative methodologies. The simulated annual decrease in groundwater levels verify the groundwater depletion and shows that the groundwater system is still under intense pressure even though farming activity in recent years has turned to less water-intensive crops. Moreover, the enduring water balance deficit indicates the necessity for the construction of collective irrigation networks that will eliminate the extravagant number of private own wells, and the need for sustainable water resources management. Unless appropriate measures and policies are implemented, the Lake Karla's aquifer system may face the problem of water scarcity.

This article presents a method that uses a configurable number of pilot points to calibrate and subsequently analyze the groundwater flow system. The use of the Newton smoothing technique in solving the inverse problem with the use of pilot points ensures stability in hydraulic conductivity parameter values and the other hydrogeological variables for Lake Karla's aquifer. Furthermore, the use of tools like PLPROC and mod2obs contribute significantly to the creation of a realistic image of the hydrogeological environment of the Karla aquifer. The pilot points approach is particularly important in data sparse areas of hydrogeological information. Hence, the proposed approach could be used in conjunction with stochastic field generation techniques for assessing model uncertainty and especially in contaminant fate and transport modelling studies.

DATA AVAILABILITY STATEMENT

All relevant data are included in the paper or its Supplementary Information.

REFERENCES

Alamanos
A.
,
Latinopoulos
D.
,
Loukas
A.
&
Mylopoulos
N.
2020
Comparing two hydro-economic approaches for multi-objective agricultural water resources planning
.
Water Resour. Manage.
34
(
14
),
4511
4526
.
http://doi.org/10.1007/s11269-020-02690-6
.
Anderson
M. P.
,
Woessner
W. W.
&
Hunt
R. J.
2015
Applied Groundwater Modeling: Simulation of Flow and Advective Transport
, 2nd edn.
Academic Press
,
London
,
UK
;
San Diego, CA, USA
.
Aquaveo
2019
Groundwater Modeling System Version 10.4 (Version 10.4)
.
Aquaveo
,
Provo
,
UT
. .
Baalousha
H. M.
,
Fahs
M.
,
Ramasomanana
F.
&
Younes
A.
2019
Effect of pilot-points location on model calibration: application to the northern karst aquifer of Qatar
.
Water
11
(
4
),
679
.
https://doi.org/10.3390/w11040679
.
Doherty
J. E.
2016
PLPROC: A Parameter List Processor
.
Watermark
,
Brisbane
,
Australia
, p.
284
.
Doherty
J. E.
2018
Groundwater Data Utilities – Part B, Program Descriptions
.
Watermark Numerical Computing
,
Brisbane
,
Australia
, p.
387
.
Doherty
J. E.
2020a
Model-Independent Parameter Estimation User Manual Part II: PEST Utility Support Software
.
Watermark
,
Brisbane
,
Australia
, p.
268
.
Doherty
J. E.
2020b
PEST: Model-Independent Parameter Estimation and Uncertainty Analysis—User Manual Part I
.
Watermark
,
Brisbane
,
Australia
, p.
394
.
Doherty
J. E.
&
Hunt
R. J.
2010
Approaches to Highly Parameterized Inversion – A Guide to Using PEST for Groundwater-Model Calibration: U.S. Geological Survey Scientific Investigations Report 2010–5169
, p.
59
.
https://doi.org/10.3133/sir20105169
.
Famiglietti
J. S.
2014
The global groundwater crisis
.
Nat. Clim. Change
4
(
11
),
945
948
.
https://doi.org/10.1038/nclimate2425
.
Hayes
P.
,
Nicol
C.
,
La Croix
A. D.
,
Pearce
J.
,
Gonzalez
S.
,
Wang
J.
,
Ahmed
H.
,
He
J.
,
Moser
A.
,
Helm
L.
,
Morris
R.
&
Gornall
D.
2020
Enhancing geological and hydrogeological understanding of the Precipice Sandstone aquifer of the Surat Basin, Great Artesian Basin, Australia, through model inversion of managed aquifer recharge datasets
.
Hydrogeol. J.
28
,
175
192
.
https://doi.org/10.1007/s10040-019-02079-9
.
Johnson
A. I.
1967
Specific Yield – Compilation of Specific Yields for Various Materials
.
U.S. Geological Survey, Water Supply Paper 1662-D
,
Washington, DC
.
https://doi.org/10.3133/wsp1662D
.
Konikow
L. F.
&
Kendy
E.
2005
Groundwater depletion: a global problem
.
Hydrogeol. J.
13
,
317
320
.
https://doi.org/10.1007/s10040-004-0411-8
.
Loukas
A.
,
Mylopoulos
N.
&
Vasiliades
L.
2007
A modeling system for the evaluation of water resources management strategies in Thessaly, Greece
.
Water Resour. Manage.
21
,
1673
1702
.
https://doi.org/10.1007/s11269-006-9120-5
.
Mylopoulos
N.
&
Sidiropoulos
P.
2016
A stochastic optimization framework for the restoration of an over-exploited aquifer
.
Hydrol. Sci. J.
61
(
9
),
1691
1706
.
https://doi.org/10.1080/02626667.2014.993646
.
Niswonger
R. G.
,
Panday
S.
&
Ibaraki
M.
2011
MODFLOW-NWT, A Newton Formulation for MODFLOW-2005: U.S. Geological Survey Techniques and Methods 6–A37
, p.
44
.
OPEKEPE. Integrated Administration and Control System (IACS)
.
Crop Data 2001–2016
.
Ministry of Rural Development and Food
,
Athens
,
Greece
.
Available from http://www.opekepe.gr.
Pisinaras
V.
,
Petalas
C.
,
Tsihrintzis
V. A.
&
Karatzas
G. P.
2013
Integrated modeling as a decision-aiding tool for groundwater management in a Mediterranean agricultural watershed
.
Hydrol. Processes
27
(
14
),
1973
1987
.
https://doi.org/10.1002/hyp.9331
.
Rajabi
M. M.
,
Ataie-Ashtiani
B.
&
Simmons
C. T.
2018
Model-data interaction in groundwater studies: review of methods, applications, and future directions
.
J. Hydrol.
567
,
457
477
.
https://doi.org/10.1016/j.jhydrol.2018.09.053
.
Sidiropoulos
P.
,
Mylopoulos
N.
&
Loukas
A.
2015
Stochastic simulation and management of an over-exploited aquifer using an integrated modeling system
.
Water Resour. Manage.
29
,
929
943
.
https://doi.org10.1007/s11269-014-0852-3
.
Sidiropoulos
P.
,
Tziatzios
G.
,
Vasiliades
L.
,
Mylopoulos
N.
&
Loukas
A.
2019
Groundwater nitrate contamination integrated modeling for climate and water resources scenarios: the case of Lake Karla over-exploited aquifer
.
Water
11
(
6
),
1201
.
https://doi.org/10.3390/w11061201
.
Sidiropoulos
P.
,
Dalezios
N. R.
,
Loukas
A.
,
Mylopoulos
N.
,
Spiliotopoulos
M.
,
Faraslis
I. N.
,
Alpanakis
N.
&
Sakellariou
S.
2021
Quantitative classification of desertification severity for degraded aquifer based on remotely sensed drought assessment
.
Hydrology
8
(
1
),
47
.
https://doi.org/10.3390/hydrology8010047
.
Singh
A.
2014
Groundwater resources management through the applications of simulation modeling: a review
.
Sci. Total Environ.
499
,
414
423
.
https://doi.org/10.1016/j.scitotenv.2014.05.048
.
Sundell
J.
,
Norberg
T.
,
Haaf
E.
&
Rósen
L.
2019
Economic valuation of hydrogeological information when managing groundwater drawdown
.
Hydrogeol. J.
27
,
1111
1130
.
https://doi.org/10.1007/s10040-018-1906-z
.
Upton
K. A.
,
Jackson
C. R.
,
Butler
A. P.
&
Jones
M. A.
2020
An integrated modelling approach for assessing the effect of multiscale complexity on groundwater source yields
.
J. Hydrol.
588
,
125113
.
https://doi.org/10.1016/j.jhydrol.2020.125113
.
Zhou
H.
,
Gómez-Hernández
J.
&
Li
L.
2014
Inverse methods in hydrogeology: evolution and recent trends
.
Adv. Water Resour.
63
,
22
37
.
https://doi.org/10.1016/j.advwatres.2013.10.014
.
This is an Open Access article distributed under the terms of the Creative Commons Attribution Licence (CC BY 4.0), which permits copying, adaptation and redistribution, provided the original work is properly cited (http://creativecommons.org/licenses/by/4.0/).