Groundwater extraction may drown mega-delta: projections of extraction-induced subsidence and elevation of the Mekong delta for the 21st century

The low-lying and populous Vietnamese Mekong delta is rapidly losing elevation due to accelerating subsidence rates, primarily caused by increasing groundwater extraction. This strongly increases the delta’s vulnerability to flooding, salinization, coastal erosion and, ultimately, threatens its nearly 18 million inhabitants with permanent inundation. We present projections of extraction-induced subsidence and consequent delta elevation loss for this century following six mitigation and non-mitigation extraction scenarios using a 3D hydrogeological model with a coupled geotechnical module. Our results reveal the long-term physically response of the aquifer system following different groundwater extraction pathways and show the potential of the hydrogeological system to recover. When groundwater extraction is allowed to increase continuously, as it did over the past decades, extraction-induced subsidence has the potential to drown the Mekong delta single-handedly before the end of the century. Our quantifications also disclose the mitigation potential to reduce subsidence by limiting groundwater exploitation and hereby limiting future elevation loss. However, the window to mitigate is rapidly closing as large parts of the lowly elevated delta plain may already fall below sea level in the coming decades. Failure to mitigate groundwater extraction-induced subsidence may result in mass displacement of millions of people and could severely affect regional food security as the food producing capacity of the delta may collapse.


Introduction
The world's third largest delta, the populous and low-lying Mekong delta in Vietnam is facing increased river . On top of that, like many other deltas in the world (Syvitski et al 2009, Nicholls andCazenave 2010), the Mekong delta experiences accelerating rates of relative sea-level rise, the combined effect of absolute sea-level rise and land subsidence. As the Mekong delta has one of the lowest delta plains in the world, on average only ∼0.8 m above local m.s.l. (Minderhoud et al 2019a), relative sea-level rise threatens the delta's nearly 18 million inhabitants and its important economic function as an environment for agri-and aquaculture production vital to South-East Asia's food production. Land subsidence is the main source of relative sea-level rise in the Mekong delta and is caused by various driving processes: i.e. natural processes like tectonics and natural compaction of the Holocene sediments (Zoccarato et al 2018) and human-induced processes driven by amongst others groundwater extraction (Erban et al 2014, Minderhoud et al 2017, drainage of shallow sediments and loading by buildings and infrastructure . The change of Vietnam to an open-market economy in 1986 (Seto 2011) was the onset of large-scale groundwater extraction in the Mekong delta that triggered extraction-induced subsidence (Minderhoud et al 2017). In the decades that followed, groundwater extraction steadily increased, providing high-quality fresh water to meet the growing agricultural, industrial and domestic demand that fueled the rapid growing economy (Wagner et al 2012). As a result, extraction-induced subsidence has accelerated over the past decades to rates exceeding 25 mm yr −1 in certain areas, making it at present the main contributor of delta-wide subsidence in the Mekong delta (Minderhoud et al 2017).
With its low elevation and accelerating sea-level rise, the Mekong delta is on the verge of a tipping point and delta management choices, including policy decisions on groundwater use, will shape the future of the delta. Developing and implementation of integrated delta policies aimed to safeguard the delta for future generations are more important than ever before and require detailed quantifications of potential future subsidence. In contrast to subsidence caused by natural processes, e.g. natural sediment compaction or tectonics, which cannot be mitigated, extraction-induced subsidence can be targeted for mitigation in order to reduce subsidence rate and future delta-elevation loss. Until now, future projections of extraction-induced subsidence were solely based on linear extrapolations of past or present subsidence rates, however, these are unable to capture, for example, the influence of mitigation measures. Furthermore, subsidence induced by groundwater extraction is a nonlinear process that can exhibit delayed responses. For example, it can take years until the effect of a change in hydraulic head (i.e. water level in an aquifer) in the aquifer-system, which drives extraction-induced subsidence, is fully expressed (Galloway andBurbey 2011, Isotton et al 2015). This means that subsidence can still happen well after groundwater extraction has stopped and hydraulic heads are rising again. It is vital to gain insights in these longer-term processes and effects to assess the impacts of future groundwater extraction on subsidence of deltas and coastal systems facing high rates of relative sea-level rise as a result of groundwater overexploitation. In this paper, we aim to study long-term, hydrogeological and geotechnical behavior of the Mekong delta's aquifer system and to provide the first process-based quantification of possible future groundwater extraction on delta-wide subsidence to support informed decision-making in the delta.
Recent research resulted in the first delta-wide 3D hydrogeological model, coupled to a geotechnical module of the Mekong delta (Minderhoud et al 2017). This model enabled process-based modelling of groundwater flow and compaction of the delta's multi-aquifer system. To enable modeling of future groundwater extractioninduced subsidence, we advanced this model further and developed six mitigation and non-mitigation scenarios focused on hydraulic head development following different groundwater extraction pathways until 2100. These scenarios reveal the long-term physical behavior of the entire delta and provide valuable insights on recharge and recovery potential of the aquifer system. Furthermore, our results provide the first non-linear, process-based quantitative estimates of potential future extraction-induced subsidence for the Mekong delta. Their potential consequences for future elevation of the delta are presented in spatially explicit maps using a vertically high resolution elevation model of the Mekong delta (Minderhoud et al 2019b). The results quantify the considerable impact that mitigation efforts, focused on the reduction of groundwater extraction, may have on reducing future elevation loss. We found that the window for mitigation to keep the delta elevated above sea level is rapidly closing, and the effectiveness of policy implementation on groundwater-extraction will determine the future of this low-lying mega-delta. These new insights are relevant for decision makers in the Mekong delta and demonstrate the impact of groundwater extraction and mitigation potential on delta subsidence for deltas and coastal plains in other places on Earth.

