Groundwater sustainability: a review of the interactions between science and policy

A wide array of numerical models can be used to simulate various elements of the hydrologic cycle, ranging from groundwater models to fully-integrated, physically-based SW-GW models. For example, several groundwater sustainability studies use MODFLOW to simulate groundwater flow (Sakiyan and Yazicigil 2004, Ramesh and Mahesha 2008, Hu et al 2010, Hsu et al 2012, Huang et al 2012, Cao et al 2013, Sarma and Xu 2014, Chung et al 2015, Yihdego and Drury 2016, Sahoo and Jha 2017, Hyndman et al 2017, Al-amin et al 2018, Piscopo et al 2019). If relevant, streamflow is simulated using the MODFLOW river package (Sakiyan and Yazicigil 2004, Hu et al 2010, Huang et al 2012, Cao et al 2013, Chung et al 2015, Sahoo and Jha 2017, Al-amin et al 2018), or accounted for as excess rainfall (Hsu et al 2012). The limitation of this approach is that MODFLOW simulates rivers as either fully connected or fully disconnected, while the transition stages between the two flow regimes exist in nature (Brunner et al 2010). A more physically-based approach includes studies (Nastev et al 2006, Stigter et al 2009) that use FEFLOW (Trefry and Muffels 2007); and the study of Calderhead et al (2012) that uses HydroGeoSphere (Brunner and Simmons 2012), which are capable of simulating streamflow, and unsaturated and saturated flow. Calderhead et al (2012) use HydroGeoSphere to additionally account for land subsidence as a sustainability factor. Using a more physically-based approach could be advantageous. For example, (Brunner et al 2010) compare the relative accuracy of MODFLOW and HydroGeoSphere in simulating SW-GW interaction and show that MODFLOW cannot simulate negative pressures beneath disconnected streams, resulting in an underestimation of the infiltration flux. In addition, the discretization of both the river and aquifer in MODFLOW can cause errors in estimating the position of the water table under the river and in simulating the height of the groundwater mound (Brunner et al 2010).

In numerical models, groundwater recharge can be estimated as a calibration parameter (Nastev et al 2006, Cao et al 2013, Mustafa et al 2018), from a pre-processing step using a numerical steady-state model (Sarma and Xu 2014), as a fraction of precipitation (Stigter et al 2009, Yihdego and Drury 2016), or as the difference between precipitation and actual evapotranspiration (von Brömssen et al 2014, Piscopo et al 2019). Alternatively, a process-based surface hydrology model can be used to simulate recharge, ranging from a lumped representation, such as using an index of effective infiltration (Hsu et al 2012) to a more detailed description of surface water processes. Models in this class include SALUS (Basso and Ritchie 2015) used by (Hyndman et al 2017), HELP (US EPA 1994) used by (Calderhead et al 2012, Sahoo and Jha 2017) and SWAT (Arnold et al 1998) used by (Huang et al 2012). SALUS simulates crop, soil, and water interaction; HELP simulates rainfall, runoff, infiltration and other water pathways, accounting for vegetation, and soil types, among other land-surface aspects; SWAT is a complex watershed model that simulates spatial–temporal impacts of land-use and climate on the hydrologic cycle at the river basin scale. A combination of a cropping system model with SWAT is also used to estimate recharge in agricultural and natural areas (Hu et al 2010). This sequential approach of using a surface water model to estimate recharge for a groundwater model is loose-coupling because it does not allow for an iterative feedback between the surface water and the groundwater model.

To better describe the interaction of surface and subsurface processes, coupled and integrated SW-GW models are emerging. Generally, a fully coupled or integrated SW-GW model would be more desirable. A simple approach is to use a surface water model that has a groundwater module. Zhang et al (2016) use SWAT to evaluate the sustainability of shallow groundwater in an irrigated plain river basin in China. Acero Triana et al (2019) compare MODFLOW and SWAT results, and show that SWAT results will not be always as accurate as using a groundwater model. The study concludes that, although SWAT has a groundwater module, such an approach will not provide a detailed enough description of the groundwater system, leading to wrong conclusions that could misinform policy. Shi et al (2012) present a coupled SWAT and MODFLOW approach with two-way iterative feedback between the two models to evaluate groundwater sustainability in a semi-humid basin in China.

Alternatively, other studies use a fully integrated SW-GW model with detailed description of the surface and subsurface processes, including the unsaturated zone. A few studies in this class implement an integrated hydrological modeling approach using the MIKE-SHE model for sustainable groundwater management in Denmark (Henriksen et al 2008) and China (Shu et al 2012, Qin et al 2013). MIKE SHE (Refsgaard and Storm 1995, DHI 2003) is an integrated model that simulates the main land-surface processes in the hydrologic cycle, and simulates the unsaturated zone using the Richards equation. However, for computational efficiency, (Henriksen et al 2008) disable the unsaturated zone component. Instead of solving the full Richards equation that is computationally expensive, a more practical approach is to simulate the unsaturated zone using a one-dimensional column. This takes advantage of the dominant vertical flow direction within the unsaturated zone when averaged over large area, as in the GSFLOW model (Regan et al 2018). Several versions of GSFLOW are applied to examine the impact of unsustainable groundwater management practices on ecosystem services (Tian et al 2015, Wu et al 2015) and the interlinked impact of climate change and human activities on groundwater sustainability (Feng et al 2018). In a study that is related to groundwater sustainability with respect to land-use and climate change, Ferguson and Maxwell (2012) use ParFlow to compare effects of climate change with pumping and irrigation on terrestrial water and energy budgets of an agricultural watershed in the semi-arid Southern Great Plains, USA. ParFlow (Kollet and Maxwell 2008) is a fully-integrated numerical model that simulates the hydrologic cycle from the bedrock to the top of the plant canopy. Integrated SW-GW models have also been tailored for specific regions, such as California's C2VSim groundwater-SW simulation model (Brush et al 2013). Macewan et al (2017) develop surrogate response functions from C2VSim to evaluate groundwater sustainability in the Central Valley agricultural region in California.

Coupled and integrated hydrological models have several advantages. First, they are more accurate in determining groundwater recharge through infiltration than external water balance models for water resource management in areas with complex SW-GW interaction (Feng et al 2018). This reduction in uncertainty is due to accounting for the dynamics of evaporation and irrigation over space and time during the simulation (Qin et al 2013). In addition, coupled and integrated hydrological models provide additional observation data types other than hydraulic head and concentration data to constrain the model (Wu et al 2014). However, these complex models could involve substantial simulation uncertainty when data for model development and calibration is insufficient, therefore hindering their wider application (Wu et al 2014, Feng et al 2018).

It is not uncommon to use phenomenological models to estimate groundwater sustainability, which is the case, for example, in Kansas (Butler et al 2016), Hawaii (Mink 1981, Liu 2007), California (Miro and Famiglietti 2018), and elsewhere. Phenomenological models include analytical models as well as other types of models as presented at the end of this section. We focus on analytical models because they are frequently used to evaluate groundwater sustainability. This class of models generally utilizes a simple closed-form equation for estimating groundwater sustainability, based mostly on the water budget. Elements of such a budget are simplified by using empirical equations. They are generally termed analytical models (Kalf and Woolley 2005) owing to the way the expressions are derived, in contrast to the spatially distributed numerical models. This class of models is appealing due to its simplicity and is preferred where the use of site-specific detailed modeling approaches is not feasible, mainly due to the lack of appropriate data. In general, at least one criterion, such as the maximum allowed water level decline, is used to constrain water use, satisfying the aquifer water budget. Studies that simply assess sustainability of the aquifer without setting such a criterion were excluded from this review.

Virtually all groundwater sustainability analytical models apply a form of the water budget equation. The models are generally lumped (e.g. Mink 1981, Zhang and Kennedy 2006, Liu 2007, Loaiciga 2008, Aksever et al 2015, Alcala et al 2015, 2018, Bailey et al 2015, Benini et al 2016, Butler et al 2016, 2018, Bailey and Tavakoli Kivi 2017, Miro and Famiglietti 2018), which is consistent with the nature of this class of models. As described by Butler et al (2016, 2018), a typical lumped model (also known as a single-cell aquifer model) considers groundwater recharge as an input and pumping and natural discharge as outputs. Water height is a parameter that is related to pumping through a linear relationship and can be translated to water volume given the surface area and storage coefficient of the aquifer. Baseline groundwater sustainability generally corresponds to the average level of pumping that causes zero average groundwater-level change, or any other sustainability criterion.

Analytical models can differ regarding processes included, approaches to quantify budget elements, and the criteria set for estimating the groundwater sustainability. Some models include more processes than others in their water budget equation and model formulation. For example, Zhang and Kennedy (2006) account for artificial sources of groundwater recharge and alteration of groundwater systems due to development in an urban context. Similarly, the model of Voudouris (2006) accounts for exploitable dynamic groundwater reserve, artificial recharge, and irrigation return. Loaiciga (2008) uses a graphical method to account for surface reservoir sizing. Miro and Famiglietti (2018) use a modified water-budget approach for confined aquifers. Benini et al (2016) utilize an index that accounts for salinization to assess aquifer vulnerability to climate and land-use changes. Other models use various techniques to account for fresh water lens in coastal aquifers (Liu 2007, Bailey et al 2015, Bailey and Tavakoli Kivi 2017). In addition, to improve the representation of the groundwater system, a few studies extend analytical models beyond the lumped or single cell aquifer. For example, Jang et al (2012) treats the shallow and deep aquifers as two compartments and perform the water balance sequentially. Also, Jafari et al (2018) develop a distributed model that utilizes the area of a Thiessen's polygon surrounding each piezometer in the aquifer and the respective value of water level decline.

Hawaii is a good example of using an analytical model to evaluate groundwater sustainability. Mink (1981) developed the Robust Analytical Model (RAM) based on a steady state assumption for an aquifer where inflows (recharge) equal outflows (leakage plus pumping), with all pumping lumped as a single value. The RAM models have been used as the primary tool to implement sustainable yield water policy in Hawaii (figure 6). Estimating sustainable yield is required to update the Water Resources Protection Plan, which is a component of the Hawaii Water Plan (CWRM 2008, 2019b). The model provides a parabolic relationship between average hydraulic head and pumping. The curve can be used to estimate optimal pumping based on known recharge and the ratio between equilibrium and initial heads. Difficulties in applying RAM are mainly related to identifying a value for the equilibrium head. Values suggested by water management agencies range between 50%–75% of the initial aquifer head. RAM2 (Liu 2007, Liu and Dai 2012) improved the original RAM model by accounting for salinity. The improvements also include analytically estimating the equilibrium head to prevent salinity of pumped water from exceeding a certain acceptable level. Limitations of RAM and RAM2 models include the inability to account for groundwater leakage. Thus, the models cannot address societal preferences with respect to ecosystem services by not addressing leakage components, such as submarine groundwater discharge, spring discharge, base-flow, evapotranspiration, and drains. Simulating leakage is also important in accounting for change in lateral inflow between different aquifer administrative units.

Although analytical models, which are the most commonly used class of phenomenological models, are widely used to evaluate groundwater sustainability, they are subject to criticism (e.g. Kalf and Woolley 2005, Henriksen et al 2008, Mulligan et al 2014). A first limitation as stated by Kalf and Woolley (2005), is that the models are unsuitable for basin groundwater sustainability estimation, mainly because they lack the ability to simulate leakage, and their treatment of inflow and outflow component interaction cannot be rigorous. Not being able to simulate groundwater leakage is critical. As explained by Pereau and Pryet (2018), if leakage is not considered, then over a long period the single-cell aquifer will dry out if pumping rates exceed recharge, and will overflow if pumping remain below recharge (Pereau and Pryet 2018). Seward et al (2006) also emphasize model limitations if only natural recharge and pumping are considered. The implications of this limitation on groundwater dependent ecosystems is discussed in section 5.2. More recently, Pereau and Pryet (2018) addressed this limitation by developing an analytical model in which leakage explicitly appears in the water budget equation as a function of the groundwater level using a linear conductance model.

Another limitation of analytical models is the lack of consideration for spatiotemporal relationships, which has several implications. First, analytical models assume instantaneous change in the water-level of the model domain, while in reality pressure diffuses gradually over space from pumping wells, and the impacts of pumping can be delayed in time (Pereau and Pryet 2018). Second, not accounting for the aquifer spatial distribution of hydrogeological data can yield inaccurate results. As indicated by Thomas and Famiglietti (2015), groundwater sustainability is critically dependent on accurate estimation of the temporal and spatial variability of groundwater behaviors as a response to both natural and anthropogenic influences. For example, Oki and Meyer (2001) compare results obtained by both the numerical model SHARP (Essaid 1990) and the analytical model RAM (Mink 1981) for a major aquifer in Hawaii. The results indicate that the field measured declines of water levels are larger than those predicted by RAM, which are consistent with the results of the numerical model analysis. The model RAM underestimates water-level declines in areas where a low-permeability confining unit exists, and in the vicinity of pumping wells. Also, Elshall et al (2018) show that the RAM2 model (Liu 2007) overestimates sustainable yield by up to 30%, seriously underestimating the chance of aquifer salinization. Third, distributed spatiotemporal response of soil moisture and SW-GW interaction to variations of precipitation and pumping can be vital when incorporating ecosystem services (Henriksen et al 2008) and human activities modeling (Mulligan et al 2014) as detailed in section 5.2, and section 5.3, respectively. Fourth, considering the spatial location of pumping wells is important. For example, in the presence of saltwater intrusion, turning off certain pumping wells near the coast can increase sustainable yield (Elshall et al 2018). In addition, a pumping well located close to a discharge zone could be critical for groundwater dependent ecosystems (Seward et al 2006).

Despite these limitations, analytical models have several advantages and usages. First, they can be applicable to specific cases, such as when groundwater sustainability evaluation is based on water-level response alone. This includes cases where groundwater levels are decoupled from environmental and river flow effects (Macewan et al 2017), and for fossil groundwater resources where natural recharge and discharge are negligible (Pereau and Pryet 2018). Second, analytical models are useful as first-stage assessment tools (Kalf and Woolley 2005) and for conceptual and generic discussions on groundwater management strategies (Pereau and Pryet 2018). Third, analytical models are especially useful when there are insufficient site-specific data to develop a high-fidelity numerical model with more mathematical and geological realism. In many locations, especially in developing countries, site-specific data needed for numerical models are generally limited, which reduces the usefulness of numerical models. To overcome some of the limitations of analytical models and the high computational cost of numerical models, some studies use more elaborate phenomenological models to evaluate groundwater sustainability. This includes the economic-engineering optimization model of California CALVIN (Harou and Lund 2008), and machine learning models (Salem et al 2017). Due to the computational efficiency of the hydrological module, these models encourage accounting for multiple aquifer governance factors as discussed in section 5.3.

In summary, the choice of using a phenomenological versus numerical model is case-specific and depends on data availability, aquifer type, and sustainability factors of interest. As discussed by Hill (2006), simpler models are often preferred, given that they are characterized by fewer parameters and shorter execution times. A recommended approach is to start with simple models and slowly build model fidelity and model complexity to reach the best fit with measurements, while avoiding underfitting or overfitting the observation data (Hill 2006, Elshall et al 2019). The challenge, especially with numerical models, is to develop the most parsimonious model given groundwater sustainability factors of interest (figure 3) and plausible data. To avoid what Voss (2011b) calls 'groundwater modeling fantasies,' in which the modeler finds herself 'adrift in the details,' the focus should be on 'down to earth' model development with reasonable data requirements (Voss 2011b). Data collection and model development can be time consuming and expensive, and require specialized training. The shared responsibility of stakeholders and the scientific community is to work together to define the state of knowledge and hydrological modeling tools required to address societal preferences and answer groundwater sustainability questions.

We use coastal aquifers as an example to illustrate the role of hydrological modeling in the scientific evolution of groundwater sustainability concepts. Several studies review management aspects of coastal aquifers. These include dealing with groundwater scientific challenges to meet societal needs in coastal areas (Michael et al 2017), modeling tools and challenges to characterize and manage coastal aquifers (Werner et al 2013), computational and conceptual issues in the calibration of seawater intrusion models (Carrera et al 2010), the use of simulation optimization method for coastal aquifer management (Singh 2014, Ketabchi and Ataie-Ashtiani 2015, Singh 2015, Sreekanth and Datta 2015), and the impact of sea-level rise and climate change on coastal aquifers (Holding et al 2016, Ferguson and Gleeson 2012, Ketabchi et al 2016, Werner et al 2017, 2013). Groundwater sustainability is discussed by a number of these articles (Ketabchi and Ataie-Ashtiani 2015, Werner et al 2017, 2013). Werner et al (2017) review approximate approaches for estimating groundwater sustainability in atoll islands, which include determining a percentage of recharge or rainfall as the groundwater sustainability, or using analytical solutions with the Dupuit approximation and the Ghijben-Herzberg principle, along with field measurement and monitoring techniques. For example, Benini et al 2016 use the Dupuit-Ghijben-Herzberg relationship (Fetter 2001) and a salinization vulnerability index based on water table level to evaluate the longterm impact of land-use, climate changes, and sea-level rise on groundwater sustainability. These approaches are very useful for first assessment, and when data is limited, as discussed in section 5.1.3. In addition, field measurements techniques, such as hydrogeophysical approaches (e.g. Vouillamoz et al 2012), and hydrogeochemical approaches (e.g. Bryan et al 2016), are used to quantify groundwater resources in coastal areas. In general, numerical modeling approaches to evaluate groundwater sustainability can better project outcomes under a wide array of natural and anthropogenic scenarios. As documented by Post et al (2018), a numerical modeling approach better enables the analysis of freshwater volumes and fluxes of submarine groundwater discharge in comparison to other techniques. Hence, this section mainly focuses on numerical modeling studies, considering that phenomenological modeling studies are discussed in section 5.1.2. We conducted full-text assessments for 45 research articles about groundwater sustainability evaluation in coastal aquifers, and 32 articles were eligible for inclusion. We exclude most of the studies that only use synthetic case studies, as well as articles that do not provide sufficient justifications for not accounting for density-dependent flow. However, we include studies that use a constant density model, while discussing the limitation of this approach. These justifications include fewer data requirements, easier model development, and less computational cost.

Groundwater sustainability evaluation is challenging in coastal aquifers as they are more vulnerable to additional threats than inland aquifers. The combination of natural and anthropogenic threats significantly affects water supplies from nearshore fresh groundwater lenses. Threats can be summarized from a number of coastal aquifer studies, including Ferguson and Gleeson (2012), Michael et al (2017), Rotzoll and Fletcher (2013), and Werner et al (2017). They include (1) drought due to the often large inter-annual variability of natural rainfall cycles, (2) seawater inundation due to storm surges, extremely high tides, and tsunamis, (3) groundwater inundation in coastal-plains due to sea level, (4) inundation by flooding of low-lying areas by seawater, (5) thinning of the fresh groundwater lenses and widening of the mixing zone due to pumping-induced drawdown and decrease of spring and submarine groundwater discharge, (6) upconing at individual wells due to over-pumping, (7) local salt storage in salt-making facilities, and (8) nutrient and antibiotics use in aquaculture. These factors adversely affect the thickness of the fresh groundwater lenses, the width and shifting of the mixing zones, the resilience of pumping wells to upconing, the resilience of fresh groundwater lenses to seawater inundation and drought, and quantity and quality of groundwater discharge to groundwater dependent ecosystems. Accordingly, groundwater sustainability evaluation in coastal aquifers should ideally take into account density dependent flow, mixing processes, recharge rates, drought periodicity, and alternative sources of freshwater (Werner et al 2017). Additional issues include accounting for submarine groundwater discharges and related water quality, including nutrient levels that are important for terrestrial ecosystem services (e.g. Hugman et al 2015, Robertson et al 2019, Burnett et al 2020) and specific coastal groundwater dependent ecosystems such as coastal lagoons (e.g. Menció et al 2017), coral reefs (e.g. Delevaux et al 2018), and nearshore marine environments (e.g. Duarte et al 2010, Hugman et al 2015).

In coastal aquifers, the volume of freshwater lens depends on the location of the freshwater-saltwater interface and the recharge-discharge mechanisms of the aquifer. Budget analysis alone cannot be sufficient for groundwater sustainability evaluation. More accurately, groundwater sustainability can be evaluated by observing the changes in the volume of freshwater lens under different pumping schemes and considering recharge and sea-level rise scenarios. In addition, the analysis would be enhanced by utilizing different aquifer conceptual models and including different aquifer performance and governance factors. The most commonly used aquifer performance factors when evaluating groundwater sustainability are drawdown values (e.g. Liu et al 2006, Hugman et al 2015, Renau-Prunonosa et al 2016, Zhou et al 2017), salinity (e.g. Post et al 2018), or both (e.g. Charalambous and Garratt 2009, Praveena and Aris 2010, Praveena et al 2012, El-Kadi et al 2014, Lathashri and Mahesha 2016, Pholkern et al 2019). Salinity thresholds can be set to regulatory values or as determined by local conditions. For example, Zhao et al (2016) study considers a monitoring point to be affected by seawater intrusion when the detected chloride concentration is higher than the limit for freshwater use of 250 mg l⁻¹ as set by U.S. Environmental Protection Agency (1973). Alternatively, Pholkern et al (2019) set groundwater salinity threshold to 1000 mg l⁻¹ TDS. When salinity is not explicitly considered, changes in head-level can be taken as a surrogate for salinity and for the location of the saltwater–freshwater interface. For example, Sedki and Ouazar (2011) define an altitude value from mean sea level for hydraulic head along the coastal boundary as a threshold to prevent saltwater intrusion. Similarly, Liu et al (2006) use an adaptive management approach with a groundwater management index based on groundwater levels in wet and dry seasons. In addition to drawdown and salinity, other aquifer performance factors of concern include spring discharge (El-Kadi et al 2014, Hugman et al 2015, 2017, Burnett et al 2020), and submarine groundwater discharge (El-Kadi et al 2014, Hugman et al 2015, Post et al 2018). Groundwater discharge volumes and salinity are particularly important for groundwater dependent ecosystems and various human activities, such as in agriculture zones and culturally significant sites (Burnett et al 2020). These additional aquifer performance factors add more restrictive constraints on fresh groundwater utilization.

Fresh groundwater lenses are highly vulnerable to salinization due to excessive groundwater pumping rates, coupled with natural recharge variations and sea-level rise. Determining the groundwater sustainability from freshwater lenses is challenging because the lens response during drought periods and the long-term effects of pumping are both difficult to predict; Post et al (2018) show that the lens contraction caused by pumping is nearly a linear function of the total pumping, yet this might not be the case under different recharge scenarios. This is mainly due to aquifer heterogeneity, including the existence of preferential flow paths (Elshall et al 2013, Werner et al 2013). In addition, using multiple geological conceptualization of the subsurface can have a large impact on groundwater sustainability estimates similar to the impact of climate change (Pholkern et al 2019). A major gap in groundwater sustainability literature in coastal aquifer is the absence of adequate characterization of subsurface heterogeneities, especially under different conceptual geological models.