Groundwater and subsidence model
We used an extended and updated version of the 3D hydrogeological groundwater model of the Mekong delta developed by Minderhoud et al (2017) to model scenarios of future groundwater extraction. Originally, this hydrogeological model was used to model groundwater flow for the period 1991-2015 using the MODFLOWbased environment iMOD (Vermeulen 2006, Vermeulen et al 2018. The 3D hydrogeological representation of the delta's multi-aquifer system was created using the iMOD SolidTool (Vermeulen et al 2018) by interpolating 95 borehole logs in ten cross-sections. The model distinguishes the seven main hydrogeological units determined by the Division of Geological Mapping for the South of Vietnam (DGMS 2004), each consisting of an aquifer and an overlying aquitard. Each aquifer and aquitard was modeled explicitly with a horizonal resolution of 1×1 km 2 . The phreatic top layer at the delta system was represented by an two meter thick model layer overlying the multi-aquifer system. The surface elevation of the delta was based on a digital elevation model derived from topographical elevation data for the Mekong delta (Minderhoud et al 2019b) and supplemented with Shuttle Radar Topography Mission (SRTM) elevation data for areas outside the delta (For additional details on the 3D hydrogeological model see Supporting Information (SI), figure S.1 and table S.1 is available online at stacks.iop.org/ERC/2/011005/mmedia). The hydrogeological model parameters were calibrated using measurements of hydraulic head throughout the delta. Variable density groundwater flow (e.g. the effect of saline water) is not included in the model. We extended the simulation period of the model to 1991-2100 to evaluate of future groundwater extraction scenarios and corresponding extraction-induced subsidence. Furthermore, we improved the model by explicitly including the surface water network of the Mekong delta.
Adding surface water network The surface water system of the Mekong delta consists of a dense network of natural river branches, canals and tidal creeks. In the previous version of the model the interaction with the surface water system was modeled through a constant drain at the delta surface, but the individual Mekong river channels and the extensive canal system were not explicitly included in the model. As a result, the modelled hydraulic heads in confined aquifers close to large rivers showed a mismatch with the observed hydraulic heads, as the modeled recharge at these locations was too low (also discussed in Minderhoud et al 2017). To use the model for future predictions and simulate groundwater dynamics over more than a century, accurate modeling of groundwater recharge is essential as small mismatches cumulate to large quantities over time. Therefore, we updated the model by explicitly including the surface water system of the Mekong delta as boundary condition to improve modeled recharge of river water to the aquifer system. We divided input data on the surface water system into four categories which are (1) main river channels, (2) secondary river channels, (3) main canals and (4) secondary canals. Average width and depth estimates for each category were used to determine river depth and to estimate bed conductance (figure S.2). River stage measurements from 1999 to 2010 were supplied by the Division of Water Resources Planning and Investigation for the South of Vietnam (DWRPIS) from 39 locations in the main and secondary rivers (Bui et al 2013). These measurements were interpolated to derive average annual river stage for the entire delta (figure S.2). Implementation of the surface water system into the model resulted in an improved modeling of surface water-groundwater interaction. Moreover, it increased the modeled recharge of confined aquifers in places where river channels cut through the Holocene aquitard, a phenomenon also observed in hydraulic head measurements. As a result, the overall mean correlation coefficient (r 2 ) between observed and modeled hydraulic heads in the Mekong delta for the period 2000-2015 improved from 0.69 to 0.73 by modeling the surface water system explicitly.

Subsidence calculations
We calculated subsidence as a result of aquifer-system compaction following decreases in hydraulic head (i.e. decreasing pressure) using a one-way coupled, geotechnical subsidence module called SUB-CR (i.e. the hydrogeological schematization does not change or compress during the modeling period). We applied the socalled abc model, in which a (recompression or swelling constant) accounts for the elastic compression, b (compression constant) and c (secondary compression constant) for the viscoplastic compression (Den Haan 1994). This model determines natural strain (i.e. degree of compression) based on the isotach concept (Šuklje 1957(Šuklje , Bjerrum 1967) as a function of effective stress and intrinsic time using the abc constants. The model only considers vertical deformation. The hydrological effect of viscous compression, which tends to increase pore pressure and therefore hydraulic head as water is squeezed from a compressing laying, was set to zero as the model was calibrated without this budget term. The parametrization of the abc constants was adopted from Minderhoud et al (2017) (table S.2). Apart from the abc constants, the modeling of secondary compression also depends on the overconsolidation ratio (OCR) which is described as follows: where σ′ p is the initial pre-consolidation stress and σ′ the momentary effective stress. A lower OCR value results in higher rates of secondary compression. For the Mekong delta only limited data is available to constrain OCRs. Hoang et al (2016) found an average OCR value of 1.6 for clayey deposits in the Mekong delta and Thoang and Giao (2015) reported a similar value for medium to stiff clays in HCMC province. Validation of modeled subsidence in the Mekong delta with InSAR-derived subsidence rates from 2006-2010 (Erban et al 2014) revealed an almost similar OCR value of 1.63 to provide the best fit (Minderhoud et al 2017). We applied the same validation after updating the model (section S1.1.1 and S1.2) and the OCR value that provided the best fit with the InSAR-derived subsidence rates remains 1.63. Therefore, we also use this OCR value to calculate the best estimate subsidence values for this study. To provide a range of subsidence calculations from least conservative (very weak sediments) to most conservative estimates (rigid sediment properties) we vary the OCR value by 0.1 (i.e. 1.53-1.73). We selected a slightly narrower OCR range than used by Minderhoud et al (2017), as it provides modeling results that were better supported by the InSAR-based subsidence rate observations providing assumably a more realistic range (section S.3.1). The presented results of average modeled cumulative subsidence and subsidence rates for the different extraction scenarios are based on the best estimate model (OCR=1.63) and specifically calculated for the Mekong delta (delineated in figure 1). Furthermore, as the extraction data does not cover the entire Mekong delta, results of average delta-wide subsidence presented in this paper were calculated based on areas in the Mekong delta within a 5 km radius of a modeled extraction (figure S.4).

Updating past groundwater extraction
In the previous model version the growth of unregistered extraction volume was estimated to be similar to the growth of the registered wells volume (Minderhoud et al 2017). The hydrogeological modeling results showed that throughout the delta the modeled hydraulic head decline during the start of the modeling period slightly underestimated measured hydraulic head declines and overestimated hydraulic head declines at the end of the original modeling period towards the present. For modelling future development of hydraulic head declines and subsidence, the overestimation of hydraulic head decline at the end of the modeling period (∼2010-2015) is problematic as extrapolation into the future increases such initially small overestimations further. To address this limit, we based the annual growth in unregistered extraction volume directly on measured hydraulic head decline in the confined aquifers. This alternative approach follows the assumption that an increase in measured hydraulic head drop in an aquifer is linked to an increase in extracted volume. It was applied for all unregistered extractions in the Mekong delta and Ho Chi Minh city (HCMC) province during the modeling period 1991-2011 (section S.5 for a detailed description of the approach and modelled annual groundwater extraction).
For the period of 2011-2018, an annual growth of 2.5% was simulated for the Mekong delta based on estimates of extraction increase by the DWRPIS. For HCMC province, extraction gradually stabilized following the water act in 2007 (HCMC 2007). The update of the annual growth of unregistered extractions resulted in an improved fit with the model and the observed hydraulic head decline (figure S.6) and the model no longer overestimates the decline towards the present. Compared to the previous model version, the overall mean correlation coefficient (r 2 ) between observed and modeled hydraulic heads in the Mekong delta for the period 2000-2015 improved from 0.69 to 0.75 with updating past extraction. Combined with implementing the surface water system in the model (section S1.1.1), the performance of the model increased further (r 2 =0.78). By improving the modeled extraction volume for the past, the modeled extracted volume provides a more realistic starting value for the modeling of future extraction pathways.

Groundwater extraction scenarios
The updated and extended model enables to evaluate the future hydraulic head evolution and consequent aquifer-system compaction of the Mekong delta. We developed six scenarios to simulate different possible pathways of future evolution of total extracted volume in the Mekong delta from 2019 to 2100 under various decrees of mitigation (table 1, figure 2). The pathways were developed either focused directly on the amount of groundwater extracted in the delta or, indirectly, on the effect of extraction, i.e. hydraulic head development over time. This second pathway type was developed because maintaining certain hydraulic heads can be implemented as mitigation measure. This is, for example, the case in HCMC where, following the water act in 2007, groundwater extraction is no longer allowed when hydraulic head in an aquifer falls below a certain predetermined level (HCMC 2007). Two non-mitigation scenarios follow pathways in which the amount of groundwater extraction continues to grow: scenario B1 represents a future in which groundwater extraction continues to increase moderately (annual increase of 2% of the 2018 volume; 55×10 3 m 3 daily extraction) and scenario B2 represents a worst case scenario in which extraction increases double this amount (annual increase of 4% of the 2018 volume; 110×10 3 m 3 daily extraction).
Four mitigation scenarios follow groundwater extraction pathways aiming to limit extraction growth and/ or reduce total extracted volume. Scenario M1, M2 and M3 have been developed in a way that they could represent realistic cases by incorporating a transition period during which the extracted volume is gradually stabilized or reduced. This period is needed to reduce groundwater use and to create the infrastructure needed to provide for alternative fresh water sources to meet the fresh water demand. Mitigation scenario M1 represents a Table 1. Modeled scenarios of groundwater extraction pathways from 2018 until 2100. The percentual annual change in extracted volume was applied to all wells included in the model.

Scenario
Extraction pathways stabilization of the extracted volume, allowing limited volume growth until 2020 to realistically incorporate the effects of new wells that are already licensed and currently being constructed. After 2020 there is no further increase and the total extracted volume (∼2.8 million m 3 daily) remains stable until 2100. Mitigation scenario M2 focuses on the stabilization of hydraulic heads and aims to maintain present hydraulic heads until 2100. In this scenario the modeled natural recharge of the delta system equals the amount of groundwater extracted from the subsurface and this requires a reduction of groundwater extraction of 50% of the 2018 volume (∼1.4 million m 3 daily). Scenario M3 investigates a situation in which the recovery of groundwater levels is the main focus but groundwater is still extracted in small quantities. In this scenarios groundwater extraction is reduced with 75% since 2018 (∼0.7 million m 3 daily), to allow recovery of the hydraulic head through natural recharge. Scenario M4 is a theoretical case to investigate the response of the Mekong delta system after an abrupt stop of all extraction after 2018. This scenario was developed to quantify the maximum recovery rate of the hydraulic heads in the aquifer system and the minimum amount of future subsidence that will inevitably happen as inheritance of three decades of hydraulic head lowering.
In each scenario the groundwater extraction pathway was described as a percentual change in total extraction volume compared to the total volume of groundwater extraction of 2018 (table 1). This annual change was simulated in all existing wells in the model located in both the Mekong delta and HCMC province. No new wells were added nor was there spatial variability in extraction growth included. For each scenario, the hydrogeological response of the groundwater system was modeled and the amount of extraction-induced subsidence quantified until 2100.

Projections of future elevation
We created spatially explicit projections of future delta elevation to local mean sea level by combining the estimates of potential extraction-induced subsidence for different scenarios with climate-change driven sealevel rise and a vertically high resolution Digital Elevation Model ( . As both processes have not yet been quantified for all locations in the Mekong delta, for this analysis we follow the assumption that they counterbalance each other in term of net elevation change over time, which also seems to have been the case in the past when the Holocene delta plain was formed (Zoccarato et al 2018).

Hydrogeological evolution following groundwater extraction pathways
We developed six scenarios to simulate different possible pathways of future evolution of total extracted volume in the Mekong delta from 2019 to 2100 under various decrees of mitigation (figure 2; table 1). Future average delta-wide hydraulic head evolution for each scenario is presented in figure 3. The average hydraulic head is an equally weighted average of the six confined aquifers (Upper Pleistocene to Upper Miocene) in the subsurface. For both non-mitigation scenarios (B1 and B2), in which groundwater extraction continues to increase, the average hydraulic head in the Mekong delta continues to drop throughout the modelling period to −30 m in scenario B1 and to almost −50 m in scenario B2 by 2100. When the extraction volume stabilizes after 2020 (scenario M1), the average hydraulic head continues to drop until it gradually stabilizes towards the end of the century at a level twice the present-day hydraulic head level below 0 meter (−13.5 m). Scenario M2 and M3 show both a short decreasing drop in average hydraulic head, after 2018 followed by a recovery of the heads in the decades afterwards. Scenario M2 represents the extraction pathway aiming to stabilize the average hydraulic head at its 2018 level of −6.5 m, which it does after a short dip and recovery. In the recovery scenario (M3), the hydraulic heads recover towards a level half the average hydraulic head in 2018 (−3.2 m in 2100, SI table S.3). The full extraction stop scenario (M4) shows a recovery of the head to levels slightly above mean sea level by the end of the modeling period, comparable to the hydrogeological situation in the delta at the end of the 20th century before the excessive exploitation of the groundwater system had started. The average hydraulic head recovered to 0 meter by the year 2077, comparable to the situation in the delta in 1997.
The above hydraulic head values are delta-wide, multi-aquifer averages. Hydraulic head evolution in each separate aquifer is dependent on its specific balance between extraction and recharge, which is also spatially variable. The maps of average hydraulic head evolution for the different extraction pathways (figure 4) show the spatial variability in the Mekong delta. In general, aquifer recharge is higher in the central and northern part of the Mekong delta, where the presence of large rivers in combination with a thinner Holocene aquitard (less lowpermeable clays at the delta surface) enable more recharge when compared to the southern part of the delta. This becomes especially visible in scenario M2 as hydraulic heads in the central and northern part of the delta slightly  increase towards the end of the century, while in the Ca Mau peninsula in the south, the heads slightly decrease compared to the 2018 levels.

Modeled extraction-induced subsidence
Future extraction-induced subsidence The evolution of subsidence rates and cumulative subsidence since 2018 is modeled until the end of the 21th century following the six groundwater extraction pathways. Figure 5 shows the spatial subsidence patterns of the best estimate model projections for the Mekong delta and HCMC province. The development of the average Mekong delta-wide subsidence rates and cumulative subsidence for each scenario during the period 2000-2100 are shown respectively, in figures 6 and 7.
In the two non-mitigation scenarios the total volume of groundwater extraction continues to increase and as a result subsidence rates remain high throughout the 21st century. In scenario B2, the delta-wide average subsidence rate keeps increasing until 2078 at 13.1 mm yr −1 (12.1-13.6 mm yr −1 for the most conservative and least conservative models) (figure 6). At the end of the century, the average delta-wide subsidence rate is 12.9 (12.3-13.1) mm yr −1 with certain areas in the delta experiencing rates up to 45 mm yr −1 (see SI table S.5 for highest modeled subsidence rates). Cumulatively, the delta subsides on average 100 (86-110) cm during the period 2018-2100 in scenario B2 (figure 7). In scenario B1, the average delta-wide subsidence rate increases from 2018 towards 2023 to 8.9 (5.84-11.9) mm yr −1 followed by a gradual decline in average rate towards 7.6 (7.0-7.9) mm yr −1 in 2100 with highest rates going up to 20 mm yr −1 . In scenario B1, during the 21st century following 2018, the delta subsides on average a total of 68 (55-70) cm.
In all four mitigation scenarios in which groundwater extraction is kept stable or decreases in the future, the subsidence rates also decrease. The largest differences between the mitigation scenarios in term of subsidence rates can be observed in the coming first decades (figure 6). The largest abrupt change in subsidence rate happens in scenario M4, where subsidence rates drop sharply from 8.9 (5.5-12.5) to 1.6 (−1.3-4.9) mm yr −1 as a result of the abrupt stop of all groundwater extraction. In 2050 average modeled subsidence rates of the Mekong delta for scenario M1 to M4 are respectively, 5.0 (3.5-5.9) mm yr −1 , 2.9 (1.7-3.7) mm yr −1 , 1.9 (1.0-2.7) mm yr −1 , 1.6 (0.7-2.3) mm yr −1 , while maximum modeled rates for each of the respected scenarios are 11.1 (10.2-11.7) mm yr −1 , 4.6 (3.6-6.1) mm yr −1 , Figure 5. The evolution of subsidence rate (left) and cumulative subsidence (right) in the Mekong delta until 2100 following different groundwater extraction pathways. Cumulative subsidence is calculated from 2018 onwards. The quantifications are based on the best estimate model.
3.4 (1.7-4.8) mm yr −1 , 2.9 (1.4-4.3) mm yr −1 (SI table S.5). Towards the end of the century, the subsidence rates all gradually converge towards each other and rates in 2100 range from 1.2-2.5 (0.7-2.7) mm yr −1 for all mitigation scenarios, with maximum rates ranging from 2.2-4.4 (1.3-4.6) mm yr −1 . Cumulatively, the mitigation scenarios result in relatively large differences. By 2100, the Mekong delta experiences, in respectively scenario M1 to M4, on average a cumulative subsidence of 39 (28-48) cm, 25 (15-33) cm, 19 (10-27) and 12 (4-18) cm since 2018. Spatially, subsidence rates vary as a result of location and depth of extraction wells and extracted volume in relation to the subsurface architecture and composition, which determines local recharge rates and compaction potential. Higher subsidence rates are found in areas with larger extraction volumes and low recharge rates, thus experiencing larger and longer sustained drawdowns in hydraulic head ( figure 5). In general, rates show an increasing trend from the apex  of the delta in the northwest towards the coastline in the southeast as the aquifer system thickens in this direction, which increases its total compressibility potential. Figures 6 and 7 also show the range in projections of absolute sea-level rise, caused by accelerated thermal expansion of seawater and melting of ice sheets as a result of global warming, using the RCP 2.6 and 8.5 climate change scenarios (Church et al 2013). At present, modeled extraction-induced subsidence rates (delta-wide average: 9 (5-13) mm yr −1 , locally up to 39 (33-42) mm yr −1 ) are much larger than rates of present absolute sealevel rise (3-4 mm yr −1 ). During this century, when groundwater extraction continues to increase, scenario B1 and B2, the extraction-induced subsidence rates will remain higher or, at some point, equal to rates of absolute sea-level rise. In case of mitigation (scenario M1-M3) rates will gradually decrease and at some point become less than climate-change driven sea-level rise.

Discussion
Limitations on modeling the hydrology and subsidence of the Mekong delta aquifer system Future change in total groundwater extraction volume was modeled based on the delta-wide extraction database of the DWRPIS which is, at present, the best available data on groundwater extraction available for the delta. Uncertainties in extracted volume and missing extractions may affect the modeled spatial patterns of hydraulic head in the aquifers (discussed in Minderhoud et al 2017). We modeled the changes in extracted volume of each extraction pathway uniformly over the delta with the aim to investigate the potential future behavior of the aquifer system. In reality, we expect changes in future extraction to be much more location-specific and temporally variable, as they are influenced by a complex interaction between groundwater governance, legislation and law enforcement which can vary locally in the delta (Ha et al 2018). Furthermore, also local differences in access to good quality surface water, socio-economical situation, technical capability and land-use practices, which, in turn, may be influenced by physical processes and feedback-mechanisms, e.g. ongoing subsidence increasing salinization, will determine future groundwater use. Our modelling approach paves the road to investigate more complex and realistic scenarios that include abovementioned factors.
Modeled extraction-induced subsidence is affected by uncertainties in the hydrogeological model and geotechnical parameterization (discussed in Minderhoud et al 2017). The difference between modeled subsidence rates for the most and least conservative geotechnical parameterization becomes smaller towards the end of the modeling period (figure 6), which shows that the relative influence of initial overconsolidation ratio uncertainty decreases in time (Suklje 1957). While geotechnical model uncertainty decreases, the influence on modeled subsidence of the groundwater extraction pathways increases. Even though locally uncertainties in modeled extraction-induced subsidence may be considerable, on the scale of the entire delta our results do provide a first indicative, process-based quantification on how the aquifer-system may respond to different groundwater extraction pathways in the 21st century.
Response of the aquifer system to future groundwater extraction The evolution of the hydrogeological situation in the aquifer system of the Mekong delta will be determined by the extraction pathway followed. When the amount of groundwater extraction in the Mekong delta remains stable after 2020 (scenario M1) or continues to increase further into the future (scenario B1 and B2), aquifer system depletion will continue as well, resulting in continuous hydraulic head drop. In scenario B2 the average hydraulic head in the delta will drop almost to −50 m by the end of the century. Although we consider this as a most extreme scenario, drawdowns of this magnitude in confined aquifers are not uncommon and have been reported throughout the world, a.o. Bangkok . Moreover, also nearby Ho Chi Minh city drawdowns more than 30 meter have been measured between 1994 and 2015. Although drawdowns are site-specific and depend on many factors like aquifers and aquitards size and storage capacity, groundwater extraction and recharge rate, physically the hydrogeological situation modeled in scenario B2 is possible.
The modeling results of the mitigation scenarios show the potential of the Mekong delta's aquifer system to recover from past and future hydraulic head drops when extraction is reduced or stopped. They also reveal the amount of present overexploitation. Currently, the amount of groundwater extraction in the delta is about twice the amount of recharge. A reduction of extraction by 50% will result in a stabilization of average hydraulic head at its present level (scenario M2), indicating that the extracted volume equals the amount of aquifer recharge by infiltration of precipitation, surface water and intrusion of sea water. Maximum recovery of the hydraulic heads in the Mekong delta's aquifer system is achieved when extraction is completely stopped (scenario M4). Under these optimal circumstances, it will take the aquifer system of the Mekong delta ∼60 years to recover from the hydraulic head drop caused by groundwater extraction during the past ∼20 years. This demonstrates the low recharge rate of the delta's aquifer system, which is likely caused by the low permeable Holocene deposits at the delta surface (Minderhoud et al 2017). This means that even after a complete extraction stop, past groundwater overexploitation will still have a continued effect on the hydrogeological situation of the aquifer system for future decades.
Past hydraulic head declines in the aquifer system have triggered aquifer-system compaction which led to subsidence of the delta surface (Minderhoud et al 2017). This process of aquifer-system compaction is sluggish as it takes time for the decrease in pore pressure to propagate into the low-permeable and compressible aquitards, causing a delay in primary compression. Secondary compression, which is time-dependent, creates additional delayed compression of the aquifer-system which can continue for decades after it is triggered. In combination with the slow recharge rate of the delta's aquifer system, this means that even when extraction is completely stopped (scenario M4), the delta will continue to subside, amounting to an average cumulative extraction-induced subsidence of 12 (4-18) cm by the end of the 21st century. This amount of subsidence is the inevitable inheritance of past groundwater overexploitation. In reality groundwater extraction will not stop after 2018, therefore the true amount of future extraction-induced subsidence will be higher. How much higher will be determined by the extraction pathway.
Future of the Mekong delta Consequences for delta elevation and relative sea-level rise Extraction-induced subsidence results in elevation loss at the delta surface. The magnitude of the impact of elevation loss on a delta's future sustainability and the livelihood of its inhabitants depends, to a large extent, by its elevation above sea level. Figure 8 shows the area in percentage of the Mekong delta plain that will experience a certain amount of cumulative subsidence for each extraction pathway in 2050 and 2100 and the average delta plain elevation of the Mekong delta (∼0.8 m a.m.s.l., Minderhoud et al 2019aMinderhoud et al , 2019b. It becomes evident that if groundwater extraction continues to increase in the future (scenario B1 and B2), extraction-induced subsidence alone may cause large parts of the delta to lose all elevation a.m.s.l. before the end of the century. In case of groundwater extraction mitigation, cumulative extraction-induced subsidence will cause less but still considerable amounts of elevation loss. The minimum elevation loss that the delta will experience by the end of the century, as a result of inevitable aquifer-system compaction (scenario M4), equals ∼15% (5%-22%) of the present average Mekong delta plain elevation a.m.s.l. Beside extraction-induced subsidence, the delta will also lose elevation a.m.s.l. as the sea level itself is rising. When we correct future delta elevation a.m.s.l. with projections of absolute sea-level rise (Church et al 2013), the percentual area that loses on average all elevation above sea level by extraction-induced subsidence considerably larger (supplementary figure S.7). Only when future extraction-induced subsidence is mitigated by strongly reducing groundwater extraction and rates of absolute sea-level rise follow a moderate pathway, the majority of the delta may still be elevated above sea level by the end of the century.
By combining the results of three plausible groundwater extraction pathways (B1, M1 and M4) with a moderate estimate of absolute sea-level rise and a high resolution vertical elevation model of the Mekong delta (Minderhoud et al 2019), the effects of different extraction scenarios are made spatially explicitly in maps projecting future elevation to mean sea level (figure 8). The maps reveal that considerable parts of the Mekong delta will likely fall below sea level during the coming century, even under the M3 scenario. In particular the lowlying SW part of the delta will fall below sea level already in the course of this century, regardless the considered scenario. The extent of those parts of the delta that remain above sea level by the end of the century varies from 32% in B1 to 44% in M1 and 57% in M3, showing the large potential of mitigation to reduced extractioninduced subsidence.
Beside extraction-induced subsidence and absolute sea-level rise, relative sea-level rise in the delta is also determined by other subsidence processes, like tectonics, natural compaction of Holocene sediments, which can amount up to rates of several cm yr −1 in the Mekong delta (Zoccarato et al 2018), drainage of surface water and loading by buildings and infrastructure . As figure 9 shows, whether and when a certain part of the delta will fall below sea level is location-specific and dependent on local conditions of relative sea-level rise and surface elevation. On top of that, sedimentation of clastic and organic sediments on the delta surface is the natural mechanism of a delta to keep up with relative sea-level rise and can even result in elevation gain. This is the case in some areas along the Mekong delta's coastline, hosting pristine mangrove forests with abundant sediment availability (Lovelock et  Time for mitigation, mitigation to create time Given its lowly elevated topography, ongoing subsidence, absolute sea-level rise and decreased sedimentation, the Mekong delta is heading for a tipping point towards a 'collapsed' state (Renaud et al 2013). As the delta is progressively losing elevation, the window for mitigation is closing rapidly. Although relative sea-level rise for the Mekong delta is caused by more than groundwater extraction-induced subsidence, it is, at present, the dominant factor. When and whether mitigation strategies to reduce extraction-induced subsidence are successfully implemented, will determine, for a large extent, how long the delta will stay above sea level. Although part of the future extraction-induced subsidence in inevitable (Scenario M4), a large part still can still be mitigated and avoided. A strong reduction in groundwater extraction would not only allow the aquifer system to recharge, and subsequently decrease aquifer system compaction, it would also reduce other water-related issues, such as salt water intrusion in the aquifer system and decrease in water quality (Renaud et al 2015, Smajgl et al 2015, Eslami et al 2019. Although the best solution to reduce extraction-induced subsidence is to immediately stop all groundwater extraction, this is realistically not an option, as people in the delta rely on groundwater for their freshwater supply. Until alternative water sources are available, such as a piped water supply, high quality surface water or desalinization, groundwater will continue to be used to meet the fresh water demand. While investments are being made and infrastructure is being developed to provide an alternative fresh water supply, implementation of smart extraction strategies can already reduce extraction-induced subsidence while supplying fresh groundwater. For example, extraction can be relocated to areas that are less exploited, contain less compression-prone, fine-grained sediments (i.e. less clay, more sand) and experience natural recharge rates (e.g. close to rivers). Additionally, extraction could be concentrated in higher elevated parts of the Figure 8. Projections of future elevation with extraction-induced subsidence following three different extraction scenarios and absolute sea-level rise according to the mid-range projections of SLR (median) of RCP 4.5. Elevation is in meter above mean sea-level based on the 'Topo DEM' (Minderhoud et al 2019a(Minderhoud et al , 2019b. Additional elevation loss by subsidence as a result of natural compaction is assumed to be counterbalanced by elevation gain through deposition of new sediments. Blueish areas are below mean sea level. delta, where subsidence is less harmful. However, such strategies only provide temporal solutions and cannot substitute the reduction in overall groundwater use. Managed aquifer recharge (MAR) provides opportunities to artificially stimulate the aquifer system recharge, by for example using groundwater injection wells. This may also help to battle salt water intrusion and simultaneously creating strategic freshwater reserves in the subsurface. All the above mentioned solutions will require fundamental changes in groundwater management and law enforcement, availability of alternative water sources, changes in agricultural practices, and investments in infrastructure to distribute fresh water over the delta. In November 2017, the Vietnamese Government issued 'resolution 120' which describes the ambition for a prosperous, sustainable and climate-resilient future of the Mekong delta, including the aspiration to end groundwater use by the year 2100. In December 2018, the Vietnamese Government issued degree 167 (167/2018/ND-CP) to restrict extraction in areas where it is causing subsidence, pollution, salinization and depletion of groundwater resources and delineate groundwater protection zones. While these intentions provide an optimistic outlook for the delta, our results reveal that the available time for implementing the required changes to reduce groundwater extraction and effectively mitigate subsidence is very limited and shorter than currently foreseen in governmental plans. The race against the clock has started and any delay in implementation will cause the delta to lose more of its elevation above sea level, which will increase the exposure of the nearly 18 million inhabitants to flooding, storm surges and groundwater salinization. The food producing capacity will decrease and the costs of flood protection to prevent permanent submersion of the delta will increase.