Significant uncertainties also exist with respect to recharge evaluation, especially regarding climate modeling and difficulties in the downscaling processes necessary for local aquifer scale simulations. Typical studies evaluate recharge scenarios utilizing data for climate change, land-use, artificial recharge, among others. Different top-down and bottom-up approaches have been adopted to evaluate natural recharge due to climate change. A top-down approach evaluates a specific future climate change scenario, while bottom-up approach assesses the impact of systematic increase or decrease in recharge. Top-down approach includes recharge estimation under different Intergovernmental Panel on Climate Change (IPCC) scenarios for future climate projections, which includes using the Special Report on Emissions Scenarios (SRES) scenarios (Rasmussen et al 2013, El-Kadi et al 2014, Benini et al 2016, Hugman et al 2017), and the Representative Concentration Pathways (RCP) scenarios (Pholkern et al 2019). The SRES scenarios are designed around a set of development assumptions of regional and global patterns of economic growth and environmental sustainability. They were superseded in IPCC Fifth Assessment Report released in 2014 by RCP scenarios that are based on the change in radiative forcing and other forcing agents. Example studies include Rasmussen et al (2013), who assume an increase in groundwater recharge of 15% in the period from 2010 to 2100 for a coastal aquifer in Denmark, based on outputs from regional climate models representing an IPCC A2 scenario of regionally oriented economic development, and a B2 scenario of local environmental sustainability. Alternatively, Pholkern et al (2019) evaluate three RCP scenarios of low-emissions mitigation scenario leading to a very low forcing level (RCP2.6), a midrange mitigation emission scenario (RCP4.5), and a very high baseline emission scenario (RCP8.5). The annual precipitation of the study basin at North East Thailand indicates that average annual rainfall tends to increase by 6%, 7%, and 11%, under the RCP2.6, RCP4.5, and RCP8.5 scenarios, respectively. The study shows that groundwater sustainability is consistently higher under RCP8.5, RCP4.5 and RCP2.6, respectively, under each of the four evaluated conceptual models. The study found that uncertainty about the model geological structure and recharge due to future climate change are significant in comparison to uncertainty about boundary conditions. Other than IPCC scenarios, El-Kadi et al (2014) evaluate a drought scenario that is developed from historical drought events. Similarly, Zhao et al (2016) conduct a frequency analysis using historic regional rainfall data to develop three recharge scenarios representing wet, normal, or dry weather conditions. On the other hand, bottom-up recharge perturbs current or pre-development recharge by a percentage decrease (Liu et al 2008, Rasmussen et al 2013, Unsal et al 2014, Lathashri and Mahesha 2016), or increase (e.g. Praveena and Aris 2010). For example, Renau-Prunonosa et al (2016) use a bottom up approach to identify critical saltwater intrusion areas under wet, intermediate, and dry conditions.

Change in land-use is another important factor that controls recharge. To study the impact of human activities on future recharge, land-use scenarios are developed to evaluate increased urbanization (El-Kadi et al 2014), stakeholder driven socioeconomic development (Benini et al 2016), and irrigation recharge return due to different agricultural practices (Charalambous and Garratt 2009). Land-use change includes developing adaptation scenarios. For example, Benini et al (2016) show that climate change is the most important driving force in comparison to land-use change dynamics. However, watershed restoration can help in decreasing water deficit.

A number of studies emphasize recharge related factors in assessing changes in fresh groundwater lenses for groundwater sustainability evaluation. These include developing river boundary condition scenarios (Hugman et al 2015, Pholkern et al 2019), and inland vertical boundary condition scenarios (Pholkern et al 2019). Pholkern et al (2019) show that the impact of top-recharge due to climate change on groundwater sustainability evaluation is more critical than the impact of uncertainty regarding inland vertical boundary recharge. Another important factor is the increase in recharge through aquifer management (Praveena et al 2012, Hugman et al 2017). Several studies examine increasing groundwater sustainability by aquifer recharge using treated wastewater or water transfer (Vandenbohede et al 2009, Praveena et al 2012, Al-Maktoumi et al 2016, Vijay and Mohapatra 2016, Zhou et al 2017). Another management strategy is to passively enhance infiltration from streams by using small dams or weirs and increasing the infiltration capacity of the streambed (Hugman et al 2017). Hugman et al (2017) show that using existing large diameter wells to infiltrate surplus water from the large surface water reservoirs that are used for public water supply is effective for reducing saltwater intrusion and increasing groundwater sustainability. In general, all studies related to climate, land-use, and managed aquifer recharge show that the thickness of the fresh groundwater lenses can generally be directly proportional to recharge. However, some areas in the aquifer can be more vulnerable to decrease in recharge than others. Also, the increase in groundwater sustainability with respect to the increases in recharge volume can be a nonlinear process under different recharge scenarios.

While increase in recharge will generally increase sustainable groundwater the impact of sea-level rise is not uniquely defined. This is mainly due to the intricate relation between several factors that include land-surface inundation, flux-control or head-control inland boundary conditions, recharge variations, aquifer bed slope, and aquifer thickness to aquifer length ratio. For example, the study by Rasmussen et al (2013) shows that only changes in groundwater recharge had an effect on the saltwater intrusion for a system with flux controlled boundary conditions, but for a system with head-controlled boundary conditions, changes in recharge, sea level, and the boundary value (e.g. the stage a drainage canal) are important for saltwater intrusion. The impact of sea-level rise on land-surface inundation and aquifer salinization have been studied before to a large extent (Ferguson and Gleeson 2012, Rotzoll and Fletcher 2013, Ketabchi et al 2016). The Ketabchi et al (2016) review article discusses and highlights the gaps in literature about factors that control the impact of sea-level rise on saltwater intrusion. Here we focus on the few studies that assesses the impact of sea-level rise on groundwater sustainability evaluation. Among these, the study of Unsal et al (2014) assesses the impact of a worst case, maximum global average sea-level rise estimated by IPCC corresponding to an increase of 0.88 m at the end of 100 years, showing a reduction in groundwater sustainability by about 25% from the base-line case. Similar studies show that sea level rise will increase the risk of seawater intrusion and a related decrease in groundwater sustainability (Zhou et al 2017). On the other hand, Lathashri and Mahesha (2016) illustrate that saltwater intrusion is intensified in the area adjoining the tidal rivers rather than that due to the sea alone, and concludes that regional sea level rise of 1 mm/year has minimal impact on groundwater sustainability. Similarly, Rasmussen et al (2013) show that the recharge variation has a much more dominant impact on saltwater intrusion in comparison to a sea level rise of 0.75 m by 2100, under flux-controlled inland boundary conditions that did not result in significant impacts on saltwater intrusion.

Tidal effect is another factor that can influence groundwater sustainability in coastal aquifers. Tidal dynamics induce periodic fluctuations of the water table in coastal aquifers resulting in further saltwater intrusion. For example, Zhao et al (2016) show that saltwater intrusion risk increases with high coastal tides. This study simulates tidal fluctuations by setting the specified head on the coastline boundary to the sea level, and evaluates tidal impact on the groundwater sustainability given nine future recharge scenarios due to climate change. The study shows that without considering tidal effects, the extent of seawater intrusion would be underestimated because tidal effects can influence the model boundary conditions, and accordingly, the groundwater levels in the aquifers. Lu et al (2013) evaluate high, medium, and low tide levels with different groundwater pumping schemes. The study shows that tide-induced seawater intrusion can significantly affect the groundwater levels and salinities when groundwater pumping exceeds a certain bound.

In the literature, diverse approaches are adopted to evaluate the groundwater sustainability in coastal aquifers. One approach is to use the Hill method (Freeze and Cherry 1979) that plots simulated pumping rates against drawdowns. For example, Liu et al (2006) use this method to obtain a linear relationship between annual pumping and drawdown, and to define a groundwater management index based on the status of groundwater level changes to allow local government officials to dynamically adjust the pumping scheme. A second approach is the zero water-level change method (e.g. Pholkern et al 2019) that defines safe yield as the average pumping over a time period in which the groundwater storage level is maintained at the beginning and end of this period. The main disadvantage of these two methods is that they do not explicitly account for groundwater discharge. However, their benefit is that they have a defined timeframe. A third approach is to define a long timeframe for the evaluation, which could be challenging, and must be based on timescales of several decades (Post et al 2018). One option is to define a multi-generational timeframe, for example, until the end of the current century (Rasmussen et al 2013, Benini et al 2016, Hugman et al 2017). Rasmussen et al (2013) further aim to capture the long-term effects of the imposed climate change by treating recharge and sea level as unchanged for an additional 200 years as post-processing steps. Another methodology is to evaluate groundwater sustainability when the numerical simulation reaches a new equilibrium condition, which can be a long period. For example, Post et al (2018) show that, even for a small permeable island like Bonriki, it takes three decades to realize new equilibrium conditions that reflect pumping stresses.

Another commonly used approach for groundwater sustainability evaluation is to assess different pumping scenarios against a set of aquifer performance and governance factors (e.g. Liu et al 2008, Rasmussen et al 2013, El-Kadi et al 2014, Hugman et al 2017). For example, Hugman et al (2017) evaluate sustainability by assessing two pumping scenarios in which groundwater use for irrigation is decreased by 25% and 50% due to decreasing irrigated land, increasing agricultural efficiency, or the promotion of the use of alternative water sources. In addition, the study evaluates the impact of spatiotemporal distribution of public-supply pumping given different users in wet and dry seasons. The aquifer factors considered comprise spring discharge and salinity levels for groundwater dependent ecosystems (Hugman et al 2017, Burnett et al 2020). Pumping schemes could also include alternative groundwater pumping options. For example, to reduce salinization, Alam and Olsthoorn (2014) evaluate the use of scavenging wells to pump shallow brackish water to replace deep fresh groundwater for irrigation.

Using a simulation optimization method to evaluate the optimal and sustainable pumping of groundwater from coastal aquifers is another common approach. This can include a manual search. For example, Pholkern et al (2019) divide the pumping scheme into eight clusters, and decrease pumping in each cluster until a feasible solution that does not violate groundwater sustainability factors is obtained. The groundwater sustainability factors considered in this study comprise drawdown change that is set to zero, reflecting the achievement of a balance between water inflow and outflow, and groundwater salinity levels at existing pumping wells that do not exceed the salinity threshold. Alternatively, a more common approach is to use an automatic search that determines the optimal pumping that satisfies the aquifer performance and governance factors (e.g. Qahman et al 2005, Schoups et al 2006, Rejani et al 2009, Javadi et al 2012, Nocchi and Salleolini 2013, Uddameri et al 2014, Renau-Prunonosa et al 2016, Kamali and Niksokhan 2017, Burnett et al 2020). For example, the study by Qahman et al (2005) develops a simulation optimization method to maximize pumping and the profit of selling water, minimize salinization and the operational and water treatment costs, and satisfy drawdown thresholds. The decision variables are pumping rates, while the constraints are head-levels and concentrations. To reduce the high computational cost of density dependent groundwater flow models, the study by Renau-Prunonosa et al (2016) uses a similar approach, but employs a constant density groundwater flow model, which assumes a direct relationship between drawdown and seawater intrusion.

The finite difference MODFLOW model is commonly used to simulate groundwater flow for groundwater sustainability evaluation in coastal and other aquifers (e.g. Liu et al 2006, Sedki and Ouazar 2011, Renau-Prunonosa et al 2016, Zhou et al 2017). However, density dependent models, such as SEAWAT (Langevin 2009), are more appropriate for coastal aquifers. Depending on the nature of the site, SW-GW interaction can be considered in SEAWAT by using the MODFLOW river package (Lathashri and Mahesha 2016), or applying the river boundary condition (Pholkern et al 2019). Alternatively, other studies (Charalambous and Garratt 2009, Hugman et al 2015, 2017) use the finite element model FEFLOW to simulate SW-GW interaction with density dependent flow. Other finite element codes used to evaluate groundwater sustainability include CODSEA-3D (Gambolati et al 1999, Lecca 2000), FEMWATER (Lin et al 1997), MOCDENS3D (Oude Essink 1998), and SUTRA, that are respectively used by Qahman et al (2005), Lal and Datta (2019), Vandenbohede et al (2009), and Burnett et al (2020). Recharge for the groundwater sustainability evaluation for coastal aquifers is generally estimated using an approximate method or a loose-coupling approach. Approximate methods include utilizing a percentage of rainfall as recharge (Sedki and Ouazar 2011, Zhou et al 2017), using water budget method (Liu et al 2006, Charalambous and Garratt 2009, Hugman et al 2015, Post et al 2018), or setting recharge as a calibration parameter (e.g. Charalambous and Garratt 2009). Alternatively, recharge can be estimated using a surface water model, such as the physically-based surface water model used by El-Kadi et al 2014) or the HELP3 model used by Pholkern et al 2019) with loose coupling to the groundwater model. Examples for public domain, constant density models include GSFLOW (Regan et al 2018) (fully coupled) and MODFLOW-SWAT (Kim et al 2008) (integrated). For density dependent flow, and to the best of our knowledge, only the coupled FEFLOW and MIKE-SHE model (Yamagata et al 2012) is available, which are not public-domain software. There is certainly a gap in coastal aquifer studies towards developing similar public domain models.

We selected karst aquifer as a second example for three main reasons. First, karst aquifers are an important water resource, because ∼15% of the global ice-free continental surface is characterized by the presence of karstifiable carbonate rock, and ∼16% of the global population lives on karst (Goldscheider et al 2020). Secondly, karst terrain has a high vulnerability in terms of groundwater quality and subsidence. Thus, the development of methods to address data uncertainty and unavailability, and the development of mathematical models that fully capture the karst complexity, and accordingly the estimation of groundwater sustainability is not a trivial task. The third reason is that karst aquifers have generally received less attention in groundwater hydrology research than porous media and fractured rock aquifers.

The model simulation cannot be more accurate than the input data. Modeling flow and transport in subsurface media has a distinctive position among all other fluid mechanics applications, since hydrogeological data is uncertain especially when it comes to karst aquifers. While porous media aquifers are often heterogeneous, conductive fractures are generally expected to have highly erratic heterogeneity, directional dependence, a dual or multicomponent nature, and multiscale behavior (Neuman 2005). Yet on top of the matrix and fracture porosities, karst media is commonly characterized by a tertiary porosity of interconnected conduits. In comparison to porous and fractured media, data unavailability is a peculiar feature that karst mediums exhibit. As these conduits occupy small fraction of the aquifer area and volume, the probability of intercepting one of the high permeability channels during a random drill is very low (Worthington and Ford 2009). In other words, the existing geophysical methods, or other site investigation methods, will not guarantee the location of active conduits, but will only increase the probability of finding them. Thus, while these conduits are controlling most of the flow, they can remain undiscovered. Even a common definition of karst is unclear. As noted by Bakalowicz (2005) 'because of the complexity of karst medium, many think that every karst aquifer is currently considered as being representative of only itself'. This data uncertainty and unavailability problem makes karst the most challenging topic in hydrogeology.

Apart from the challenging data issue, the development and application of mathematical models to simulate flow and transport in karstified rocks are still in their infancy. Karst aquifers did not receive much attention in the past in comparison to porous media and fractured aquifers. In part, many of the methods such as the stochastic flow and transport equation that were developed for porous and fractured media did not extend to karst media. On the other hand, the methods that were developed or adopted to address the special complexity of karst are still under development in a sense that none of the existing codes can fully capture the full complexity of karst robustly. The complexity arises due to the presence of (1) diffuse, fracture, and conduit recharge, (2) laminar flow in porous matrix, fracture, and conduits filled with sediments, (3) turbulent flow in vuggy matrix, fracture, and conduits, (4) saturated or unsaturated flow for both matrix and fractures, (5) pipe or open channel flow for conduits, (6) solute transport under laminar and turbulent flow conditions, and (7) diffuse and spring discharge.

Unlike porous media aquifers, phenomenological models are generally more commonly used in karst aquifers than numerical models due to the aforementioned difficulties. Phenomenological models in karst, which are termed global models, can be classified as black box or grey box models. These models are represented functions that relate the system inputs, outputs, and responses. These do not provide specific spatial information, but temporal information of the global system responses to different inputs and thus they can be used for predication. Also, by analyzing the system inputs and outputs, the system components or functions can be resolved. Thus, these methods can be extended in application to characterize vulnerability at a spring catchment (Butscher and Huggenberger 2008), nitrate transport (Pinault et al 2001), coastal aquifer (Pinault et al 2004), and surface water and groundwater exchange (Pinault and Schomburgk 2006), among other applications related to groundwater sustainability. The use of these methods is common in karst hydrology, since they can bypass the data problems as previously highlighted, yet the range of their applications is limited. On the other hand, the applications of numerical models can be wider with respect to assessment and management of water resources and mitigation of hazards such as subsidence.

Studies of groundwater sustainability for karst aquifers face a number of challenges in numerical modeling, aquifer characterization, and groundwater and surface water interaction. The major challenge is that the dual-media (matrix and conduit) characteristics of karst aquifers make groundwater sustainability definition complicated because groundwater flow is laminar in a matrix, but turbulent in conduits. Therefore, groundwater sustainability should be defined for matrix and conduits separately. However, since the spatial distribution of conduits is always unknown, groundwater sustainability cannot be defined specifically for conduits. One solution to this dilemma is to use the concept of equivalent porous medium by assuming that conduit flow is also laminar (or more specifically that the representative elementary volume of a karst aquifer is large enough to average out local influence of karst conduits), and thus can be studied by the tools developed for porous and fractured media (Scanlon et al 2003, Ghasemizadeh et al 2012). The field-scale modeling of Davis et al (2010) demonstrates that this concept is valid when using MODFLOW and MT3DMS to simulate groundwater flow and nitrogen transport in the Woodville Karst Plain of northern Florida, USA. However, Kuniansky (2016) points out that the concept of equivalent porous medium cannot simulate the response of the groundwater system to short-term hydrological events such as heavy rainfalls. This is not surprising because the concept of equivalent porous medium is merely an approximation and cannot fully describe turbulent flow in conduits.

In the last decade, building on the concept of discrete conduit continuum, MODFLOW-CFP (conduit flow process) (Shoemaker et al 2007) and CFPv2 (Reimann et al 2014) have been developed to couple laminar flow in porous media and turbulent flow in conduits. While the concept of discrete conduit continuum is theoretically more advanced than the concept of equivalent porous medium and the computer codes have been used for field-scale studies (e.g. Gallegos et al 2013, Xu et al 2015), the research of coupling laminar and turbulent flows in karst aquifers is still in its infant stage, and many mechanisms of the coupling are still largely unknown. For example, while MODFLOW-CFP and CFPv2 assume that the flow exchange between matrix and conduits is a linear function of the difference between the heads of matrix and conduits, this assumption may not be valid (Pacheco Castro et al 2020). In addition, field applications of the computer codes require knowing spatial distribution of karst features (e.g. conduits and sinkholes), where the distribution is the only large unknown at any field sites. Lacking quantitative tools of groundwater modeling is the major challenge of evaluating groundwater sustainability in karst aquifers.

Efforts have been spent to develop quantitative tools without requiring detailed aquifer characterization. For example, models based on the concept of lumped-parameter equivalent porous medium have been developed (e.g. Scanlon et al 2003, Santos and Andreu 2010, Jódar et al 2014). Such models treat the matrix as a reservoir and conduits as another reservoir, and thus do not require knowing the spatial distribution of karst features. However, many parameters (e.g. reservoir areas and water depth) of these models are fudge factors, making the model not useful in practice. These kinds of problems also exist in statistical models developed for groundwater sustainability studies. For example, Yin et al (2012) develop an artificial neural network for an assessment of sustainable yield of karst water in Huaibei, China. Using machine learning techniques for groundwater sustainability studies may be a future direction if available data are sufficient to develop a statistical model (e.g. an artificial neural network), but insufficient to develop a mechanistic model (e.g. a MODFLOW-CFP model).

Studying groundwater sustainability of karst aquifers requires advanced understanding of groundwater and surface water interactions in the aquifers, because groundwater and surface water are strongly and dynamically connected through karst features such as sinkholes. For example, a stream can entirely disappear when it falls into a sinkhole and becomes groundwater. This phenomenon is well described in the review article of Tihansky (1999). The groundwater and surface water interactions make the groundwater sustainability studies challenging because of the dynamic nature of the interactions. For example, a productive karst aquifer may become unsustainable quickly if its recharge from surface water is reduced substantially over a drought period. Studying groundwater and surface water interactions again faces the difficulty of lacking knowledge and data. For example, it is rarely known how much surface water directly drains into groundwater aquifers. Resolving this problem requires investigating resources on aquifer characterization, and monitoring program development.

Practical applications of the ecosystem services concept in groundwater management are still scarce (Tuinstra and van Wensem 2014) and specifics on groundwater dependent ecosystem management are lacking (Rohde et al 2017). Groundwater dependent ecosystems (GDE) are ecosystems that are maintained by direct or indirect access to groundwater. Few review articles (Klove et al 2014, Tuinstra and van Wensem 2014, De La Hera et al 2016, Rohde et al 2017) have discussed the potential contribution of the ecosystem services concept to groundwater sustainability. Here we extend this discussion by conducting a systemic review of articles that incorporated the ecosystem services concept in groundwater sustainability evaluation. Given the pool of articles from the screening stage, we identified 40 research articles related to GDE, with 26 articles including a groundwater modeling component. The eligibility criteria for inclusion in this review article is groundwater sustainability evaluation relevance, and all the articles were eligible to be included. In addition, we conducted a full-text assessment of eight review articles about GDE. In this section we first provide a brief overview about the concept of ecosystem services and its importance for groundwater sustainability, threats to ecosystem services due to unsustainable groundwater management practices, and the state of science in incorporating the ecosystem services concept in groundwater sustainability from an ecohydrology perspective. This is followed by a summary of our aforementioned systematic review, which mainly focuses on how the ecosystem services concept has been incorporated in groundwater sustainability evaluation.

Layering GDE on top of an inherently complicated SW-GW system adds another layer of complexity, which needs to be considered when evaluating a groundwater sustainability policy. GDE include springs, wetlands, peatlands, wet forests, lakes, rivers, streams, hyporheic zone, estuary, coastal lagoons, nearshore marine environments, and aquifers. GDE provide ecosystem services, which are benefits that people obtain from ecosystems. These include provisioning services (e.g. supply of water, food, water, and timber), regulating services (e.g. water purification, waste treatment, erosion regulation, flood control, drought buffer, and flora and fauna habitat), cultural services (e.g. provision of recreational, educational, research, aesthetic and spiritual benefits), and supporting services that are necessary for the production of all other ecosystem services such as water and nutrient cycling (Millennium Ecosystem Assessment 2005, Knüppe et al 2016). The reader is referred to other studies (Tuinstra and van Wensem 2014, Knüppe et al 2016) that provide a detailed description of ecosystem services that depend on groundwater. Groundwater is important for safeguarding ecosystem functions and services because groundwater helps to maintain water levels, temperature, oxygen content, dissolved ions, nutrients, organic matters, and (bio)geochemical conditions required by plants and animals in GDE (Kløve et al 2011, Gleeson and Richter 2018). For detailed discussion about the importance of various groundwater attributes to different classes of GDE the reader is referred to Eamus et al (2016).

The aforementioned ecosystem functions and services will be altered by groundwater pumping. Although this relation is poorly understood, a few articles reviewed the impact of groundwater pumping on wetlands (De La Hera et al 2016), streamflow (Gleeson and Richter 2018), and multiple GDE under climate change (Klove et al 2014). Yet there is a general understanding that groundwater capture can detrimentally impact ecosystems. Groundwater pumping and capture translate into (1) increased recharge rates, (2) decreased areas of rejected recharge where the water table is close to the ground surface that is vital for wetlands, (3) changes in groundwater pressure, (4) increases in drawdown, (5) decreases in streamflow, baseflow, spring discharge, and submarine groundwater discharge, and (6) changes in water quality and physio-chemical properties (Kalf and Woolley 2005, Konikow and Leake 2014, De La Hera et al 2016). For example, the impacts of increased drawdown are that phreatophytes will no longer able to reach the water-table (De La Hera et al 2016), and soil salinization may increase (Shi et al 2012). For details about the impacts of groundwater pumping and quality on GDE the reader is refereed to Eamus et al (2016).

Ecosystem services is an integral policy component in sustainable groundwater management. The aforementioned external pressures from the groundwater system can change the structure and functioning of GDE, cause habitat destruction, result in loss of sensitive species and biodiversity, and lead to gradual loss of ecosystem services as discussed in details by Klove et al (2014). This conflict between ecosystems and human needs is increasing, as in the past century human population quadrupled, irrigated agricultural land increased by six-fold, and water pumping rates from fresh water ecosystems increased eightfold (Sophocleous 2007). Thus, safeguarding ecosystem services is becoming mandatory under different groundwater sustainability policies. Rohde et al (2017)'s review article provides a global synthesis of managing GDE under groundwater sustainability policies with special focus on Australia, California, the European Union, and South Africa.

The science to quantify ecological water needs to support GDE is in its infancy (Wheater and Gober 2015, De La Hera et al 2016, Rohde et al 2017, Gleeson and Richter 2018). In Hawaii, for example, the impacts of submarine groundwater discharge on marine ecosystems (which are especially important to the indigenous Hawaiian culture and Hawaii's tourism-based economy) have been studied (e.g. Duarte et al 2010, Delevaux et al 2018). Yet we do not know the ecological water requirements, which is how much water the GDE needs. Thus, the groundwater sustainability estimates that would safeguard these ecosystem services have never been quantified, except as a trade-off curve between pumping and available fraction of the pre-development ecological water (Burnett et al 2020). In other words, a significant challenge is to determine the groundwater needs of GDE, and the acceptable amount of water that can be withdrawn from ecosystems without severally degrading the natural functioning and productivity of GDE (Sophocleous 2007, De La Hera et al 2016). This is a challenging task because how ecosystems depend on hydrological drivers, and how they respond to predicted changes in hydrology at different spatial and temporal scale issues is currently not well described, largely due to knowledge gaps in ecohydrology (Klove et al 2014, Rohde et al 2017). This is also due to sparsity of monitoring data that includes information on geomorphology, ecology and hydrology of ecosystems (Klove et al 2014). In addition, the diversity of groundwater dependent ecosystems (e.g. Howard and Merrifield 2010) makes it difficult to provide a one-size-fits-all management solution (Rohde et al 2017). Although most tools for managing groundwater- ecological systems are not yet well defined, De La Hera et al (2016) conclude that the integration of groundwater processes and ecohydrological concepts is currently possible, and may result in more sustainable outcomes.

How to incorporate ecosystem services to groundwater sustainability evaluation is unclear. To bridge the gap between scientific research and regulatory needs, Klove et al (2014) suggest establishment of a list of measurable indicators that describe GDE vulnerability at several spatial and temporal scales with multiple factors (physical, chemical and biotic); this is followed by identification of the linkages between these indicators and groundwater hydrology, climate change impacts, and human activities. In practice, to account for ecosystem services in groundwater sustainability evaluation, there are two related questions that can be assessed. These are: 'How much ecological water will be available given a certain pumping scheme? What are the ecological water needs of the GDE?' While the first question is relatively easy to answer, the second question is much more involved. Hydrological modeling can be used to answer the first question using multiple hydraulic, water budget and water quality indicators. For example, water-level is a commonly used indicator to account for ecological water demand (Stigter et al 2009, Shi et al 2012, Esteban and Dinar 2013, Li et al 2013, Pereau and Pryet 2018, Burnett et al 2020). Groundwater quality is another important indicator (Hinsby et al 2008, Burnett et al 2020). Rohde et al (2017) discussed the scientific underpinning for establishing indicators with thresholds, highlighting that the inherent diversity of groundwater dependent ecosystems requires indicators and thresholds to be locally determined, due to differences in reliance on groundwater, species composition, and adaptive capacities to varying threats. The challenge is how to define the threshold of these indicators. Multiple approaches are used to address this question.

One straightforward approach is to conduct scenario analysis to show the impact of different indicator thresholds. For example, Hu et al (2010) tested the impact of different agricultural and irrigation scenarios on ecological water, showing that improving agricultural practices can help groundwater recovery and increase ecological water supply. Similarly, Stigter et al 2009 used FEFLOW to test future scenarios of recharge due to climate change and water demand increases, showing the difficulty of meeting ecological water for GDE as the lowering of hydraulic heads would significantly reduce or stop spring discharge, and increase saltwater intrusion. Ronayne et al (2017) used MODFLOW with streamflow-routing package to evaluate groundwater sustainability for the natural base-case scenario, and managed aquifer recharge scenario, with the purpose of increasing streamflow to improve riverine species habitat in the South Platte River, Colorado. The study uses capture concepts to assess how pumping and recharge alter the head-dependent flows to and from the aquifer that in turn can help explain impacts on streamflow. A similar study (Scherberg et al 2014) shows that a managed aquifer recharge scenario increased groundwater discharge through springs and stream beds, benefiting aquatic habitat rather than building long-term aquifer storage.

A second approach to is to use insights gained from hydrological modeling to determine the indicator threshold. For example, Hsu et al (2012) in the Chih-Ben watershed, Taiwan used changes in groundwater levels that can lead to variation in soil-water content, and accordingly the deterioration of the GDE as an indicator. The indicator threshold is determined using mechanistic SW-GW model simulations that showed that groundwater sustainability evaluation is sensitive to any small drawdown in the downstream region. This causes a significant decrease of groundwater levels in the mountain areas impacting biodiversity and other ecosystem services. Also, Mulligan et al (2014) use average stream flow as an indicator and the threshold is set based on a percentage of the historic flow that would not cause infeasible solutions in the model runs.

A third approach is to use a simulation optimization method with trade-off evaluation of different indicator thresholds. Shi et al (2012) used a trade-off approach for meeting ecological environment water demand while balancing the integrated benefits of the society, economy, environment, and resources. This study evaluated groundwater sustainability under twelve strategies of four surface reservoir supply scenarios with low to high proportions of water crossed by three river ecological water demand thresholds of minimum (10% of multi-annual runoff), suitable (30% of multi-annual runoff), and optimal (60% of multi-annual runoff). By conducting trade-off evaluation and ranking these solutions, Shi et al (2012) show that the optimized groundwater yield could be sustained by high reservoir supplies to maintain suitable ecological water demand, while meeting human needs. Similarly, Burnett et al (2020) use a simulation optimization method to develop a trade-off curve of optimal groundwater pumping at different levels of spring discharge, which is important for terrestrial GDE.

A fourth approach to determine the indicator threshold is to estimate the water needs of GDE. The goal is to define measurable ecological indicators that assess GDE; and link these indicators to hydrological models to evaluate the groundwater sustainability accordingly. Monitoring and targeted research can identify which aspects of the hydrological regime are most important in supporting the structure and function of an ecosystem, and accordingly the ecological water requirements can be quantified using this best available scientific information (Rohde et al 2017). For example, following the guidelines from the Danish EPA with a participatory evaluation process, Henriksen et al (2008) formulate an ecological-based indicator focusing on baseflow and the ecological objectives of river reach. This indicator was used with other hydraulic and aquatic indicators to evaluate groundwater sustainability in Denmark using an integrated SW-GW via the model MIKE-SHE. Establishing these indicators with thresholds requires studying multiple ecological factors. For example, Aldous and Bach (2014) estimate groundwater sustainability by establishing hydro-ecological relationships that can be used to develop quantitative, measurable thresholds that are sensitive to changes in groundwater quantity and drawdown. Also, in response to the Massachusetts Water Management Act, Levangie (2008) shows that successive monitoring and research efforts resulted in moving from safe yield to groundwater sustainability policy through defining an indicator that focuses on habitat and streamflow requirements of the basin. Ecohydrological modeling can also be used to formulate these indicators. For example, Yang et al (2016) use field survey data and an eco-environment assessment model to develop ecohydrological indicators for the Tuwei River watershed in Shaanxi Province, China. The indicators are used to evaluate and rank multiple pumping schemes using an ensemble of mechanistic and phenomenological hydrological models that simulate water budget, groundwater levels, and water quality. Generally, selecting and formulating these indicators to maximize the amount and quality of information about the ecological integrity of a system, while minimizing the time and expense involved, remains a challenge (Klove et al 2014, De La Hera et al 2016).

Fifth, a pragmatic approach is to evaluate the indicator thresholds from a hydroeconomic modeling perspective. As explained in detail in section 4.3, hydroeconomic models can be used to find optimal groundwater pumping schemes that would maximize a utility function such as increasing ecosystem services, and to evaluate trade-offs between human and ecological water needs. In doing so, the optimum path of pumping rate does not necessarily converge to the recharge rate, but depends on the costs associated with GDE damages (Pereau and Pryet 2018). As an example, Esteban and Dinar (2013) used an linear utility function to estimate damages to GDE given drawdown and reduction in the flows that feed GDE. Commonly, hydroeconomic models use lumped phenomenological hydrological models. Pereau and Pryet (2018) reviewed hydroeconomic modeling studies highlighting that most of these studies can be considered an illustration of the water budget myth as they do not account for natural discharge in the water budget, which is important for sustainable groundwater management in the presence of ecosystem damages. Thus, the study of Pereau and Pryet (2018) developed a hydro-economic model that takes into account environmental flows both in the water budget, and in the utility function. Kahil et al (2016) note that a major gap in hydroeconomic modeling is the weak integration of physically-based representations of SW-GW interactions to inform complex basin scale policy choices. Henriksen et al (2008) emphasize the need for a mechanistic SW-GW model to allow for discretization and explicit representation in space and time of exchange flow components between saturated zone flows, unsaturated zone, overland flows, river reaches, and drains. On the other hand, other studies (Seward et al 2006, Seward 2010) emphasize the importance of phenomenological models particularly when an adaptive management approach is adopted.

Modeling is only one tool to evaluate ecosystem services indicators, and other practical approaches include adaptive management, risk assessment, and participation. The concept of adaptive management has a long history of application in ecosystem management (Knuppe and Pahl-Wostl 2011). This is mainly because we recognize that GDE are complex systems; our knowledge about GDE is imperfect; our ability to predict the characteristics and responses of GDE to groundwater pumping is limited; the outcomes of our groundwater management practices on GDE are uncertain; human activities and societal attitudes towards GDE are diverse and variable; our knowledge about GDE management needs to be developed incrementally; and more ecohydrological research and modeling studies alone will not help us to improve groundwater sustainability evaluation. The concept of adaptive management acknowledges these concerns, and adopts a pragmatic, flexible and experimental 'learning by doing' approach to management. Adaptive management enables water managers to allocate groundwater based on routine monitoring (e.g. groundwater levels, water quality metrics, instream flow criteria, and vegetation growth), targeted scientific investigations (e.g. plant water-use modeling, numerical groundwater modeling, and environmental or isotopic tracers), participation, and exploration programs to test potential pumping schemes (Seward et al 2006, Rohde et al 2017). This is to determine the hydrological conditions and thresholds required to maintain a GDE (Rohde et al 2017), to learn about the resilience and buffering capacity of the systems (Henriksen et al 2008), and to select prefered options (Seward et al 2006).

Several case studies show the effectiveness of the adaptive management approach when successfully implemented. For example, Rohde et al (2017) state that adaptive management is the core of Australia's approach for managing GDE. As such the Australian government commissioned the development of a practical 'GDE Toolbox' to assist Australian state agencies in revising conceptual models, identifying threshold responses, and identifying and managing GDE for water plans. The study of Rohde et al (2017) shows that New South Wales adopted the most specific and comprehensive statewide approach for managing GDE by combining the adaptive management strategies required in the Australia National Water Initiative with ecological valuation and risk assessment. This is to strategically target GDE with high value and high risk that require the greatest attention. This is especially needed when large uncertainty exists in the early management years, and limited financial resources are available for monitoring (Rohde et al 2017). Similar regulatory programs were undertaken by the state of Michigan, USA, the province of Ontario, Canada, and under the EU Water Framework Directive as discussed by Gleeson and Richter (2018). For example, to meet the California Sustainable Groundwater Management Act requirements for GDE, a guidance document (Rohde et al 2018) was developed based on best available science to help agencies, consultants, and stakeholders to efficiently incorporate GDE into groundwater sustainability plans. Participation is an integral part of adaptive management because indicator thresholds, and trade-offs cannot be decided on by science alone and community values need to be integrated (Seward 2010), and to enable stakeholders and decision makers to make risk‐based consensus decisions (Gleeson and Richter 2018).

While the concept of adaptive management may appear to be in conflict with the stance that sustainable groundwater volume is a fixed number or range that science can determine with more research, this is not the case. Seward et al (2006) clarify that the role of scientists is to identify a range of sustainability options with expected consequences, while the role of groundwater users and managers is to select the preferred option; then scientists can monitor the outcomes of that option and revise the sustainability scenarios as required. This is further elaborated in section 6.

Human activities add a third layer of complexity to the modeling framework. Sustainable water resources management requires consideration of human-groundwater systems interaction for several reasons. From the policy side, accounting for human activities is imperative for groundwater sustainability policy-making and implementation. For example, most of the articles related to human activities modeling, that were eligible to be included after our full-text review, address a policy issue. These articles addressed policies related to regulation of flow, quality, allocation, or pricing of groundwater (Li et al 2013, Cabello et al 2015, Gohar et al 2019), policies related to improving groundwater use efficiency and overdraft reduction (Li et al 2013, Davidsen et al 2016, Al-amin et al 2018, Feng et al 2018), adaptation policies (Kahil et al 2016, Wurl et al 2018), and groundwater replacement policies (Badiuzzaman et al 2017). These policy components also encompass evaluation and comparisons of the impact of different groundwater policy options and mechanisms (Aistrup et al 2017, Esteban and Dinar 2013, Farhadi et al 2016, Harou and Lund 2008, Hyndman et al 2017, Zellner and Reeves 2012), groundwater policy guidance and development (Zellner 2007, Hu et al 2010, Shi et al 2012, Fernald et al 2015, Salem et al 2017), and evaluation of modeling approaches related to groundwater sustainability policy (Mulligan et al 2014, Macewan et al 2017).

From the science side, limited characterization of the human-groundwater systems interaction has resulted in inadequate groundwater sustainability evaluation. This is mainly because the groundwater sustainability is a function of both aquifer performance and governance factors (figure 3). As these factors are changing through space and time, it is becoming increasingly clear that improved accounting of changes, interactions, and feedback within coupled human-water resources systems is important (Montanari et al 2013, Vogel et al 2015). Moreover, human activities have a big impact on model simulation and predictions. Vogel et al (2015) challenge the emphasis on natural processes as the main predictor of hydrological response, given the integration of human systems into hydrological design. Accounting for human activities is important for studying the improvement, efficiency, and sustainability of groundwater resources systems because societal and land-use changes can cause fundamental shifts in basin functions that may be of the same order of magnitude, if not greater, than those predicted by climate change (Ceola et al 2016, Ferguson and Maxwell 2012). While mathematical models about natural aspects of groundwater resources systems have been a matter of research for centuries and are generally well developed, human activities and behavioral models in hydrology are much newer and less developed (Giuliani et al 2016).

In this section, we review human activities modeling in the groundwater sustainability literature, which is relatively limited in comparison to natural aspects, although rapidly expanding (figure 10). We conducted a full-text review of 46 articles about groundwater sustainability at different geographic locations that describe groundwater systems and human activities. After our full-text review, 36 articles were eligible to be included in the review article based on the relevance of the work. For example, we excluded articles that are not directly related to groundwater sustainability, and articles that are mainly economics-driven rather than hydrology-driven. We classified these articles under the three research areas of integrated water resources management, hydroeconomic modeling, and sociohydrology resulting in 8, 11, and 17 articles respectively. This section summarizes human activities modeling approaches that are used in these studies. These include land-use scenarios and backcasting, simulation optimization methods, hydroeconomic optimization, regional economic-engineering models, system dynamics modeling, holistic natural-human models, and agent-based modeling.

Several studies in hydrosociology and integrated water resources management address human activities modeling for sustainable groundwater development. Hydrosociology is the study of the impacts of human activities on water resources systems and the societal impacts of water resources projects (Falkenmark 1979). Integrated water resources management extends hydrosociology. The Global Water Partnership (2000) defines integrated water resource management as 'a process which promotes the coordinated development and management of water, land, and related resources in order to maximize the resultant economic and social welfare in an equitable manner without compromising the sustainability of vital eco-systems.' The eight reviewed articles under these two research sub-fields use mechanistic groundwater models to relate human activities to groundwater systems. While MODFLOW is commonly used (Hu et al 2010, Li et al 2013, Hyndman et al 2017, Wurl et al 2018), Shi et al (2012) use coupled SWAT-MODFLOW, Feng et al (2018) use GSFLOW, and Yin et al (2012) use MIKE-SHE.

In these studies, sustainable groundwater development is evaluated such that human activities are modeled using scenario analysis and simulation optimization. The scenario analysis approach estimates how the groundwater system responds to different scenarios related to human activities. These scenarios represent changes in agricultural and domestic water demand (Ramesh and Mahesha 2008), regulations to reduce overexploitation (Li et al 2013, Feng et al 2018), human activities in relation to climate conditions (Pulido-Velazquez et al 2015, Wurl et al 2018, Feng et al 2018), land-use (Pulido-Velazquez et al 2015, Hyndman et al 2017), and irrigation and soil management practices (Hu et al 2010, Shu et al 2012). These scenarios were developed using an array of approaches, from simple to complex. As a simple approach, Feng et al (2018) use an attribution approach based on historical trend analysis of coupled SW-GW model simulations to develop scenarios of upstream inflow and groundwater pumping that are a function of climate conditions and human activities. While agricultural water use in this study was determined externally and assumed unchanged under future climate conditions, Ramesh and Mahesha (2008) estimate increase in agricultural and domestic water demand based on population forecasting. This was used to evaluate sustainable groundwater development under different development scenarios that are a function of agricultural and domestic water demand and recharge. Similarly, Li et al (2013) define sustainable management practices by using a ground-settlement recovery scenario and groundwater exploitation planning scenarios generated from local development plans. Wurl et al (2018) adopt a more elaborate stakeholder-driven approach to conduct PESTE (political, economic, social, technological, and environmental) and SWOT (strengths, weaknesses, opportunities, and threats) analyses to define water insecurity scenarios and sustainable pumping scenarios. More elaborate dynamic scenarios can also be developed. For example, Shu et al (2012) adopt an integrated modeling approach using MIKE-SHE to conduct a scenario analysis of the effect of different crop rotations, irrigation intensity, water transfer projects, and other water management options. Also, Hu et al (2010) use a crop-growth model to develop agricultural water-saving scenarios with agronomic factors such as water use efficiency and crop productivity to ensure groundwater recovery. Similarly, Hyndman et al (2017) use SALUS (System Approach to Land-use Sustainability) to develop dynamic land-use scenarios to study groundwater sustainability in a megacity. Pulido-Velazquez et al (2015) sequentially couple SWAT and MODFLOW with a nitrate mass-transport model in MT3DMS to evaluate the impact of climate and land-use change on groundwater sustainability strategies such that land-use scenarios were developed based on historical analysis and by accounting for multiple socioeconomic factors.

Whereas scenario analysis provides a solution given a scenario, the purpose of simulation optimization is to maximize a target objective. Shi et al (2012) evaluate sustainable groundwater development for a river basin using an integrated simulation optimization method that assesses the integral benefit of water resources development and utilization by maximizing water use efficiency, optimizing environmental water demand, and minimizing anthropogenic influence on groundwater systems. This is carried out using an integrated hydrological model coupled with a benefit-cost-loss model, which accounts for multiple aquifer performance and governance factors using cost and trade-off functions (e.g. substitute expense method, shadow engineering method, allocation coefficient method, opportunity cost method, and market value method). These factors include benefits (industrial, irrigation and domestic water supply, flood control, inter-generational equity, tourism, and culture), cost (surface water and groundwater), and loss (groundwater and surface pollution, phreatic evaporation, soil salinization, increase exploitation cost, river sedimentation, impact on farmland and people, reservoir surface evaporation).

We reviewed 21 hydroeconomic studies and 11 were eligible to be included based on relevance. We only included hydroeconomic modeling studies and not hydroeconomic indicator studies. For example, Livingston and Garrido (2004) suggest a set of hydroeconomic indicators for assessing the relative success or failure of groundwater policies that address emerging physical, economic, and environmental stress on groundwater resources. Other similar studies that are not modeling-related were excluded. We also excluded studies that utilize global hydrological models to study groundwater sustainability at a large scale (Zaveri et al 2016).

Hydroeconomic modeling (Noel and Howitt 1982, Harou et al 2009) aims to optimize the economic objectives of a water system subject to natural or societal constraints, or to evaluate optimal water conservation options and development projects (Pande and Sivapalan 2017). Hydroeconomic models attempt to capture the complexity of water users' economic decisions and the biophysical constraints of the water system (Macewan et al 2017). Generally, the overarching goal of an economically efficient (optimized) water policy is to ensure that water can be sustainably provided to meet users' needs at an affordable price (Kahil et al 2016). Depending on the groundwater sustainability objectives, the hydroeconomic model that is linked to the hydrological model could be a simple crop or water demand model addressing one or a number of human activities related to groundwater sustainability such as agriculture (Harou and Lund 2008, Valderrama et al 2011, Esteban and Dinar 2013, Alcala et al 2015, Davidsen et al 2016, Salem et al 2017, Aistrup et al 2017, Gohar et al 2019), urban water demand (Harou and Lund 2008, Gohar et al 2019) and tourism (Alcala et al 2015). Hydroeconomic models could use analytical functions to account for ecosystem services (Harou and Lund 2008, Esteban and Dinar 2013, Mulligan et al 2014, Kahil et al 2016, Aistrup et al 2017). Some of these models could be more comprehensive model accounting for multiple aquifer performance and governance components (Mulligan et al 2014, Kahil et al 2016, Macewan et al 2017).

Existing economic models could take different forms such as economic optimization objective functions with analytical constraints (Esteban and Dinar 2013, Davidsen et al 2016, Kahil et al 2016, Salem et al 2017, Gohar et al 2019), dynamics modeling causal diagrams with analytical functions representing hydrological systems and human activities (Valderrama et al 2011, Alcala et al 2015), economic optimization agent-based model (Mulligan et al 2014), and Open Modeling Interface Standard (Aistrup et al 2017). Note that a system dynamics model uses a causal loop diagram to represent the constituent components of the system and their interactions including feedback loops to determine the system's behavior over a period of time. An agent-based model simulates the behaviors, actions and interaction of individual and collective agents to assess their effects on the system as a whole. The Open Modeling Interface Standard (Gregersen et al 2007) is an interface that facilitates model integration by allowing existing models to run simultaneously and share information. In addition, more comprehensive region-specific economic models are used such as California's CALVIN regional economic-engineering model (Harou and Lund 2008) and California's Statewide Agricultural Production (SWAP) economic model (Macewan et al 2017).

Some of the aforementioned hydroeconomic modeling studies use phenomenological groundwater models such as lumped water budget models (Valderrama et al 2011, Esteban and Dinar 2013, Alcala et al 2015, Davidsen et al 2016, Aistrup et al 2017, Gohar et al 2019), the economic-engineering optimization model of California CALVIN (Harou and Lund 2008), response functions generated from California's C2VSim groundwater-surface water simulation model (Macewan et al 2017), and machine learning models (Salem et al 2017). Alternatively, Kahil et al (2016) and Mulligan et al (2014) use the distributed mechanistic model of MODFLOW. This approach more accurately simulates the spatial and temporal heterogeneity aquifer and interaction with river flow (Kahil et al 2016). In addition, Mulligan et al (2014) compare lumped and distributed models to show that a detailed distributed mechanistic model with spatially heterogeneous representation (rather than a lumped groundwater model) allows for spatially disaggregated human activities representation, which is necessary to properly investigate policy instruments in a groundwater basin. Such a fully-coupled approach between human-groundwater system is generally a more robust framework to support water management decisions when optimal groundwater policy has multiple objectives with complex interactions with environmental variables (Macewan et al 2017). However, phenomenological models could still be advantageous when mechanistic understanding of the aquifer system is not possible due to data scarcity. In addition, calibration of fully-coupled human-groundwater hydroeconomic model could be difficult, and data exchange between two systems is not generally straightforward. Data exchange between hydrological and socioeconomic systems can be bidirectional with system dynamics modeling (Valderrama et al 2011, Alcala et al 2015), agent-base modeling (Mulligan et al 2014), and Open Modeling Interface Standard (Aistrup et al 2017). This bidirectional feedback, which allows for developing coevolving scenarios between human-water system, is one of the emerging ideas and ways forward for the progress of hydrology (Montanari et al 2015).

With respect to sustainable groundwater development, while hydroeconomic modeling operationalizes economic concepts and incorporates them into water management to find feasible and optimal solutions, sociohydrology integrates the human-dimension as an endogenous component of water systems to develop an understanding of the dynamics and bidirectional feedback of the coupled human-water system (Sivapalan and Blöschl 2015, Wheater and Gober 2015).

Hydrology is relevant to society, and thus requires integration with the natural and social sciences to better support water management decisions (Wagener et al 2010, Lund 2015, Vogel et al 2015). With respect to sustainable groundwater development, hydroeconomic models are mainly normative models aiming at operationalizing and incorporating economic concepts into water management to maximize a utility function to find feasible and optimal solutions (Giuliani et al 2016). On the other hand, sociohydrology is more descriptive, aiming at integrating the human-dimension as an endogenous component of the water system. This is to develop an understanding about the dynamics, bidirectional feedback, and coevolving scenarios of the coupled human-water system in the past and when exposed to altered boundary conditions such as climate or other socioeconomic drivers (Sivapalan and Blöschl 2015, Wheater and Gober 2015). As noted by Wesselink et al (2017), the hydrological science community has recently launched sociohydrology as the research theme for the current decade (2013–2022) to advance hydrological science for the benefit of society. As such, the water system acts as a changing interface between environment and society (Montanari et al 2013). This poses a necessity and a challenge for water sustainability analysis to explore the endogenization of human-dimension (e.g. values, societal preferences, technology, economics, etc) in space and in space-time (Pande and Sivapalan 2017). A sociohydrological system is generally composed of four sub-systems that are the hydrological, ecological, economic, and societal sub-systems (Liu et al 2015). For details about sociohydrological system conceptualization and modeling the reader is referred to recent reviews (Sivapalan and Blöschl 2015, Blair and Buytaert 2016, Roobavannan et al 2018).

As the science of sociohydrological modeling is emerging, so are sustainable groundwater development studies from a sociohydrology prospective. We reviewed 17 studies that we divide into two categories. The first category includes articles that provide a sociohydrological perspective about the human-water system of the study area over a long period of time. These articles do not necessarily include a modeling component, but rather present a comprehensive overview based on reviewing and synthesizing existing data and literature. By considering the changes in values and norms governing the aquifer system along with natural events such as increase or decrease of precipitation, the aim of these studies is to describe changes in the structure, dynamics and different stages of sustainable groundwater development thresholds and tipping points. Thus, these studies provide insights on how groundwater resources react to pressures from human activities and in turn how the society reacts to threats to groundwater resources, with examples from the Indus basin in Pakistan (Archer et al 2010), the irrigated valleys of northern New Mexico (Fernald et al 2015), the North China Plain (Han et al 2017), the Mediterranean region (Leduc et al 2017), the agricultural region in Central California (Rudestam et al 2018), and the Murray-Darling Basin in Australia (Sanderson and Curtis 2016). This could also include using a human-ecosystem-water modeling framework with a water accounting method such as the concept of water metabolism (Cabello et al 2015).

The second category includes articles about coupled human-water system modeling to study sustainable groundwater development. Diverse coupled human-water system modeling approaches are used in sociohydrology, including agent-based modeling, coupled component modeling, system dynamics, game theory, Bayesian network, and pattern-oriented modeling (among others) as reviewed by Blair and Buytaert (2016). The reviewed groundwater sustainable development studies (Zellner 2007, Zellner and Reeves 2012, Farhadi et al 2016, Al-amin et al 2018) use agent-based modeling coupled with MODFLOW mechanistic groundwater flow models. Castilla-Rho et al (2019) use agent-based modeling with FlowLogo (Castilla-Rho et al 2015), which is finite-difference solution of the governing equations of groundwater flow with agent-based simulation. Two studies in Michigan (Zellner 2007, Zellner and Reeves 2012) use an agent-based land-use model to study sustainable groundwater development with respect to urban growth by suggesting alternative forms of development and water-use to understand how these processes interact to create the observed patterns of water resource depletion and sustainability. Al-amin et al (2018) develop an agent-based model for Arizona that simulates water demands to capture the dynamic interactions among household-level consumers and policy makers to implement mandatory restriction policies for groundwater sustainability. Farhadi et al (2016) use an agent-based model in Iran with a multi-objective simulation optimization method that is informed by stakeholders. This is to evaluate groundwater sustainability policy mechanisms that encourage agents to cooperate with the management decisions. Castilla-Rho et al (2019) investigate case studies in the Murray–Darling Basin (Australia), the California Central Valley (USA), and the transboundary Punjab aquifer (India and Pakistan) to show that effective groundwater sustainability regulations and their implementation need to account for culture.

From the above analysis we can conclude that unlike classic Integrated Water Resource Management (IWRM) approaches that mainly focus on the one-way relationship between human-water systems, sociohydrology can provide more insights through a two-way dynamic relationship. This might not be useful in all cases. Yet such involved analysis is particularly important when human factors could have a larger impact on groundwater sustainability solutions in comparison to other factors such as climate and geology. In that case, descriptive sociohydrological models can provide insights into the dynamic interaction between human-water systems, and help to find better sustainable groundwater development solutions (Zellner 2007, Zellner and Reeves 2012, Farhadi et al 2016, Al-amin et al 2018). In addition, a few studies turn descriptive sociohydrological tools to normative approaches whose ultimate goal is to predict the optimal human decisions toward sustainable groundwater development using both phenomenological models (Valderrama et al 2011, Alcala et al 2015, Aistrup et al 2017) and mechanistic models (Mulligan et al 2014, Farhadi et al 2016). Using a mechanistic model can provide site-specific details, while using a phenomenological model can provide a regional economic context and an understanding about macro-level properties, functioning, dynamic interactions, pathways, and sustainability tipping points of human-water systems. However, major challenges to human activity modeling studies stem from the profound uncertainties associated with both the natural and human aspects of the complex and dynamic water systems (Wheater and Gober 2015, Blair and Buytaert 2016) as discussed in the following section.

Modeling is a process by which we communicate (what we think to be) our knowledge about the system, and uncertainty analysis is a process by which we communicate our incomplete knowledge about the system. In groundwater sustainability literature, several studies identified and discussed multiple sources of uncertainty without addressing all or any of them (e.g. El-Kadi et al 2014, Hu et al 2010, Lathashri and Mahesha 2016, Urrutia et al 2018, Piscopo et al 2019). As such, uncertainty analysis is an excellent tool for the modelers to understand and communicate model uncertainty and limitations. For example, Piscopo et al (2019) identify sources of uncertainty to caution that a proposed groundwater management plan is preliminary since uncertainties of the system are still hardly quantifiable. Also, communicating to end-users what the model can and cannot do reduces the risk of model misuse. For example, Hugman et al (2013) caution that conceptualizing karstic aquifer systems with a single continuum equivalent porous media may result in significant uncertainty when simulating smaller scale effects such as locations of well fields. Even if uncertainty is not quantified, a qualitative analysis of uncertainty can be informative (Gillespie et al 2012) and helpful to mitigate uncertainty (Gallardo et al 2009). For example, Gillespie et al (2012) propose two conceptual models representing local recharge or intra-basin flow conceptualizations and discuss the limitations of the two models with implications on aquifer sustainability. As another example, Gallardo et al (2009) discuss the model limitations and suggest a safety factor of 20% to account for uncertainty.

Generating different predictive, explorative, and normative scenarios (Börjeson et al 2006) are the simplest and most commonly used approaches for uncertainty analysis. In the context of groundwater sustainability, these include pumping scenarios (e.g. Stigter et al 2009, Calderhead et al 2012, Hugman et al 2013, Lathashri and Mahesha 2016, Macewan et al 2017, Feng et al 2018, Urrutia et al 2018, Gohar et al 2019), recharge scenarios (e.g. Stigter et al 2009, Calderhead et al 2012, El-Kadi et al 2014, Lathashri and Mahesha 2016, Feng et al 2018, Urrutia et al 2018, Wurl et al 2018, Gohar et al 2019, Pholkern et al 2019), sea-level rise scenarios (e.g. Lathashri and Mahesha 2016), water export scenarios (Calderhead et al 2012), land-use scenarios (e.g. Passarello et al 2014), 'what-if' stakeholders preferences scenarios (Timani and Peralta 2015), and policy scenarios (e.g. Aistrup et al 2017, Guillaume et al 2012), among others. While scenarios can be generated for multiple purposes, as discussed below, we first present the simplest case of generating scenarios, which is to demonstrate future possibilities and to generate discussion. For example, El-Kadi et al (2014) conduct a study in Jeju Island, South Korea to show that if an historical drought occurred in future, it would decrease sustainable yield by 16%, decrease spring discharge by 28%, and dry up 27% of springs in comparison to the baseline case. Several studies develop scenarios to advocate for the need of adaptation measures and new polices. For example, Urrutia et al (2018) present recharge scenarios due to climate change and pumping scenarios due to increase in population and mining activities in the Atacama Desert in northern Chile to demonstrate the need to use alternative water resources such as desalination to minimize the impact of the combined effects of economic growth and climate change on the aquifer. Similarly, Stigter et al (2009) demonstrate the combined impact from climate change and water demand increases on aquifer discharge and risk of ecological degradation in Algarve, Portugal; this is to point out to stakeholders the need to prepare societal and technical tools to alleviate these impacts, and to broaden the current definitions of sustainability in the study region. Other studies additionally advocate for the need for new approaches for decision support under uncertainty. For example, Passarello et al (2014) develop land-use change scenarios for an urban area in Texas to demonstrate to stakeholders the influence of urbanization and the implications of these scientific uncertainties on policy and urban water management decisions. Similarly, Gober et al (2010) integrate climate change uncertainty into formal decision analysis for water planning to offer insights into water planning in Phoenix, Arizona, and to demonstrate the need for new approaches to decision making under uncertainty.

As Kitanidis (2015) state: 'The sooner we arrive at this realization [that all models are wrong] the better, so that we can either apply a model appropriately or move on to the next step, uncertainty quantification.' To quantify the impact of parametric, model, and scenario uncertainty on model predictions is the most common task in uncertainty analysis. For example, Henriksen et al (2008) estimate the total exploitable groundwater resource of Denmark to be 1 × 10⁹ m³/year and provide an uncertainty estimate for each regional estimate, which ranges from ±10% to ±40%. Uncertainty is ubiquitous in management models to evaluate groundwater sustainability. Groundwater modelers are typically faced with (1) complex subsurface heterogeneity, (2) state variables and parameters that are scale, spatial and time dependent, (3) data scarcity about subsurface geology, (4) uncertainty about top and inland boundary recharge and boundary conditions, and (5) computationally intensive numerical models that generally hinder uncertainty quantification. Despite these challenges, systematic investigation into the uncertainty quantification and its impact on decision support has been limited in groundwater hydrology (Heße et al 2019). This is especially true in groundwater sustainability literature. Advancements in uncertainty quantification in subsurface hydrology have been reviewed and discussed for model data integration (Rajabi et al 2018), uncertainty of subsurface characterization (Scheidt et al 2018), and conceptual uncertainty (Enemark et al 2019). However, the state-of-the-art techniques and tools highlighted in these studies are not widely applied in groundwater sustainability numerical modeling studies. Deterministic groundwater models are often used, and uncertainty quantification in the model outputs is generally assessed using sensitivity analysis (Hu et al 2010, Huang et al 2012, Calderhead et al 2012, El-Kadi et al 2014, Pholkern et al 2019), and multiple deterministic conceptual models of the subsurface (Timani and Peralta 2015, Feng et al 2018, Pholkern et al 2019). Also, Lal and Datta (2019) use an ensemble surrogate model within a simulation optimization framework to address parametric uncertainty of hydraulic conductivity and porosity. Sources of uncertainty discussed in these articles regarding parametric and conceptual uncertainty include subsurface geology, hydraulic conductivity, anisotropy ratio, specific yield, recharge, boundary condition and fluxes, pumping rates and locations, riverbed hydraulic conductance, porosity, dispersivity, etc. Qin et al (2013) consider uncertainty from using a coarse model grid. Delottier et al (2017) consider uncertainty in steady versus transient solutions. Out of the reviewed numerical modeling case studies about groundwater sustainability, only two studies (Mustafa et al 2018, 2019) use a probabilistic single model ensemble with Markov chain Monte Carlo (MCMC) sampling to address parametric uncertainty, and Mustafa et al (2019) use a probabilistic multi-model ensemble with Bayesian model averaging to address parametric and conceptual uncertainty. While MCMC is not used due to its high computational cost, Monte Carlo simulation approaches with low fidelity (Hill et al 2016) or multi-fidelity (Zhang et al 2018) can be adopted to reduce the computational burdens. Other alternatives include using parallel computing (e.g. Elshall et al 2015) or surrogate models (e.g. Zhang et al 2013). More case studies that use these and similar tools are needed to advance the science of groundwater sustainability evaluation.

Groundwater models are only one layer of uncertainty in water-ecology-human models. Uncertainty analysis in water-ecology-human systems has received less attention. While several articles have discussed the uncertainty of human-water systems (Guillaume et al 2012, Di Baldassarre et al 2016), ecosystems (Hamel and Bryant 2017), water-ecosystems (Beven and Alcock 2012), future scenarios (Maier et al 2016), and decision support (Reichert et al 2015), not much has been done in practice in the groundwater sustainability literature. Uncertainty of ecosystem services tied to groundwater sustainability is discussed in section 5.1.3 within the concept of adaptive management. The challenges and solutions of uncertainty quantification in ecosystem services are discussed by Hamel and Bryant (2017) and can be helpful to water-ecosystem models in hydrology. Regarding human activities, few of the reviewed studies consider several sources of uncertainty in their human-water models. These sources include weather factors (Valderrama et al 2011, Aistrup et al 2017), economic factors (Valderrama et al 2011, Susnik et al 2013), agent productivity parameters (Mulligan et al 2014), and human behaviors (Noel and Cai 2017). Using a phenomenological model, Guillaume et al (2012) conduct an uncertainty analysis of a dynamic coupled economic-groundwater model for groundwater sustainability evaluation. Sources of uncertainty considered in this study include allocation policies in future planning, processes of interest, rainfall variability, spatial pumping distribution, transmissivity and storativity of the aquifer, irrigation choices, agricultural price models of local conditions and crop yield parameters, relevant parameters for decisions, and model output uncertainty. Guillaume et al (2012) demonstrate how a variety of uncertainties in such a model can be addressed with a number of methods including propagation of scenarios and bounds on parameters, multiple models, block bootstrap time-series sampling, and robust linear regression for model calibration. Guillaume et al (2012) also provide an uncertainty typology for coupled human-water models, which can help advance this under-researched area.

Uncertainty analysis is a model diagnostic tool. For example, Tsai and Elshall (2013) develop a hierarchical Bayesian model averaging method that segregates uncertain model components; this is to comparatively evaluate the candidate propositions of each uncertain model component, to understand the individual contribution of each uncertain model component to the model prediction and variance, and to prioritize the contribution of each uncertain model component to the overall model uncertainty. Similarly, Dai et al (2017) develop hierarchical sensitivity analysis to identify important system processes under conceptual and parametric uncertainty. These and similar methods serve as a learning tool to advance knowledge about the model. These methods mainly involve combinatorial design to represent the uncertain model components. In groundwater sustainability literature, Calderhead et al (2012), for example, use multiple scenarios with combinatorial design to represent several uncertain model components, and show that the impact of climate change on recharge plays only a minor role in the occurrence of subsidence in Mexico City in comparison to pumping scenarios and groundwater export. Similarly, Pholkern et al (2019) in Northeast Thailand show that variable depths and thicknesses of the aquifer have a higher impact on sustainable yield estimates than model boundary conditions do. A main limitation of these and other similar studies (e.g. Hugman et al 2013, Rasmussen et al 2013, Unsal et al 2014, Zhao et al 2016) is that these conclusions are based on the comparative analysis of the results of different model components without accounting for the probability of these components and their interaction (e.g. Dai et al 2017). Feng et al (2018) use an attribution approach to study interaction between climatic and human impacts on groundwater sustainability for a coastal aquifer in northern China, leading to several insights related to aquifer function that call for strict regulations on groundwater pumping.

The identified important sources of uncertainty may be used to facilitate a data-worth analysis, which aims at designing data collection plans such that the expected benefit of new information exceeds its cost. In hydrology this is generally done using a Bayesian framework to identify new data locations (e.g. Pham and Tsai 2015, Neuman et al 2012, Lu et al 2012) or types (e.g. Wöhling et al 2015). We did not identify a similar study in groundwater sustainability literature. However, Timani and Peralta (2015) in Utah use a multi-model simulation optimization approach to reconcile two disparate conceptual models that are contested among stakeholders, and show about a 25% difference of maximum perennial-yield. With the simulation optimization procedure for the two disparate models, Timani and Peralta (2015) identify field data that is most needed to resolve this conflict. Additionally, Li et al (2014) extend the data-worth to information-worth analysis using a numerical groundwater model with different representations of information about the aquifer and its risk of contamination. Li et al (2014) assess the effectiveness of aquifer monitoring information in achieving more sustainable use, showing that pumping rates differ when risk information that synthesizes data on aquifer conditions is provided to the users, and that the level of information about the state of the aquifer also effects extraction behavior. The study highlights the importance of contamination data, showing that pumping is significantly reduced in experiments where contamination is possible compared to those where the pumping cost is the only factor discouraging groundwater use.

The reliability of a plan or design can change if uncertainty is introduced to the problem. For example, Chitsazan et al (2015) combine chance‐constrained programming with Bayesian model averaging to assess the impact of geological structure uncertainty in groundwater quality control design in comparison to traditional chance‐constrained programming; the study shows that considering parametric uncertainty alone overestimates the design reliability. Uncertainty analysis can help identify robust plans given uncertainty. Guillaume et al (2012) explain that uncertainty analysis provides an answer to the question: What if a model assumption is wrong?, and hence allows stakeholders to choose a policy with an understanding of the possible adverse impacts, or that will provide the desired outcome if the 'best assumption' is changed. In groundwater sustainability literature, several studies consider uncertainty analysis to identify robust plans given uncertainty (e.g. Guillaume et al 2012, Mulligan et al 2014, Uddameri et al 2014, Gohar et al 2019). For example, Uddameri et al (2014) use a fuzzy simulation optimization approach to identify a better policy to cope with the uncertainty regarding specifying desired future conditions due to incomplete understanding of the aquifer dynamics in South Texas. This fuzzy approach yielded lower estimates of groundwater availability in comparison to the crisp optimization scheme, as it accounts for stakeholders' uncertainty. Also, Mulligan et al (2014) compare two groundwater-use policies in California under the frame of productivity uncertainty, and use two modeling approaches to explore the effect of modeling assumptions on the projected performance of these polices.

A critical model parameter realization(s), critical model(s) or critical scenario(s) refers to the ones that have the most influential effect on the solution depending on the desired reliability level (Kourakos and Mantoglou 2008). In groundwater sustainability literature, Wurl et al (2018) offer an analysis of hydrological resilience of a water-limited arid ecosystem in northern Mexico, under future pumping scenarios and changing climate conditions; the study aims to recognize water insecurity scenarios and to define appropriate actions towards more sustainable groundwater use through involvement of local stakeholders. The identified critical scenario or model can be then used for further analysis. For example, Ostad-Ali-Askari et al (2019) evaluate several pumping and agricultural practice strategies to restore the aquifer sustainability in the study area, given the identified critical scenario of climate change. Seward et al (2006) note that identifying which conceptual model to be examined must be done in consultation with all the stakeholders.

While most of the abovementioned studies use uncertainty analysis as a learning tool to learn more about the model, the model solution, and the problem addressed by the model, other uncertainty analysis methods are specifically tailored to provide deeper insights about the problem that the model addresses. For example, Susnik et al (2013) compare system dynamics modeling and object-oriented Bayesian networks modelling to support groundwater management decision in the Kairouan aquifer system, Tunisia. System dynamics modeling (e.g. Calderhead et al 2012, Susnik et al 2013, Fernald et al 2015) implicitly accounts for uncertainty (and probabilistic uncertainty can be incorporated into it), while Bayesian networks are mainly a probabilistic framework. By comparing these two modelling paradigms, Susnik et al (2013) show that system dynamics modeling is a cyclic approach that allows the user to discover potentially hidden dynamics in a system by simulating non-linear feedback processes, while Bayesian networks is an acyclic approach that incorporates the variability and uncertainty in every single variable with probabilistic outputs for key variables. Susnik et al (2013) define a hidden dynamic as a behavior that emerges due to the interaction of all model components, which is not necessarily apparent from studying each model element independently. The study finds that the analyses of both models agree, indicating current overexploitation of the aquifer, and that pumping reduction offers the best solution to end aquifer overexploitation. However, the study notes that system dynamics modeling has the potential for stakeholder collaboration, while Bayesian networks can be overly complex as understanding of probabilistic distributions may not be straightforward for stakeholders.