Articles | Volume 12, issue 3
Earth Syst. Dynam., 12, 919–938, 2021

Special issue: Understanding compound weather and climate events and related...

Earth Syst. Dynam., 12, 919–938, 2021

Research article 13 Sep 2021

Research article | 13 Sep 2021

Exploring how groundwater buffers the influence of heatwaves on vegetation function during multi-year droughts

Exploring how groundwater buffers the influence of heatwaves on vegetation function during multi-year droughts
Mengyuan Mu1, Martin G. De Kauwe1,2, Anna M. Ukkola1, Andy J. Pitman1, Weidong Guo3, Sanaa Hobeichi1, and Peter R. Briggs4 Mengyuan Mu et al.
  • 1ARC Centre of Excellence for Climate Extremes and Climate Change Research Centre, University of New South Wales, Sydney 2052, Australia
  • 2School of Biological Sciences, University of Bristol, Bristol, BS8 1TQ, United Kingdom
  • 3School of Atmospheric Sciences and Joint International Research Laboratory of Atmospheric and Earth System Sciences, Nanjing University, Nanjing 210023, China
  • 4Climate Science Centre, CSIRO Oceans and Atmosphere, Canberra 2601, ACT, Australia

Correspondence: Mengyuan Mu (


The co-occurrence of droughts and heatwaves can have significant impacts on many socioeconomic and environmental systems. Groundwater has the potential to moderate the impact of droughts and heatwaves by moistening the soil and enabling vegetation to maintain higher evaporation, thereby cooling the canopy. We use the Community Atmosphere Biosphere Land Exchange (CABLE) land surface model, coupled to a groundwater scheme, to examine how groundwater influences ecosystems under conditions of co-occurring droughts and heatwaves. We focus specifically on south-east Australia for the period 2000–2019, when two significant droughts and multiple extreme heatwave events occurred. We found groundwater plays an important role in helping vegetation maintain transpiration, particularly in the first 1–2 years of a multi-year drought. Groundwater impedes gravity-driven drainage and moistens the root zone via capillary rise. These mechanisms reduced forest canopy temperatures by up to 5 C during individual heatwaves, particularly where the water table depth is shallow. The role of groundwater diminishes as the drought lengthens beyond 2 years and soil water reserves are depleted. Further, the lack of deep roots or stomatal closure caused by high vapour pressure deficit or high temperatures can reduce the additional transpiration induced by groundwater. The capacity of groundwater to moderate both water and heat stress on ecosystems during simultaneous droughts and heatwaves is not represented in most global climate models, suggesting that model projections may overestimate the risk of these events in the future.

1 Introduction

Droughts and heatwaves are important socioeconomic and environmental phenomena, impacting regional food production (Kim et al., 2019; Lesk et al., 2016), water resources (Leblanc et al., 2009; Orth and Destouni, 2018), and the resilience of ecosystems (Ibáñez et al., 2019; Ruehr et al., 2019; Sandi et al., 2020). When droughts and heatwaves co-occur (a “compound event”) the consequences can be particularly severe, reducing the terrestrial carbon sink (Ciais et al., 2005), potentially accelerating tree die-off (Allen et al., 2010, 2015; Birami et al., 2018), and setting conditions conducive to wildfires (Jyoteeshkumar reddy et al., 2021). One region experiencing severe coincident heatwaves and drought is Australia (Mitchell et al., 2014). Drought in Australia is associated with large-scale modes of variability, including the El Niño–Southern Oscillation and the Indian Ocean Dipole (van Dijk et al., 2013), and periods of below-average rainfall can extend for multiple years (Verdon-Kidd and Kiem, 2009). Heatwaves are commonly synoptically driven, associated with blocking events that can be sustained over many days (Perkins-Kirkpatrick et al., 2016; Perkins, 2015). Modes of variability and synoptic situations are important in setting up conditions conducive to drought and heatwave. However, once a heatwave or drought has become established, land–atmosphere interactions can intensify and prolong both heatwaves and droughts (Miralles et al., 2019), affect their intensity, and influence the risk of their co-occurrence (Mukherjee et al., 2020). The role of the land surface in amplifying or dampening heatwaves and droughts is associated with the partitioning of available energy between sensible and latent heat (Fischer et al., 2007; Hirsch et al., 2019) and is regulated by subsurface water availability (Teuling et al., 2013; Zhou et al., 2019). As soil moisture becomes more limiting, more of the available energy is converted into sensible heat, reducing evaporative cooling via latent heat. Changes in the surface turbulent energy fluxes influence the humidity in the boundary layer, the formation of clouds, incoming solar radiation, and the generation of rainfall (D'Odorico and Porporato, 2004; Seneviratne et al., 2010; Zhou et al., 2019). The sensible heat fluxes warm the boundary layer, leading to heat that can accumulate over several days and exacerbate heat extremes (Miralles et al., 2014), which can in turn increase the atmospheric demand for water and intensify drought (Miralles et al., 2019; Schumacher et al., 2019).

Vegetation access to groundwater has the potential to alter these land–atmosphere feedbacks by maintaining vegetation function during extended dry periods, supporting transpiration and moderating the impact of droughts and heatwaves (Marchionni et al., 2020; Miller et al., 2010). Where the water table is relatively shallow, capillarity may bring water from the groundwater towards the surface root zone, increasing plant water availability. Where the water table is deeper, phreatophytic vegetation with tap roots can directly access groundwater (Zencich et al., 2002). The presence of groundwater and the access to groundwater by vegetation are therefore likely to buffer vegetation drought and heatwave stress. For example, groundwater may help vegetation sustain transpiration and consequently cool plant canopies via evaporation. This is particularly critical during compound events where cessation of transpiration would increase the risk of impaired physiological function and the likelihood that plants would exceed thermal limits and risk mortality (Geange et al., 2021; O'sullivan et al., 2017; Sandi et al., 2020).

Quantifying the influence of groundwater on vegetation function has remained challenging as concurrent observations of groundwater dynamics, soil moisture, and energy and water fluxes are generally lacking over most of Australia and indeed many parts of the world. Land surface models (LSMs) provide an alternative tool for studying the interactions between groundwater, vegetation, and surface fluxes in the context of heatwaves and droughts (Gilbert et al., 2017; Martinez et al., 2016a; Maxwell et al., 2011; Shrestha et al., 2014). However, there has been very little work focused on the influence of groundwater on droughts and heatwaves occurring at the same time (Keune et al., 2016; Zipper et al., 2019). Our key goal in this paper is therefore to examine the timescales and extent to which vegetation utilises groundwater during drought and heatwaves and determine the degree to which groundwater can mitigate the impacts of compound extremes. We focus on droughts and heatwaves occurring over south-eastern (SE) Australia during 2000–2019 using the Community Atmosphere Biosphere Land Exchange (CABLE) LSM. SE Australia is an ideal case study since its forest and woodland ecosystems are known to be dependent on groundwater (Eamus and Froend, 2006; Kuginis et al., 2016; Zencich et al., 2002), and it has experienced two multi-year droughts and record-breaking heatwaves over the last 2 decades. By examining the role of groundwater in influencing droughts and heatwaves and by understanding how well CABLE can capture the relevant processes, we aim to build confidence in the simulations of land–atmosphere interactions for future droughts and heatwaves.

2 Methods

2.1 Study area

The climate over SE Australia varies from humid temperate near the coast to semi-arid in the interior. In the last 20 years, SE Australia experienced the 9-year Millennium Drought during 2001–2009 (van Dijk et al., 2013), where rainfall dropped from a climatological average (1970–1999) of 542 to 449 mm yr−1, and a 3-year intense recent drought during 2017–2019 where rainfall dropped to 354 mm yr−1 (Fig. S1). It has also suffered record-breaking summer heatwaves in 2009, 2013, 2017, and 2019 (Bureau of Meteorology, 2013, 2017, 2019; National Climate Centre, 2009). Here we investigate groundwater interactions during the period 2000–2019, focusing on the Millennium Drought (MD, 2001–2009) and the recent drought (RD, 2017–2019).

2.2 Overview of CABLE

CABLE is a process-based LSM that simulates the interactions between climate, plant physiology, and hydrology (Wang et al., 2011). Above ground, CABLE simulates the exchange of carbon, energy, and water fluxes, using a single-layer, two-leaf (sunlit and shaded) canopy model (Wang and Leuning, 1998), with a treatment of within-canopy turbulence (Raupach, 1994; Raupach et al., 1997). CABLE includes a six-layer soil model (down to 4.6 m) with soil hydraulic and thermal characteristics dependent on the soil type and soil moisture content. CABLE has been extensively evaluated (e.g. Abramowitz et al., 2008; Wang et al., 2011; Zhang et al., 2013) and benchmarked (Abramowitz, 2012; Best et al., 2015) at global and regional scales.

Here we adopt a version of CABLE (Decker, 2015; Decker et al., 2017) which includes a dynamic groundwater component with aquifer water storage. This version, CABLE-GW, has been previously evaluated by Decker (2015), Ukkola et al. (2016b), and Mu et al. (2021) and shown to perform well for simulating water fluxes. CABLE code is freely available upon registration (, last access: 29 August 2021); here we use CABLE SVN revision 7765.

2.3 Hydrology in CABLE-GW

The hydrology scheme in CABLE-GW solves the vertical redistribution of soil water via a modified Richards equation (Zeng and Decker, 2009):

(1) θ t = - z K z Ψ - Ψ E - F soil ,

where θ is the volumetric water content of the soil (mm3 mm−3), K is the hydraulic conductivity (mm s−1), z is the soil depth (mm), Ψ and ΨE are the soil matric potential (mm) and the equilibrium soil matric potential (mm), and Fsoil is a sink term related to subsurface runoff and transpiration (s−1) (Zeng and Decker, 2009; Decker, 2015). To simulate groundwater dynamics, an unconfined aquifer is added to the bottom of the soil column with a simple water balance model:

(2) d W aq d t = q re - q aq , sub ,

where Waq is the mass of water in the aquifer (mm), qaq,sub is the subsurface runoff in the aquifer (mm s−1), and qre is the water flux between the aquifer and the bottom soil layer (mm s−1) computed by the modified Darcy's law:

(3) q re = ( K aq + K bot ) 2 Ψ aq - Ψ n - Ψ E , aq - Ψ E , n z wtd - z n ,

where Kaq and Kbot are the hydraulic conductivity in the aquifer and in the bottom soil layer (mm s−1), Ψaq and ΨE,aq are the soil matric potentials for the aquifer (mm), and Ψn and ΨE,n are the soil matric potentials for the bottom soil layer (mm). zwtd and zn are the depth of the water table (mm) and the lowest soil layer (mm), respectively. The positive qre refers to the downward water flow from soil column to aquifer (i.e. vertical drainage, Dr), and the negative qre is the upward water movement from the aquifer to the soil column (i.e. recharge, Qrec). CABLE-GW assumes that the groundwater aquifer sits above impermeable bedrock, giving a bottom boundary condition of

(4) q out = 0 .

CABLE-GW computes the subsurface runoff (qsub, mm s−1) using

(5) q sub = sin d z d l q ^ sub e - z wtd f p ,

where dzdl is the mean subgrid-scale slope, q^sub is the maximum rate of subsurface drainage (mm s−1), and fp is a tunable parameter. qsub is generated from the aquifer and the saturated deep soil layers (below the third soil layer).

2.4 Experiment design

To explore how groundwater influences droughts and heatwaves, we designed two experiments, with and without groundwater dynamics, driven by the same 3 h time-evolving meteorology forcing and the same land surface properties (see Sect. 2.5 for datasets) for the period 1970–2019 at a 0.05 spatial resolution with a 3 h time step. To correct a tendency for high soil evaporation, we implemented a parameterisation of soil evaporation resistance that has previously been shown to improve the model (Decker et al., 2017; Mu et al., 2021).

2.4.1 Groundwater experiment (GW)

This simulation uses the default CABLE-GW model, which includes the unconfined aquifer to hold the groundwater storage and simulates the water flux between the bottom soil layer and the aquifer. We first ran the default CABLE-GW with fixed CO2 concentrations at 1969 levels for 90 years by looping the time-evolving meteorology forcing over 1970–1999. At the end of the 90-year spin-up, moisture in both the soil column and the groundwater aquifer reached an effective equilibrium when averaged over the study area. We then ran the model from 1970 to 2019 with time varying CO2. We omit the first 30 years of this period and analyse the period 2000–2019 to allow for further equilibrium with the time-evolving CO2.

2.4.2 Free drainage experiment (FD)

Many LSMs, including those used in the Coupled Model Intercomparison Project 5 (CMIP5), still use a free drainage assumption and neglect the parameterisation of the unconfined aquifer. To test the impact of this assumption we decoupled the aquifer from the bottom soil layer and thus removed the influence of groundwater dynamics (experiment FD). In FD, at the interface between the bottom soil layer and the aquifer, soil water can only move downwards as vertical drainage at the rate defined by the bottom soil layer's hydraulic conductivity:

(6) q re = K bot .

This vertical drainage is added to the subsurface runoff flux:

(7) q sub = q sub + q re .

The simulated water table depth (WTD, i.e. zwtd) in CABLE-GW affects the water potential gradient between the soil layers via ΨE (Zeng and Decker, 2009) and impacts qsub (Eq. 5). However, in FD, decoupling the soil column from the aquifer and adding vertical drainage directly to subsurface runoff causes an artificial and unrealistic decline in WTD. To solve this problem, we assume a fixed WTD in the FD simulations at 10 m in order to remove this artefact from the simulation of ΨE and qsub.

The FD simulations are initialised from the near-equilibrated state at the end of the 90-year spin-up used in GW. The period 1970–2019 is then simulated using varying CO2 and the last 20 years are used for analysis.

2.4.3 Deep root experiment (DR)

The parameterisation of roots, including the prescription of root parameters in LSMs, is very uncertain (Arora and Boer, 2003; Drewniak, 2019) and LSMs commonly employ root distributions that are too shallow (Wang and Dickinson, 2012). The vertical distribution of roots influences the degree to which plants can utilise groundwater, and potentially the role groundwater plays in influencing droughts and heatwaves. To explore the uncertainty associated with root distribution, we added a “deep root” (DR) experiment by increasing the effective rooting depth in CABLE for tree areas. In common with many LSMs, CABLE-GW defines the root distribution following Gale and Grigal (1987):

(8) f root = 1 - β root z ,

where froot is the cumulative root fraction (between 0 and 1) from the soil surface to depth z (cm), and βroot is a fitted parameter specified for each plant functional type (PFT) (Jackson et al., 1996). In CABLE, the tree areas in our study region are simulated as evergreen broadleaf PFT (Fig. S2a) with a βroot= 0.962, implying that only 8 % of the simulated roots are located below the soil depth of 0.64 m (Fig. S2b). However, field observations (Canadell et al., 1996; Eberbach and Burrows, 2006; Fan et al., 2017; Griffith et al., 2008) suggest that the local trees tend to have a far deeper root system, possibly to help cope with the high climate variability. We therefore increased βroot for the evergreen broadleaf PFT to 0.99, which assumes 56 % of roots are located at depths below 0.64 m and 21 % of roots below 1.7 m (Fig. S2b). This enables the roots to extract larger quantities of deep soil water moisture, which is more strongly influenced by groundwater.

Given that we lack the detailed observations to set root distributions across SE Australia, we undertake the DR experiment as a simple sensitivity study. We only run this experiment during January 2019, when the record-breaking heatwaves compound with the severe recent drought. The DR experiment uses identical meteorology forcing and land surface properties as GW and FD and is initialised by the state of the land surface on the 31 December 2018 from the GW experiment.

2.5 Datasets

Our simulations are driven by the atmospheric forcing from the Australian Water Availability Project (AWAP), which provides daily gridded data covering Australia at 0.05 spatial resolution (Jones et al., 2009). This dataset has been widely used to force LSMs for analysing the water and carbon balances in Australia (Haverd et al., 2013; De Kauwe et al., 2020; Raupach et al., 2013; Trudinger et al., 2016). The AWAP forcing data include observed fields of precipitation, solar radiation, minimum and maximum daily temperatures, and vapour pressure at 09:00 and 15:00. Since AWAP forcing does not include wind and air pressure we adopted the near-surface wind speed data from McVicar et al. (2008) and assume a fixed air pressure of 1000 hPa. Due to missing observations before 1990, the solar radiation input for 1970–1989 was built from the 1990–1999 daily climatology. Similarly, wind speeds for 1970–1974 are built from the 30-year climatology from 1975 to 2004. We translated the daily data into 3-hourly resolution using a weather generator (Haverd et al., 2013).

The land surface properties for our simulations are prescribed based on observational datasets. Land cover type is derived from the National Dynamic Land Cover Data of Australia (DLCD) (, last access: 29 August 2021). We classify DLCD's land cover types to five CABLE PFTs: crop (irrigated/rainfed crop, pasture, and sugar DLCD classes), broadleaf evergreen forest (closed/open/scattered/sparse tree), shrub (closed/open/scattered/sparse shrubs and open/scattered/sparse chenopod shrubs), grassland (open/sparse herbaceous), and barren land (bare areas). We then resample the DLCD dataset from the 250 m resolution to the 0.05 resolution. The leaf area index (LAI) in CABLE is prescribed using a monthly climatology derived from the Copernicus Global Land Service product (, last access: 29 August 2021). The climatology was constructed by first creating a monthly time series by taking the maximum of the 10 daily time steps each month and then calculating a climatology from the monthly data over the period 1999–2017. The LAI data were resampled from the original 1 km resolution to the 0.05 resolution following De Kauwe et al. (2020). Soil parameters are derived from the soil texture information (sand, clay, and silt fraction) from SoilGrids (Hengl et al., 2017) via the pedotransfer functions in Cosby et al. (1984) and resampled from 250 m to 0.05 resolution.

To evaluate the model simulations, we use monthly total water storage anomaly (TWSA) at 0.5 spatial resolution from the Gravity Recovery and Climate Experiment (GRACE) and GRACE Follow On products (Landerer et al., 2020; Watkins et al., 2015; Wiese et al., 2016, 2018). The RLM06M release is used for February 2002–June 2017 and for June 2018–December 2019. We also use total land evapotranspiration from the 2000–2018 monthly Derived Optimal Linear Combination Evapotranspiration (DOLCE version 2, Hobeichi et al., 2021a) at 0.25 resolution, as well as the 2000–2019 daily Global Land Evaporation Amsterdam Model (GLEAM version 3.5,; Martens et al., 2017; Miralles et al., 2011) at 0.5 spatial resolution. For daytime land surface temperature (LST) we use the Moderate Resolution Imaging Spectroradiometer (MODIS) datasets from Terra and Aqua satellites (products MOD11A1 and MYD11A1, Wan and Li, 1997; Wan et al., 2015a, b) at 1 km spatial resolution. We only consider pixels and time steps identified as good quality (QC flags 0). Only the daytime LST values are used due to the lack of good-quality nighttime LST data. The Terra overpass occurs at 10:00 and Aqua at 14:00 local time. To analyse the compound events in January 2019, we linearly interpolate the 3-hourly model outputs to 14:00 to match the overpass time of the Aqua LST. The GRACE, GLEAM, and MODIS datasets were resampled to the AWAP resolution using bilinear interpolation.

To evaluate model performance during heatwaves, we identify heatwave events using the excess heat factor index (EHF, Nairn and Fawcett, 2014). EHF is calculated using the daily AWAP maximum temperature, as the product of the difference of the previous 3 d mean to the 90th percentile of the 1970–1999 climatology and the difference of the previous 3 d mean to the preceding 30 d mean. A heatwave occurs when the EHF index is greater than 0 for at least 3 consecutive days. We only focus on summer heatwaves occurring between December and February of the following year.

3 Results

3.1 Simulations for the Millennium Drought and the recent drought

Previous studies have shown that simulations by LSMs diverge as the soil dries (Ukkola et al., 2016a), associated with systematic biases in evaporative fluxes and soil moisture states in the models (Mu et al., 2021; Swenson and Lawrence, 2014; Trugman et al., 2018). We therefore first evaluate how well CABLE-GW captures the evolution of terrestrial water variability during two recent major droughts.

Figure 1a shows the total water storage anomaly during 2000–2019 observed by GRACE and simulated in GW and FD. Both GW and FD accurately capture the interannual variability in total water storage for SE Australia (r= 0.96 in GW and 0.90 in FD). Both model configurations simulate a decline in TWSA through the first drought period (up to 2009, see Fig. S1), the rapid increase in TWSA from 2010 associated with higher rainfall, a decline from around 2012 due to the re-emergence of drought conditions, and the rapid decline during the recent drought after conditions had eased in 2016 (Fig. S1). FD underestimates the magnitude of monthly TWSA variance (standard deviation, SD = 37.18 mm) compared to GRACE (47.74 mm) or GW (47.67 mm), particularly during the wetter periods (2000, 2011–2016) and the first ∼2 years of the droughts (2001–2002, 2017–2018) (Fig. 1a). This underestimation in FD compared to GW is linked to the lack of aquifer water storage in the FD simulations, which provides a reservoir of water that changes slowly and has a memory of previous wet/dry climate conditions (Fig. 1a).

Figure 1(a) Total water storage anomaly (TWSA) during 2000–2019 and (b) accumulated PE for the two droughts over SE Australia. In panel (a), observations from GRACE are shown in black, the GW simulation in green, FD in orange, and the aquifer water storage anomaly in GW in blue. The shading in panel (a) highlights the two drought periods. The left top corner of panel (a) displays the correlation (r) between GRACE and GW/FD, as well as the standard deviation (SD, mm) of GRACE, GW, and FD over the periods when GRACE and the simulations coincide. Panel (b) shows the accumulated PE for two periods; the dark lines show the 2001–2009 Millennium Drought (MD) and the light lines show the 2017–2019 recent drought (RD). The correlation (r) between the P and E is shown in the legend of panel (b).


Figure 1b shows the accumulated precipitation (P) minus evapotranspiration (E) over the two drought periods. GW increases the evapotranspiration relative to FD such that the accumulated PE decreases from about 786 to 455 mm during the Millennium Drought, which is within the range of DOLCE (460 mm) and GLEAM (97 mm) estimates. A similar result, although over a much shorter period, is also apparent for the recent drought (Fig. 1b). The lower PE in GW suggests that the presence of groundwater storage can alleviate the vegetation water stress during droughts and reduces the reliance of E on P, indicated by a small reduction in the correlation (r) between E and P from 0.28 in FD to 0.24 in GW for MD and a reduction from 0.42 to 0.37 for RD (Fig. 1b). Although the evapotranspiration products display some differences, the GW simulations are closer overall to both the DOLCE and the GLEAM observationally constrained estimates. The better match of GW than FD to the two evapotranspiration products implies that adding groundwater improves the simulations during droughts, whilst the remaining mismatch would tend to suggest further biases in simulated evapotranspiration arising from multiple sources (e.g. a mismatch in leaf area index or contributions from the understorey). The difference in E is also demonstrated spatially in Fig. S3. During the Millennium Drought, the GW simulations show a clear improvement over FD in two aspects. GW shows smaller biases in E along the coast where FD underestimates E strongly (Fig. S3b–c). The areas where E is underestimated are also smaller in extent in GW, suggesting that GW overall reduces the dry bias. The magnitude of the bias in GW reaches around 300 mm over small areas of SE Australia, while in the FD simulations biases are larger, reaching 400 mm over a larger area. Plant photosynthesis assimilation rates are associated with transpiration via stomata conductance. Figure S4 presents the spatial maps of gross primary productivity (GPP) during the two droughts. GW simulations increase carbon uptake by 50–300 g C yr−1 along the coasts (Fig. S4c, f). However, since CABLE uses a prescribed LAI and does not simulate any feedback between water availability and plant growth (e.g. defoliation) and its impact on GPP, we only focus on how GW influences evapotranspiration and the surface energy balance in the subsequent sections.

Overall, Figs. 1 and S3 indicates that representing groundwater in the model improves the simulation of the interannual variability in the terrestrial water cycle and storage, particularly during droughts.

3.2 The role of groundwater in sustaining evapotranspiration during droughts

We next explore the mechanisms by which including groundwater modifies the simulation of evapotranspiration. Figure 2 displays the overall influence of groundwater on water fluxes during the recent drought. GW simulates 50–200 mm yr−1 more E over coastal regions where there is high tree cover (Fig. 2a; see Fig. S2a for land cover). Adding groundwater also increases E in most other regions, although the impact is negligible in many inland and non-forested regions (i.e. west of 145 E). We identified a clear connection between E (Fig. 2a) and the simulated WTD in the GW simulations (Fig. S5). GW simulates 110 mm yr−1 more E when the WTD is shallower than 5 m deep, 22 mm yr−1 when the WTD is 5–10 m deep, but only 3 mm yr−1 more when the WTD is below 10 m. Higher transpiration (Et) in GW explains 78 % of the evapotranspiration difference between GW and FD where WTD is shallower than 5 m (Fig. 2b). This is confirmed by the change in the soil evaporation (ΔEs) (Fig. 2c) where adding groundwater increases Es by negligible amounts over most of SE Australia but by up to 25 mm yr−1 in regions underlain by shallow groundwater (Fig. S5), which is consistent with field observations that indicate that Es can be substantial under conditions of a very shallow water table (Thorburn et al., 1992). In the very shallow WTD areas, the excess Es in GW results from the capillary rise of moisture from the shallow groundwater to the surface.

Figure 2The overall influence of groundwater during the recent drought. Panels (a)(e) are the difference (GW–FD) in evapotranspiration (ΔE), transpiration (ΔEt), soil evaporation (ΔEs), vertical drainage (ΔDr) and recharge from the aquifer to soil column (ΔQrec), respectively. In the bottom right of panels (a)(e), the average of each variable over selected water table depths (WTDs) is provided. Panel (f) is the maximum nighttime water stress factor difference (Δβ) between 03:00 (i.e. predawn when the soil is relatively moist following capillary lift overnight) and 21:00 the previous day. We only include rainless nights in January 2019 to calculate Δβ to remove any influence of overnight rainfall.

A significant factor in explaining how groundwater influences E is through changes in vertical drainage and recharge from the aquifer to the soil column. Figure 2d shows that the vertical drainage (Dr) both increases and decreases depending on the location. The addition of groundwater reduces vertical drainage by 74 mm yr−1 where WTD is shallower than 5 m. In some regions, the drainage increases with the inclusion of groundwater by up to 100 mm yr−1, especially in the areas where WTD is ∼5 m. This is associated with the WTD being slightly below the bottom of the soil column (4.6 m). When the groundwater aquifer is nearly full in GW, the wetter soil in the bottom layer leads to a much higher hydraulic conductivity in GW than in FD, leading to higher vertical drainage in GW and a positive ΔDr. Inland, where the WTD tends to be much deeper there is no significant difference in Dr between GW and FD.

Figure 2e shows the difference in recharge into the upper soil column (ΔQrec) between GW and FD. The recharge from the aquifer into the bottom soil layer provides 17 mm yr−1 of extra moisture in GW, where the WTD is between 5–10 m, and 10 mm yr−1 where the WTD is < 10 m, partially explaining the changes in E and Et in areas with a deep WTD. However, there is no significant ΔQrec in regions with a shallow WTD (∼5 mm yr−1), suggesting the influence of groundwater is mainly via reduced drainage in these locations. Recharge from the aquifer to the soil column can only occur when WTD is below the soil column (bottom boundary at 4.6 m depth). If WTD is shallow and within the soil column, the interface is saturated and no recharge from the aquifer to the soil column can occur and water only moves downwards by gravity.

The combined impact of reduced drainage in GW (Fig. 2d) and recharge from the aquifer into the root zone (Fig. 2e) is an increased water potential gradient between the drier top soil layers and the wetter deep soil layers, encouraging overnight capillary rise. Taking the hot and dry January 2019 as an example, when the compound events occurred, Fig. 2f shows the maximum water stress factor difference (Δβ) overnight (between 21:00 and 03:00, i.e. predawn when soil is relatively moist following capillary lift overnight). We only consider rainless nights to exclude the impact of drainage induced by precipitation. The water stress factor (β) is based on the root distribution and moisture availability in each soil layer and represents the soil water stress on transpiration as water becomes limiting. Figure 2f implies that while the redistribution of moisture is small overall, in some locations it can reduce moisture stress by up to 4 %–6 %.

3.3 The impact of groundwater during heatwaves

We next explore whether the higher available moisture due to the inclusion of groundwater enables the canopy to cool itself via evapotranspiration during heatwaves by examining the temperature difference between the simulated canopy temperature (Tcanopy, C) and the forced air temperature (Tair, C). We focus on the forested regions (Fig. S2a) as the role of groundwater in enhancing plant water availability was shown to be largest in these regions (Fig. 2).

Figure 3a shows the average TcanopyTair (ΔT, C) over the forested regions for summer heatwaves from the GW and FD simulations, with the grey line indicating the median ΔT difference. During heatwaves, the inclusion of groundwater moistens the soil and supports higher transpiration, cooling the canopy and reducing ΔT relative to FD by up to 0.76 C (e.g. summer heatwaves in 2013). As the drought lengthens in time, the depletion of moisture gradually reduces this effect, from an average reduction of 0.52 C of the first 3 years to 0.16 C of the last 3 years in the Millennium Drought (Fig. 3a). The impact of groundwater is clear in the evaporative fraction (Fig. 3b) where in periods of higher rainfall (e.g. 2010–2011; Fig. S1) and at the beginning of a drought (2001, 2017) the evaporative fraction (EF) is higher (0.03 to 0.18). This implies that more of the available energy is exchanged with the atmosphere in the form of latent rather than sensible heat. However, the strength of the cooling effect decreases as the droughts extend and the transpiration difference (ΔEt, mm d−1) diminishes quickly (Fig. 3c) because the vegetation becomes increasingly water-stressed (Fig. 3d) which consequently limits transpiration. For all variables (ΔT, EF, Et, and β), the difference between GW and FD is greatest during the wetter periods (e.g. 2013) and in the first 1–2 years of the multi-year drought (2001–2002 for the Millennium Drought or 2017–2018 for the recent drought). After the drought becomes well established, the FD and GW simulations converge as depleting soil moisture reservoirs reduce the impact of groundwater on canopy cooling and evaporative fluxes.

Figure 3Groundwater-induced differences in (a) TcanopyTair (ΔT), (b) evaporative fraction (EF), (c) transpiration (Et), and (d) water stress factor (β) during 2000–2019 summer heatwaves over forested areas (the green region in Fig. S2a). The left y axis is the scale for boxes. The blue boxes refer to the GW experiment and the red boxes to FD. For each box, the middle line is the median, the upper border is the 75th percentile, and the lower border is 25th percentile. The right y axis is the scale for the grey lines which display the difference in the medians (GW–FD). The shadings highlight the two drought periods.


Figure 4Land response to heatwaves during the recent drought. Panels (a)(c) are the mean TcanopyTair (ΔT), evaporative fraction (EF), and soil water stress factor (β) in GW, respectively, during 2017–2019 summer heatwaves. Panels (d)(f) are the difference (GW–FD) of Tcanopy, EF, and β. In the bottom right of each plot, the average of each variable over selected water table depths (WTDs) is provided. Note that the colour bar is switched between (d) and (e)(f).

Figure 4a shows the spatial map of ΔT simulated in GW during heatwaves in the 2017–2019 drought. It indicates that both land cover type (Fig. S2a) and WTD (Fig. S5) contribute to the ΔT pattern. The evaporative cooling via transpiration is stronger over the forested areas compared to cropland or grassland and stronger in the regions with a wetter soil associated with a shallower WTD. However, EF is mainly determined by WTD (compare Figs. 4b and S5). Inland, where the WTD is deeper and the soil is drier, most of the net radiation absorbed by the land surface is partitioned into sensible rather than latent heat (Fig. 4b). However, in the coastal regions with a shallow WTD, the wetter soil reduces the water stress (Fig. 4c), enables a higher EF (Fig. 4b), and alleviates heat stress on the leaves (Fig. 4a). Along the coast where WTD is shallow, GW simulates a cooler canopy temperature due to the higher evaporative cooling (Fig. 4e), which is the consequence of a lower soil water stress (Fig. 4f) linked to the influence of groundwater (Fig. S5).

Figure 5 shows the density scatter plot of ΔT versus WTD in SE Australia forested areas during heatwaves in 2000–2019. A shallow WTD moderates the temperature difference between the canopy and the ambient air during heatwaves leading to a smaller temperature difference. Meanwhile, as the WTD increases, due to the limited rooting depth in the model, the ability of the groundwater to support transpiration and offset the impact of high air temperatures is reduced. Figure 5 shows a large amount of variations but nonetheless implies a threshold of ∼6 m, whereafter there is a decoupling and little influence from groundwater during heatwaves. However, the absolute value of the threshold is likely CABLE-specific and associated with the assumption of a 4.6 m soil depth, which also sets the maximum rooting depth (roots can only extend to the bottom of the soil and cannot directly access the groundwater aquifer in CABLE). The CABLE soil depth comes from observational evidence of most roots being situated within the top 4.6 m (Canadell et al. 1996). Since the model assumes no roots exist in the groundwater aquifer, when the water table is below this depth, the water fluxes become largely uncoupled between the soil column and the groundwater aquifer, leading to a negligible impact of GW below ∼6 m depth.

Figure 5A density scatter plot of TcanopyTair (ΔT) versus water table depth (WTD) in GW simulations over forested areas in all heatwaves during 2000–2019. Every tree pixel on each heatwave day accounts for one record, and the darker colours show higher recorded densities.


3.4 The impact of groundwater during the drought and heatwave compound events

To examine the influence of groundwater on heatwaves occurring simultaneously with drought, we focus on a case study of the record-breaking heatwaves in January 2019, which is the hottest month on record for the study region (Bureau of Meteorology, 2019). The unprecedented prolonged heatwave period started in early December 2018 and continued through January 2019 with three peaks. We select 2 d (15 and 25 January 2019), when heatwaves spread across the study region, from the second and third heatwave phases (Fig. S6).

We evaluate CABLE Tcanopy against MODIS LST observations, concentrating on forested areas where MODIS LST should more closely reflect vegetation canopy temperatures, but note that this comparison is not direct as the satellite estimate will contain contributions from the understorey and soil. Figure 6a–b show the good-quality MODIS LST minus Tair at 14:00 (ΔTMOD_14:00) over forested regions on the 15 and 25 January 2019, and Fig. 6c–d display the matching GW-simulated ΔT at 14:00 (ΔTGW_14:00). Overall, ΔTGW_14:00 increases from the coast to the interior in both heatwaves, consistent with the ΔTMOD_14:00 pattern in both heatwaves, although ΔTGW_14:00 appears to be biased high relative to ΔTMOD_14:00 along the coastal forests (Fig. S7a–b).

Figure 6The simulation of two heatwaves on 15 (left column) and 25 January 2019 (right column). The first row shows the difference between MODIS land surface temperature (LST) and Tair at 14:00 (ΔTMOD_14:00) (only forested areas with good LST quality data are displayed). The second row is the GW simulation of ΔT at 14:00 (ΔTGW_14:00). The third row is the difference of ΔT at 14:00 between GW and FD simulations (ΔTGW_14:00-ΔTFD_14:00). The last row is the same as the third row but for the difference between the DR and GW simulations (ΔTDR_14:00-ΔTGW_14:00). Note that the comparison between GW/FD/DR and MODIS LST is shown in Fig. S7.

Figure 6e–f show the ΔT14:00 difference between GW and FD. Access to groundwater can reduce canopy temperature by up to 5 C, in particular where the WTD is shallow. While reductions of 5 C are clearly limited in spatial extent, the overall pattern of cooling is quite widespread and coincident with the groundwater-induced Et increase (Fig. S8a–b), implying a reduction in heat stress along coastal regions with a shallow WTD during heatwaves. Generally, GW matches MODIS LST better than FD despite the bias in both simulations (compare Figs. S7a–b and S7c–d). Nevertheless, the temperature reduction between GW and FD is still modest (< 1 C) for most of the forested regions. This may be related to the shallow root distribution assumed in many LSMs, which prevents roots from directly accessing the moisture stored in the deeper soil (note, CABLE assumes 92 % of a forest's roots are in the top 0.64 m, Fig. S2b). To examine this possibility, we performed the deep root (DR) sensitivity experiment which prescribed more roots in the deeper soil for forests (56 % below 0.64 m depth). Figure 6g–h illustrate the difference between ΔT14:00 in DR and ΔTGW_14:00. By enabling access to moisture in the deeper soil, the LSM simulates further cooling by 0.5–5 C across the forests associated with an Et increase of 25–250 W m−2 (Fig. S8c–d). The prescribed deeper roots also lead to an overall better simulation of ΔT at 20:00 relative to the MODIS LST (Figs. S7e–f vs. Fig. S7a–b).

Figure 7 shows the diurnal cycles of ΔT for the two selected regions (red boxes in Fig. 6) compared with the MODIS LST estimates. The region highlighted for the 15 January (Fig. 7a) has a WTD of 4–7 m, while the region highlighted for the 25 January (Fig. 7b) has a WTD < 4 m (Fig. S5). In both regions, the simulated ΔT is highest in FD, lower in GW and lowest in DR. Where the WTD is 4–7 m (Fig. 7a), the three simulated ΔT are slightly lower than ΔT calculated by MODIS LST (red squares). However, in the shallower WTD region (Fig. 7b), the simulated ΔT between experiments is more dispersed across experiments and exceeds the MODIS ΔT at both time points, implying that neglecting groundwater dynamics and deep roots is more likely to cause an overestimation of heat stress in the shallower WTD region. The shallower WTD region (Fig. 7b) tends to have a high LAI coverage, implying that the MODIS LST represents a good approximation of the canopy temperature over this region. Consequently, the lower MODIS ΔT implies that CABLE is likely underestimating transpiration, leading to an overestimation of ΔT in all three simulations.

Figure 7Diurnal cycle of ΔT on 15 (left column) and 25 January 2019 (right column) over the selected regions shown in Fig. 6. The shadings show the uncertainty in every simulation defined as one standard deviation (SD) among the selected pixels. The red dots are MODIS LST minus Tair with the uncertainty shown by the red error bars. For both regions, only pixels available in MODIS are shown.


3.5 Constraints on groundwater mediation during the compound events

We finally probe the reasons for the apparent contradiction between the large impact of groundwater on E during drought (Fig. 2a) but a smaller impact on ΔT during the compound events (Fig. 7). Figure 8 shows three factors (β, vapour pressure deficit (D), and Tair) that constrain the impact of groundwater on ΔT in CABLE during heatwaves in January 2019. Figure 8a shows the difference in ΔT between GW and FD as a function of Δβ, suggesting that the inclusion of groundwater has a large impact on ΔT when there is a coincidental and large difference in β between the GW and FD simulations. Figure 8b indicates a clear threshold at D= 3 kPa where GW and FD converge, while Fig. 8c shows a convergence threshold when Tair exceeds 32 C. Above these two thresholds, access to groundwater seemingly becomes less important in mitigating plant heat stress. There are two mechanisms in CABLE that explain this behaviour. First, as D increases, CABLE predicts that stomata begin to close following a square root dependence (De Kauwe et al., 2015; Medlyn et al., 2011). Second, as Tair increases, photosynthesis becomes inhibited as the temperature exceeds the optimum for photosynthesis. In both instances, evaporative cooling is reduced, regardless of the root zone moisture state dictated by groundwater access. That is to say, access to groundwater has limited capacity to directly mediate the heat stress on plants during a compound event when the air is very dry or very hot.

Figure 8Density scatter plots showing the three factors that influence the difference in TcanopyTair between GW and FD (ΔT, expressed as GW–FD difference). Panel (a) is ΔT difference against the β difference (GW–FD) (Δβ), (b) is ΔT difference against vapour pressure deficit (D), and (c) is ΔT difference against Tair. Each point corresponds to a tree pixel on a heatwave day in January 2019. The darker colours illustrate where the records are denser. The correlation (r) between the x and y axes is shown in the bottom left of each panel.


4 Discussion

In the absence of direct measurements, we used the CABLE-GW LSM, constrained by satellite observations to investigate how groundwater influences ecosystems under conditions of co-occurring droughts and heatwaves. We found that the influence of groundwater was most important during the wetter periods and the first  2 years of a multi-year drought (∼2001–2002 and 2017–2018; Figs. 1 and 3). This primarily occurred via impedance of gravity-driven drainage (Fig. 2d) but also via capillary rise from the groundwater aquifer (Fig. 2e). This moistening enabled the vegetation to sustain higher E for at least a year (Fig. S9). As the droughts progressed into multi-year events, the impact of groundwater diminished due to a depletion of soil moisture stores regardless of whether groundwater dynamics were simulated.

When a heatwave occurs during a drought, and in particular early in a drought, the extra transpiration enabled by representing groundwater dynamics helps reduce the heat stress on vegetation (e.g. the reduction of 0.64 C of ΔT over the forests in 2002, Fig. 3a). This effect is particularly pronounced in regions with a shallower WTD (e.g. where the groundwater was within the first 5 m, there was a 0.5 C mean reduction in ΔT in the recent drought; Fig. 4d). Importantly, the role played by groundwater diminishes as the drought lengthens beyond 2 years (Fig. 3). Additionally, either the lack of deep roots or stomatal closure caused by high D/Tair can reduce the additional transpiration induced by groundwater. The latter plant physiology feedback dominates during heatwaves co-occurring with drought, even if the groundwater's influence has increased root zone water availability.

Our results highlight the impact of groundwater on both land surface states (e.g. soil moisture) and on surface fluxes and how this impact varies with the length and intensity of droughts and heatwaves. The results imply that the dominant mechanism by which groundwater buffered transpiration was through impeding gravity-driven drainage. We found a limited role for upward water movement from the aquifer due to simulated shallow WTD (which was broadly consistent with the observations in Fan et al., 2013). Further work will be necessary to understand how groundwater interacts with droughts and heatwaves and what these interactions mean for terrestrial ecosystems and the occurrence of the compound extreme events, particularly under the projection of intensifying droughts (Ukkola et al., 2020) and heatwaves (Cowan et al., 2014).

4.1 Changes in the role of groundwater in multi-year droughts

Groundwater is the slowest part of the terrestrial water cycle to change (Condon et al., 2020) and can have a memory of multi-year variations in rainfall (Martínez-de la Torre and Miguez-Macho, 2019; Martinez et al., 2016a). Our results show that the effect of groundwater on the partitioning of available energy between latent and sensible heat fluxes is influenced by the length of drought. As the drought extends in time, the extra E sustained by groundwater decreases (e.g. during the Millennium Drought, Fig. S9). The role of a drying landscape in modifying the partitioning of available energy between latent and sensible heat fluxes is well known and has been extensively studied (Fan, 2015; Miralles et al., 2019; Seneviratne et al., 2010). Our results add to the knowledge by quantifying the extent of the groundwater control, and eliciting the timescales of influence and the mechanisms at play. The importance of vegetation–groundwater interactions on multi-year timescales has been identified previously. Humphrey et al. (2018) hypothesised that climate models may underestimate the amplitude of global net ecosystem exchange because of a lack of deep-water access. Our regionally based results support this hypothesis and in particular highlight the importance of groundwater for explaining the amplitude of fluxes in wet periods as well as sustaining evapotranspiration during drought (Fig. 1).

4.2 Implications for land–atmosphere feedbacks during compound events

Our results show that during drought–heatwave compound events, the existence of groundwater eases the heat stress on the forest canopy and reduces the sensible heat flux to the atmosphere. This has the potential to reduce heat accumulating in the boundary layer and help ameliorate the intensity of a heatwave (Keune et al., 2016; Zipper et al., 2019). The presence of groundwater helps dampen a positive feedback loop whereby during drought–heatwave compound events, the high exchange of sensible and low exchange of latent heat can heat the atmosphere and increase the atmospheric demand for water (De Boeck et al., 2010; Massmann et al., 2019), intensifying drying (Miralles et al., 2014). The lack of groundwater in many LSMs suggests a lack of this moderating process and consequently a risk of overestimating the positive feedback on the boundary layer in coupled climate simulations. Our results show that neglecting groundwater leads to an average overestimate in canopy temperature by 0.2–1 C where the WTD is shallow (Fig. 4d) but as much as 5 C in single heatwave events (Fig. 6e–f), leading to an increase in the sensible heat flux (Fig. 4e).

The capacity of groundwater to moderate this positive land–atmosphere feedback is via modifying soil water availability. Firstly, soil water availability influenced by WTD affects how much water is available for E. In the shallow WTD regions, the higher soil water is likely to suppress the mutual enhancement of droughts and heatwaves (Keune et al., 2016; Zipper et al., 2019), particularly early in a drought. However, this suppression becomes weaker as the WTD deepens, in particular at depths beneath the root zone (e.g. 4.6 m in CABLE-GW) or as a drought lengthens. Our results imply that the land amplification of heatwaves is likely stronger in the inland regions (Hirsch et al., 2019), where the WTD is lower than 5 m and the influence of groundwater diminishes (Fig. S5), and once a drought has intensified significantly.

On a dry and hot heatwave afternoon, plant physiology feedbacks to high D and high Tair dominate transpiration and reduce the influence of groundwater in moderating heatwaves. In CABLE, stomatal closure occurs either directly due to high D (> 3 kPa) (De Kauwe et al. 2015) or indirectly due to biochemical feedbacks on photosynthesis at high Tair (> 32 C) (Kowalczyk et al., 2006); both processes reduce transpiration to near zero, eliminating the buffering effect of groundwater on canopy temperatures. While the timing of the onset of these physiology feedbacks varies across LSMs due to different parameterised sensitivities of stomatal conductance to atmospheric demand (Ball et al., 1987; Leuning et al., 1995) and different temperature dependence parameterisations (Badger and Collatz, 1977; Bernacchi et al., 2001; Crous et al., 2013), importantly, stomatal closure during heat extremes would be model invariant.

4.3 Uncertainties and future directions

Our study uses a single LSM, and consequently the parameterisations included in CABLE-GW influence the quantification of the role of groundwater on droughts and heatwaves. We note CABLE-GW has been extensively evaluated for water cycle processes (Decker, 2015; Decker et al., 2017; Mu et al., 2021; Ukkola et al., 2016b), but evaluation for groundwater interactions remains limited due to the lack of suitable observations (e.g. regional WTD monitoring or detailed knowledge of the distribution of root depths). Figure 1 gives us confidence that CABLE-GW is performing well, based on the evaluation against the GRACE, DOLCE, and GLEAM products, as well as previous work that showed the capacity of CABLE-GW to simulate E well (Decker, 2015; Decker et al., 2017). However, we also note that key model parameterisations that may influence the role of groundwater are particularly uncertain.

We need to be cautious about the “small” groundwater impact on the canopy temperature and associated turbulent energy fluxes during high D or high Tair (Figs. 3, 4, 6 and 7). The thresholds of D and Tair currently assumed by LSMs are in fact likely to be species specific. Australian trees in particular have evolved a series of physiological adaptations to reduce the negative impact of heat extremes. It is important to note that most LSMs parameterise their stomatal response to D for moderate ranges (< 2 kPa), which leads to significant biases at high D (Yang et al., 2019), a feature common in Australia and during heatwaves in general. New theory is needed to ensure that models adequately capture the full range of stomatal response to variability in D (low and high ranges). Similarly, while there is strong evidence to suggest that the optimum temperature for photosynthesis does not vary predictably with the climate of species origin (Kumarathunge et al., 2019) (implying model parameterisations do not need to vary with species), findings from studies do vary (Cunningham and Reed, 2002; Reich et al., 2015). Moreover, evidence that plants acclimate their photosynthetic temperature response is strong (Kattge and Knorr, 2007; Kumarathunge et al., 2019; Mercado et al., 2018; Smith et al., 2016; Smith and Dukes, 2013). As a result, it is likely that LSMs currently underestimate groundwater influence during heatwaves due to the interaction with plant physiology feedbacks. This is a key area requiring further investigation. For example, Drake et al. (2018) demonstrated that during a 4 d heatwave > 43 C, Australian Eucalyptus parramattensis trees did not reduce transpiration to zero as models would commonly predict, allowing the trees to persist unharmed in a whole-tree chamber experiment. Although De Kauwe et al. (2019) did not find strong support for this phenomenon across eddy covariance sites, if this physiological response is common across Australian woodlands, it would change our view on the importance of soil water availability (therefore groundwater) for the evolution of heatwave or even compound events. Coupled model sensitivity experiments may be important to determine the magnitude that such a physiological feedback would present and could guide the direction of future field/manipulation experiments.

Root distribution and root function and thereby how roots utilise groundwater are uncertain in models (Arora and Boer, 2003; Drewniak, 2019; Wang et al., 2018; Warren et al., 2015) and indeed in observations (Fan et al., 2017; Jackson et al., 1996; Schenk and Jackson, 2002). Models often ignore how roots forage for water and respond to moisture heterogeneity, limiting the model's ability to accurately reflect the plant usage of groundwater (Warren et al., 2015). In LSMs, roots are typically parameterised using a fixed distribution and normally ignore water uptake from deep roots. This assumption neglects any climatological impact of root distribution and the differentiation in root morphology and function (fine roots vs. tap roots), leading to a potential underestimation of groundwater utilisation in LSMs (see our deep root experiment, Fig. 6g–h). This assumption may be particularly problematic in Australia where vegetation has developed significant adaptation strategies to cope with both extreme heat and drought, including deeply rooted vegetation that can access groundwater (Bartle et al., 1980; Dawson and Pate, 1996; Eamus et al., 2015; Eberbach and Burrows, 2006; Fan et al., 2017). We also note that CABLE does not directly consider hydraulic redistribution, defined as the passive water movement via plant roots from moister to drier soil layers (Burgess et al., 1998; Richards and Caldwell, 1987). Neglecting hydraulic redistribution has the potential to underestimate the groundwater transported upwards and understate the importance of groundwater on ecosystems.

On the atmosphere side, the existence of groundwater increases the water flux from the land to the atmosphere, particularly in regions of shallow WTD, during the first 1–2 years of a drought. This has the potential to moisten the lower atmosphere and may encourage precipitation (Anyah et al., 2008; Jiang et al., 2009; Martinez et al., 2016b; Maxwell et al., 2011). However, our experiments are uncoupled from the atmosphere so while there is the potential for the higher E to affect the boundary layer moisture (Bonetti et al., 2015; Gilbert et al., 2017; Maxwell et al., 2007), clouds, and precipitation, we cannot conclude that it would until we undertake future coupled simulations.

Finally, we note we have focused on the role of groundwater in a natural environment. Humans extract large quantities of groundwater in many regions (Döll et al., 2014; Wada, 2016). Adding human management of groundwater into LSMs enables an examination of how this affects the vulnerability of ecosystems to heatwaves and drought and may ultimately identify those vulnerable ecosystems close to tipping points that are priorities for protection.

5 Conclusions

In conclusion, we used the CABLE LSM, constrained by satellite observations, to explore the timescales and extent to which groundwater influences vegetation function and turbulent energy fluxes during multi-year droughts. We showed that groundwater moistened the soil during the first  2 years of a multi-year drought, which enabled the vegetation to sustain higher evaporation (50–200 mm yr−1 over the coastal forest regions) during drought onset. This cooled the forest canopy on average by 0.03–0.76 C in heatwaves during 2001–2019 and by as much as 5 C in regions of shallow water table depths in the heatwave in January 2019, helping to moderate the heat stress on vegetation during heatwaves. However, the ability of groundwater to buffer vegetation function varied with the length and intensity of droughts and heatwaves, with its influence decreasing with prolonged drought conditions. Importantly, we also demonstrated that the capacity of the groundwater to buffer evaporative fluxes during heatwaves is constrained by plant physiology feedbacks which regulate stomatal control, irrespective of soil water status. Given the increased risk of regional heatwaves and droughts in the future, the role of groundwater on land–atmosphere feedbacks and on terrestrial ecosystems needs to be better understood in order to constrain future projections.

Code and data availability

The CABLE code is available at (NCI, 2021) after registration. Here, we use CABLE revision r7765. Scripts for plotting and processing model outputs are available at (Mu, 2021). The DOLCE version 2 dataset is available from the NCI data catalogue at (Hobeichi et al., 2021b). The GRACE dataset is available at (Wiese et al., 2019). The GLEAM version 3.5 dataset is available at (GLEAM, 2021). The datasets of MOD11A1 (, Wan et al., 2015a) and MYD11A1 (, Wan et al., 2015b) were acquired from the NASA Land Processed Distributed Active Archive Center (LP DAAC), located in the USGS Earth Resources Observation and Science (EROS) Center in Sioux Falls, South Dakota, USA (, LP DAAC, 2021).


The supplement related to this article is available online at:

Author contributions

MM, MGDK, AJP, and AMU conceived the study, designed the model experiments, investigated the simulations, and drafted the paper. SH and PRB provided the evaluation and the meteorology forcing datasets. All authors participated in the discussion and revision of the paper.

Competing interests

The contact author has declared that neither they nor their co-authors have any competing interests.


Publisher’s note: Copernicus Publications remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Special issue statement

This article is part of the special issue “Understanding compound weather and climate events and related impacts (BG/ESD/HESS/NHESS inter-journal SI)”. It is not associated with a conference.


Mengyuan Mu, Martin G. De Kauwe, Andy J. Pitman, Anna M. Ukkola, and Sanaa Hobeichi acknowledge support from the Australian Research Council (ARC) Centre of Excellence for Climate Extremes (CE170100023). Mengyuan Mu acknowledges support from the UNSW University International Postgraduate Award (UIPA) scheme. Martin G. De Kauwe and Andy J. Pitman acknowledge support from the ARC Discovery Grant (DP190101823) and Anna M. Ukkola support from the ARC Discovery Early Career Researcher Award (DE200100086). Martin G. De Kauwe acknowledges support from the NSW Research Attraction and Acceleration Program (RAAP). We thank the National Computational Infrastructure at the Australian National University, an initiative of the Australian Government, for access to supercomputer resources. Mengyuan Mu thanks the University of Nanjing for hosting her research through 2020.

Financial support

This research has been supported by the Australian Research Council (grant nos. CE170100023, DP190101823, and DE200100086) and the NSW Research Attraction and Acceleration Program.

Review statement

This paper was edited by Jakob Zscheischler and reviewed by three anonymous referees.


Abramowitz, G., Leuning, R., Clark, M., and Pitman, A.: Evaluating the performance of land surface models, J. Climate, 21, 5468–5481,, 2008. 

Abramowitz, G.: Towards a public, standardized, diagnostic benchmarking system for land surface models, Geosci. Model Dev., 5, 819–827,, 2012. 

Allen, C. D., Macalady, A. K., Chenchouni, H., Bachelet, D., McDowell, N., Vennetier, M., Kitzberger, T., Rigling, A., Breshears, D. D., Hogg, E. H. (Ted), Gonzalez, P., Fensham, R., Zhang, Z., Castro, J., Demidova, N., Lim, J.-H., Allard, G., Running, S. W., Semerci, A., and Cobb, N.: A global overview of drought and heat-induced tree mortality reveals emerging climate change risks for forests, Forest Ecol. Manag., 259, 660–684,, 2010. 

Allen, C. D., Breshears, D. D., and McDowell, N. G.: On underestimation of global vulnerability to tree mortality and forest die-off from hotter drought in the Anthropocene, Ecosphere, 6, 129,, 2015. 

Anyah, R. O., Weaver, C. P., Miguez-Macho, G., Fan, Y., and Robock, A.: Incorporating water table dynamics in climate modeling: 3. Simulated groundwater influence on coupled land-atmosphere variability, J. Geophys. Res., 113, D07103,, 2008. 

Arora, V. K. and Boer, G. J.: A representation of variable root distribution in dynamic vegetation models, Earth Interact., 7, 1–19,<0001:AROVRD>2.0.CO;2, 2003. 

Badger, M. R. and Collatz, G. J.: Studies on the kinetic mechanism of RuBP carboxylase and oxygenase reactions, with particular reference to the effect of temperature on kinetic parameters, Year book–Carnegie Institution of Washington, Baltimore, Maryland, USA, 355–361 pp., 1977. 

Ball, J. T., Woodrow, I. E., and Berry, J. A.: A model predicting stomatal conductance and its contribution to the control of photosynthesis under different environmental conditions, in Progress in Photosynthesis Research, Springer Netherlands, Dordrecht, the Netherlands, 221–224,, 1987. 

Bartle, G. A., Murray, A. M., and Macpherson, D. K.: The distribution of root length, and the limits to flow of soil water to roots in a dry sclerophyll forest, Forest Sci., 26, 656–664, 1980. 

Bernacchi, C. J., Singsaas, E. L., Pimentel, C., Portis Jr., A. R., and Long, S. P.: Improved temperature response functions for models of Rubisco-limited photosynthesis, Plant. Cell Environ., 24, 253–259,, 2001. 

Best, M. J., Abramowitz, G., Johnson, H. R., Pitman, A. J., Balsamo, G., Boone, A., Cuntz, M., Decharme, B., Dirmeyer, P. A., Dong, J., Ek, M., Guo, Z., Haverd, V., van den Hurk, B. J. J., Nearing, G. S., Pak, B., Peters-Lidard, C., Santanello, J. A., Stevens, L., and Vuichard, N.: The plumbing of land surface models: benchmarking model performance, J. Hydrometeorol., 16, 1425–1442,, 2015. 

Birami, B., Gattmann, M., Heyer, A. G., Grote, R., Arneth, A., and Ruehr, N. K.: Heat waves alter carbon allocation and increase mortality of Aleppo Pine under dry conditions, Front. For. Glob. Chang., 1, 8,, 2018. 

Bonetti, S., Manoli, G., Domec, J. C., Putti, M., Marani, M., and Katul, G. G.: The influence of water table depth and the free atmospheric state on convective rainfall predisposition, Water Resour. Res., 51, 2283–2297,, 2015. 

Bureau of Meteorology: Special Climate Statement 43-extreme heat in January 2013, available at: (last access: 29 August 2021), 2013. 

Bureau of Meteorology: Special Climate Statement 61-exceptional heat in southeast Australia in early 2017, available at: (last access: 29 August 2021), 2017. 

Bureau of Meteorology: Special Climate Statement 68 – widespread heatwaves during December 2018 and January 2019, available at: (last access: 29 August 2021), 2019. 

Burgess, S. S. O., Adams, M. A., Turner, N. C., and Ong, C. K.: The redistribution of soil water by tree root systems, Oecologia, 115, 306–311,, 1998. 

Canadell, J., Jackson, R. B., Ehleringer, J. B., Mooney, H. A., Sala, O. E., and Schulze, E.-D.: Maximum rooting depth of vegetation types at the global scale, Oecologia, 108, 583–595,, 1996. 

Ciais, P., Reichstein, M., Viovy, N., Granier, A., Ogée, J., Allard, V., Aubinet, M., Buchmann, N., Bernhofer, C., Carrara, A., Chevallier, F., De Noblet, N., Friend, A. D., Friedlingstein, P., Grünwald, T., Heinesch, B., Keronen, P., Knohl, A., Krinner, G., Loustau, D., Manca, G., Matteucci, G., Miglietta, F., Ourcival, J. M., Papale, D., Pilegaard, K., Rambal, S., Seufert, G., Soussana, J. F., Sanz, M. J., Schulze, E. D., Vesala, T., and Valentini, R.: Europe-wide reduction in primary productivity caused by the heat and drought in 2003, Nature, 437, 529–533,, 2005. 

Condon, L. E., Atchley, A. L., and Maxwell, R. M.: Evapotranspiration depletes groundwater under warming over the contiguous United States, Nat. Commun., 11, 873,, 2020. 

Cosby, B. J., Hornberger, G. M., Clapp, R. B., and Ginn, T. R.: A statistical exploration of the relationships of soil moisture characteristics to the physical properties of soils, Water Resour. Res., 20, 682–690,, 1984. 

Cowan, T., Purich, A., Perkins, S., Pezza, A., Boschat, G., and Sadler, K.: More frequent, longer, and hotter heat waves for Australia in the twenty-first century, J. Climate, 27, 5851–5871,, 2014. 

Crous, K. Y., Quentin, A. G., Lin, Y.-S., Medlyn, B. E., Williams, D. G., Barton, C. V. M., and Ellsworth, D. S.: Photosynthesis of temperate Eucalyptus globulus trees outside their native range has limited adjustment to elevated CO2 and climate warming, Glob. Change Biol., 19, 3790–3807,, 2013. 

Cunningham, S. C. and Read, J.: Do temperate rainforest trees have a greater ability to acclimate to changing temperatures than tropical rainforest trees?, New Phytol., 157, 55–64,, 2003. 

Dawson, T. E. and Pate, J. S.: Seasonal water uptake and movement in root systems of Australian phraeatophytic plants of dimorphic root morphology: a stable isotope investigation, Oecologia, 107, 13–20,, 1996. 

De Boeck, H. J., Dreesen, F. E., Janssens, I. A., and Nijs, I.: Climatic characteristics of heat waves and their simulation in plant experiments, Glob. Change Biol., 16, 1992–2000,, 2010. 

Decker, M.: Development and evaluation of a new soil moisture and runoff parameterization for the CABLE LSM including subgrid-scale processes, J. Adv. Model. Earth Syst., 7, 1788–1809,, 2015. 

Decker, M., Or, D., Pitman, A. and Ukkola, A.: New turbulent resistance parameterization for soil evaporation based on a pore-scale model: Impact on surface fluxes in CABLE, J. Adv. Model. Earth Sy., 9, 220–238,, 2017. 

De Kauwe, M. G., Kala, J., Lin, Y.-S., Pitman, A. J., Medlyn, B. E., Duursma, R. A., Abramowitz, G., Wang, Y.-P., and Miralles, D. G.: A test of an optimal stomatal conductance scheme within the CABLE land surface model, Geosci. Model Dev., 8, 431–452,, 2015. 

De Kauwe, M. G., Medlyn, B. E., Pitman, A. J., Drake, J. E., Ukkola, A., Griebel, A., Pendall, E., Prober, S., and Roderick, M.: Examining the evidence for decoupling between photosynthesis and transpiration during heat extremes, Biogeosciences, 16, 903–916,, 2019. 

De Kauwe, M. G., Medlyn, B. E., Ukkola, A. M., Mu, M., Sabot, M. E. B. B., Pitman, A. J., Meir, P., Cernusak, L. A., Rifai, S. W., Choat, B., Tissue, D. T., Blackman, C. J., Li, X., Roderick, M., and Briggs, P. R.: Identifying areas at risk of drought-induced tree mortality across South-Eastern Australia, Glob. Change Biol., 26, 5716–5733,, 2020. 

D'Odorico, P. and Porporato, A.: Preferential states in soil moisture and climate dynamics, P. Natl. Acad. Sci. USA, 101, 8848–8851,, 2004. 

Döll, P., Müller Schmied, H., Schuh, C., Portmann, F. T., and Eicker, A.: Global-scale assessment of groundwater depletion and related groundwater abstractions: Combining hydrological modeling with information from well observations and GRACE satellites, Water Resour. Res., 50, 5698–5720,, 2014. 

Drake, J. E., Tjoelker, M. G., Vårhammar, A., Medlyn, B. E., Reich, P. B., Leigh, A., Pfautsch, S., Blackman, C. J., López, R., Aspinwall, M. J., Crous, K. Y., Duursma, R. A., Kumarathunge, D., De Kauwe, M. G., Jiang, M., Nicotra, A. B., Tissue, D. T., Choat, B., Atkin, O. K. and Barton, C. V. M. M.: Trees tolerate an extreme heatwave via sustained transpirational cooling and increased leaf thermal tolerance, Glob. Change Biol., 24, 2390–2402,, 2018. 

Drewniak, B. A.: Simulating dynamic roots in the energy exascale earth system land model, J. Adv. Model. Earth Sy., 11, 338–359,, 2019. 

Eamus, D. and Froend, R.: Groundwater-dependent ecosystems: the where, what and why of GDEs, Aust. J. Bot., 54, 91–96,,, 2006. 

Eamus, D., Zolfaghar, S., Villalobos-Vega, R., Cleverly, J., and Huete, A.: Groundwater-dependent ecosystems: recent insights from satellite and field-based studies, Hydrol. Earth Syst. Sci., 19, 4229–4256,, 2015. 

Eberbach, P. L. and Burrows, G. E.: The transpiration response by four topographically distributed Eucalyptus species, to rainfall occurring during drought in south eastern Australia, Physiol. Plant., 127, 483–493,, 2006. 

Fan, Y.: Groundwater in the Earth's critical zone: Relevance to large-scale patterns and processes, Water Resour. Res., 51, 3052–3069,, 2015. 

Fan, Y., Li, H., and Miguez-Macho, G.: Global patterns of groundwater table depth, Science, 339, 940–943,, 2013. 

Fan, Y., Miguez-Macho, G., Jobbágy, E. G., Jackson, R. B., and Otero-Casal, C.: Hydrologic regulation of plant rooting depth, P. Natl. Acad. Sci. USA, 114, 10572–10577,, 2017. 

Fischer, E. M., Seneviratne, S. I., Vidale, P. L., Lüthi, D., and Schär, C.: Soil moisture–atmosphere interactions during the 2003 European summer heat wave, J. Climate, 20, 5081–5099,, 2007. 

Gale, M. R. and Grigal, D. F.: Vertical root distributions of northern tree species in relation to successional status, Can. J. Forest Res., 17, 829–834,, 1987. 

Geange, S. R., Arnold, P. A., Catling, A. A., Coast, O., Cook, A. M., Gowland, K. M., Leigh, A., Notarnicola, R. F., Posch, B. C., Venn, S. E., Zhu, L., and Nicotra, A. B.: The thermal tolerance of photosynthetic tissues: a global systematic review and agenda for future research, New Phytol., 229, 2497–2513,, 2021. 

Gilbert, J. M., Maxwell, R. M., and Gochis, D. J.: Effects of water-table configuration on the planetary boundary layer over the San Joaquin River Watershed, California, J. Hydrometeorol., 18, 1471–1488,, 2017. 

GLEAM: Global Land Evaporation Amsterdam Model (GLEAM) Version 3.5,, last access: 29 August 2021. 

Griffith, S. J., Bale, C., and Adam, P.: Environmental correlates of coastal heathland and allied vegetation, Aust. J. Bot., 56, 512,, 2008. 

Haverd, V., Raupach, M. R., Briggs, P. R., Canadell, J. G., Isaac, P., Pickett-Heaps, C., Roxburgh, S. H., van Gorsel, E., Viscarra Rossel, R. A., and Wang, Z.: Multiple observation types reduce uncertainty in Australia's terrestrial carbon and water cycles, Biogeosciences, 10, 2011–2040,, 2013. 

Hengl, T., Mendes de Jesus, J., Heuvelink, G. B. M., Ruiperez Gonzalez, M., Kilibarda, M., Blagotić, A., Shangguan, W., Wright, M. N., Geng, X., Bauer-Marschallinger, B., Guevara, M. A., Vargas, R., MacMillan, R. A., Batjes, N. H., Leenaars, J. G. B., Ribeiro, E., Wheeler, I., Mantel, S., and Kempen, B.: SoilGrids250m: Global gridded soil information based on machine learning, PLoS One, 12, e0169748,, 2017. 

Hirsch, A. L., Evans, J. P., Di Virgilio, G., Perkins-Kirkpatrick, S. E., Argüeso, D., Pitman, A. J., Carouge, C. C., Kala, J., Andrys, J., Petrelli, P., and Rockel, B.: Amplification of Australian heatwaves via local land-atmosphere coupling, J. Geophys. Res.-Atmos., 124, 13625–13647,, 2019. 

Hobeichi, S., Abramowitz, G., and Evans, J. P.: Robust historical evapotranspiration trends across climate regimes, Hydrol. Earth Syst. Sci., 25, 3855–3874,, 2021a. 

Hobeichi, S., Abramowitz, G., and Evans, J. P.: Derived Optimal Linear Combination Evapotranspiration – DOLCE v2.1, NCI National Research Data Collection [data set],, 2021b. 

Humphrey, V., Zscheischler, J., Ciais, P., Gudmundsson, L., Sitch, S., and Seneviratne, S. I.: Sensitivity of atmospheric CO2 growth rate to observed changes in terrestrial water storage, Nature, 560, 628–631,, 2018. 

Ibáñez, I., Acharya, K., Juno, E., Karounos, C., Lee, B. R., McCollum, C., Schaffer-Morrison, S., and Tourville, J.: Forest resilience under global environmental change: Do we have the information we need? A systematic review, edited by: Zang, R., PLoS One, 14, e0222207,, 2019. 

Jackson, R. B., Canadell, J., Ehleringer, J. R., Mooney, H. A., Sala, O. E., and Schulze, E. D.: A global analysis of root distributions for terrestrial biomes, Oecologia, 108, 389–411,, 1996. 

Jiang, X., Niu, G.-Y., and Yang, Z.-L.: Impacts of vegetation and groundwater dynamics on warm season precipitation over the Central United States, J. Geophys. Res., 114, D06109,, 2009. 

Jones, D., Wang, W., and Fawcett, R.: High-quality spatial climate data-sets for Australia, Aust. Meteorol. Oceanogr. J., 58, 233–248,, 2009. 

Jyoteeshkumar reddy, P., Sharples, J. J., Lewis, S. C., and Perkins-Kirkpatrick, S. E.: Modulating influence of drought on the synergy between heatwaves and dead fine fuel moisture content of bushfire fuels in the Southeast Australian region, Weather Clim. Extrem., 31, 100300,, 2021. 

Kattge, J. and Knorr, W.: Temperature acclimation in a biochemical model of photosynthesis: a reanalysis of data from 36 species, Plant. Cell Environ., 30, 1176–1190,, 2007. 

Keune, J., Gasper, F., Goergen, K., Hense, A., Shrestha, P., Sulis, M., and Kollet, S.: Studying the influence of groundwater representations on land surface-atmosphere feedbacks during the European heat wave in 2003, J. Geophys. Res.-Atmos., 121, 13301–13325,, 2016. 

Kim, W., Iizumi, T., and Nishimori, M.: Global Patterns of Crop Production Losses Associated with Droughts from 1983 to 2009, J. Appl. Meteorol. Clim., 58, 1233–1244,, 2019. 

Kowalczyk, E. A., Wang, Y. P., and Law, R. M.: The CSIRO Atmosphere Biosphere Land Exchange (CABLE) model for use in climate models and as an offline model, CSIRO Mar. Atmos. Res. Pap., 13, 1–42,, 2006. 

Kuginis, L., Dabovic, J., Burne, G., Raine, A., and Hemakumara, H.: Methods for the identification of high probability groundwater dependent vegetation ecosystems, DPI Water: Sydney, NSW, available at: (last access: 29 August 2021), 2016. 

Kumarathunge, D. P., Medlyn, B. E., Drake, J. E., Tjoelker, M. G., Aspinwall, M. J., Battaglia, M., Cano, F. J., Carter, K. R., Cavaleri, M. A., Cernusak, L. A., Chambers, J. Q., Crous, K. Y., De Kauwe, M. G., Dillaway, D. N., Dreyer, E., Ellsworth, D. S., Ghannoum, O., Han, Q., Hikosaka, K., Jensen, A. M., Kelly, J. W. G. G., Kruger, E. L., Mercado, L. M., Onoda, Y., Reich, P. B., Rogers, A., Slot, M., Smith, N. G., Tarvainen, L., Tissue, D. T., Togashi, H. F., Tribuzy, E. S., Uddling, J., Vårhammar, A., Wallin, G., Warren, J. M., and Way, D. A.: Acclimation and adaptation components of the temperature dependence of plant photosynthesis at the global scale, New Phytol., 222, 768–784,, 2019. 

Landerer, F. W., Flechtner, F. M., Save, H., Webb, F. H., Bandikova, T., Bertiger, W. I., Bettadpur, S. V., Byun, S. H., Dahle, C., Dobslaw, H., Fahnestock, E., Harvey, N., Kang, Z., Kruizinga, G. L. H., Loomis, B. D., McCullough, C., Murböck, M., Nagel, P., Paik, M., Pie, N., Poole, S., Strekalov, D., Tamisiea, M. E., Wang, F., Watkins, M. M., Wen, H., Wiese, D. N., and Yuan, D.: Extending the global mass change data record: GRACE Follow-On instrument and science data performance, Geophys. Res. Lett., 47, e2020GL088306,, 2020. 

Leblanc, M. J., Tregoning, P., Ramillien, G., Tweed, S. O., and Fakes, A.: Basin-scale, integrated observations of the early 21st century multiyear drought in southeast Australia, Water Resour. Res., 45, 1–10,, 2009. 

Lesk, C., Rowhani, P., and Ramankutty, N.: Influence of extreme weather disasters on global crop production, Nature, 529, 84–87,, 2016. 

Leuning, R., Kelliher, F. M., Pury, D. G. G., and Schulze, E.-D.: Leaf nitrogen, photosynthesis, conductance and transpiration: scaling from leaves to canopies, Plant, Cell Environ., 18, 1183–1200,, 1995. 

LP DAAC: The Application for Extracting and Exploring Analysis Ready Samples (AρρEEARS), available at:, last access: 1 September 2021. 

Marchionni, V., Daly, E., Manoli, G., Tapper, N. J., Walker, J. P., and Fatichi, S.: Groundwater buffers drought effects and climate variability in urban reserves, Water Resour. Res., 56, e2019WR026192,, 2020. 

Martens, B., Miralles, D. G., Lievens, H., van der Schalie, R., de Jeu, R. A. M., Fernández-Prieto, D., Beck, H. E., Dorigo, W. A., and Verhoest, N. E. C.: GLEAM v3: satellite-based land evaporation and root-zone soil moisture, Geosci. Model Dev., 10, 1903–1925,, 2017. 

Martínez-de la Torre, A. and Miguez-Macho, G.: Groundwater influence on soil moisture memory and land–atmosphere fluxes in the Iberian Peninsula, Hydrol. Earth Syst. Sci., 23, 4909–4932,, 2019. 

Martinez, J. A., Dominguez, F., and Miguez-Macho, G.: Effects of a groundwater scheme on the simulation of soil moisture and evapotranspiration over southern South America, J. Hydrometeorol., 17, 2941–2957,, 2016a. 

Martinez, J. A., Dominguez, F., and Miguez-Macho, G.: Impacts of a groundwater scheme on hydroclimatological conditions over southern South America, J. Hydrometeorol., 17, 2959–2978,, 2016b. 

Massmann, A., Gentine, P., and Lin, C.: When does vapor pressure deficit drive or reduce evapotranspiration?, J. Adv. Model. Earth Sy., 11, 3305–3320,, 2019. 

Maxwell, R. M., Chow, F. K., and Kollet, S. J.: The groundwater–land-surface–atmosphere connection: Soil moisture effects on the atmospheric boundary layer in fully-coupled simulations, Adv. Water Resour., 30, 2447–2466,, 2007. 

Maxwell, R. M., Lundquist, J. K., Mirocha, J. D., Smith, S. G., Woodward, C. S., and Tompson, A. F. B.: Development of a coupled groundwater–atmosphere model, Mon. Weather Rev., 139, 96–116,, 2011. 

McVicar, T. R., Van Niel, T. G., Li, L. T., Roderick, M. L., Rayner, D. P., Ricciardulli, L., and Donohue, R. J.: Wind speed climatology and trends for Australia, 1975–2006: Capturing the stilling phenomenon and comparison with near-surface reanalysis output, Geophys. Res. Lett., 35, 1–6,, 2008. 

Medlyn, B. E., Duursma, R. A., Eamus, D., Ellsworth, D. S., Prentice, C. I., Barton, C. V. M., Crous, K. Y., De Angelis, P., Freeman, M., and Wingate, L.: Reconciling the optimal and empirical approaches to modelling stomatal conductance, Glob. Chang. Biol., 17, 2134–2144,, 2011. 

Mercado, L. M., Medlyn, B. E., Huntingford, C., Oliver, R. J., Clark, D. B., Sitch, S., Zelazowski, P., Kattge, J., Harper, A. B., and Cox, P. M.: Large sensitivity in land carbon storage due to geographical and temporal variation in the thermal response of photosynthetic capacity, New Phytol., 218, 1462–1477,, 2018. 

Miller, G. R., Chen, X., Rubin, Y., Ma, S., and Baldocchi, D. D.: Groundwater uptake by woody vegetation in a semiarid oak savanna, Water Resour. Res., 46, 2009WR008902,, 2010. 

Miralles, D. G., Holmes, T. R. H., De Jeu, R. A. M., Gash, J. H., Meesters, A. G. C. A., and Dolman, A. J.: Global land-surface evaporation estimated from satellite-based observations, Hydrol. Earth Syst. Sci., 15, 453–469,, 2011. 

Miralles, D. G., Teuling, A. J., van Heerwaarden, C. C., and Vilà-Guerau de Arellano, J.: Mega-heatwave temperatures due to combined soil desiccation and atmospheric heat accumulation, Nat. Geosci., 7, 345–349,, 2014. 

Miralles, D. G., Gentine, P., Seneviratne, S. I., and Teuling, A. J.: Land-atmospheric feedbacks during droughts and heatwaves: state of the science and current challenges, Ann. N. Y. Acad. Sci., 1436, 19–35,, 2019. 

Mitchell, P. J., O'Grady, A. P., Hayes, K. R., and Pinkard, E. A.: Exposure of trees to drought-induced die-off is defined by a common climatic threshold across different vegetation types, Ecol. Evol., 4, 1088–1101,, 2014. 

Mu, M.: bibivking/Groundwater_Vegetation_Heatwave_Drought: Groundwater Vegetation Heatwave Drought (v1.0), Zenodo [code],, 2021. 

Mu, M., De Kauwe, M. G., Ukkola, A. M., Pitman, A. J., Gimeno, T. E., Medlyn, B. E., Or, D., Yang, J., and Ellsworth, D. S.: Evaluating a land surface model at a water-limited site: implications for land surface contributions to droughts and heatwaves, Hydrol. Earth Syst. Sci., 25, 447–471,, 2021. 

Mukherjee, S., Ashfaq, M. and Mishra, A. K.: Compound Drought and Heatwaves at a Global Scale: The role of natural climate variability-associated synoptic patterns and land-surface energy budget anomalies, J. Geophys. Res. Atmos., 125, e2019JD031943,, 2020. 

Nairn, J. R. and Fawcett, R. J. B.: The excess heat factor: A metric for heatwave intensity and its use in classifying heatwave severity, Int. J. Environ. Res. Pu., 12, 227–253,, 2014. 

National Climate Centre: The exceptional January-February 2009 heatwave in south-eastern Australia, available at: (last access: 29 August 2021), 2009. 

NCI: CABLE: The Community Atmosphere Biosphere Land Exchange Model, available at:, last access: 4 August 2021. 

Orth, R. and Destouni, G.: Drought reduces blue-water fluxes more strongly than green-water fluxes in Europe, Nat. Commun., 9, 3602,, 2018. 

O’sullivan, O. S., Heskel, M. A., Reich, P. B., Tjoelker, M. G., Weerasinghe, L. K., Penillard, A., Zhu, L., Egerton, J. J. G., Bloomfield, K. J., Creek, D., Bahar, N. H. A., Griffin, K. L., Hurry, V., Meir, P., Turnbull, M. H., and Atkin, O. K.: Thermal limits of leaf metabolism across biomes, Glob. Change Biol., 23, 209–223,, 2017. 

Perkins-Kirkpatrick, S. E., White, C. J., Alexander, L. V, Argüeso, D., Boschat, G., Cowan, T., Evans, J. P., Ekström, M., Oliver, E. C. J., Phatak, A., and Purich, A.: Natural hazards in Australia: heatwaves, Clim. Change, 139, 101–114,, 2016. 

Perkins, S. E.: A review on the scientific understanding of heatwaves – Their measurement, driving mechanisms, and changes at the global scale, Atmos. Res., 164–165, 242–267,, 2015. 

Raupach, M. R.: Simplified expressions for vegetation roughness length and zero-plane displacement as functions of canopy height and area index, Bound.-Lay. Meteorol., 71, 211–216,, 1994. 

Raupach, M. R., Finkele, K., and Zhang, L.: SCAM: a soil-canopy-atmosphere model: description and comparisons with field data, CSIRO Centre for Environmental Mechanics, Canberra, Australia,, 1997. 

Raupach, M. R., Haverd, V., and Briggs, P. R.: Sensitivities of the Australian terrestrial water and carbon balances to climate change and variability, Agr. Forest Meteorol., 182–183, 277–291,, 2013. 

Reich, P. B., Sendall, K. M., Rice, K., Rich, R. L., Stefanski, A., Hobbie, S. E., and Montgomery, R. A.: Geographic range predicts photosynthetic and growth response to warming in co-occurring tree species, Nat. Clim. Change, 5, 148–152,, 2015. 

Richards, J. H. and Caldwell, M. M.: Hydraulic lift: Substantial nocturnal water transport between soil layers by Artemisia tridentata roots, Oecologia, 73, 486–489,, 1987. 

Ruehr, N. K., Grote, R., Mayr, S., and Arneth, A.: Beyond the extreme: recovery of carbon and water relations in woody plants following heat and drought stress, Tree Physiol., 39, 1285–1299,, 2019. 

Sandi, S. G., Rodriguez, J. F., Saintilan, N., Wen, L., Kuczera, G., Riccardi, G., and Saco, P. M.: Resilience to drought of dryland wetlands threatened by climate change, Sci. Rep.-UK, 10, 13232,, 2020. 

Schenk, H. J. and Jackson, R. B.: The global biogeography of roots, Ecol. Monogr., 72, 311–328, doi:[0311:TGBOR]2.0.CO;2, 2002. 

Schumacher, D. L., Keune, J., van Heerwaarden, C. C., Vilà-Guerau de Arellano, J., Teuling, A. J., and Miralles, D. G.: Amplification of mega-heatwaves through heat torrents fuelled by upwind drought, Nat. Geosci., 12, 712–717,, 2019. 

Seneviratne, S. I., Corti, T., Davin, E. L., Hirschi, M., Jaeger, E. B., Lehner, I., Orlowsky, B., and Teuling, A. J.: Investigating soil moisture–climate interactions in a changing climate: A review, Earth-Sci. Rev., 99, 125–161,, 2010. 

Shrestha, P., Sulis, M., Masbou, M., Kollet, S., and Simmer, C.: A scale-consistent terrestrial systems modeling platform based on COSMO, CLM, and ParFlow, Mon. Weather Rev., 142, 3466–3483,, 2014. 

Smith, N. G. and Dukes, J. S.: Plant respiration and photosynthesis in global-scale models: incorporating acclimation to temperature and CO2, Glob. Change Biol., 19, 45–63,, 2013. 

Smith, N. G., Malyshev, S. L., Shevliakova, E., Kattge, J., and Dukes, J. S.: Foliar temperature acclimation reduces simulated carbon sensitivity to climate, Nat. Clim. Change, 6, 407–411,, 2016. 

Swenson, S. C. and Lawrence, D. M.: Assessing a dry surface layer-based soil resistance parameterization for the Community Land Model using GRACE and FLUXNET-MTE data, J. Geophys. Res.-Atmos., 119, 10299–10312,, 2014. 

Teuling, A. J., Van Loon, A. F., Seneviratne, S. I., Lehner, I., Aubinet, M., Heinesch, B., Bernhofer, C., Grünwald, T., Prasse, H., and Spank, U.: Evapotranspiration amplifies European summer drought, Geophys. Res. Lett., 40, 2071–2075,, 2013. 

Thorburn, P. J., Walker, G. R., and Woods, P. H.: Comparison of diffuse discharge from shallow water tables in soils and salt flats, J. Hydrol., 136, 253–274,, 1992. 

Trudinger, C. M., Haverd, V., Briggs, P. R., and Canadell, J. G.: Interannual variability in Australia's terrestrial carbon cycle constrained by multiple observation types, Biogeosciences, 13, 6363–6383,, 2016. 

Trugman, A. T., Medvigy, D., Mankin, J. S., and Anderegg, W. R. L. L.: Soil moisture stress as a major driver of carbon cycle uncertainty, Geophys. Res. Lett., 45, 6495–6503,, 2018. 

Ukkola, A. M., De Kauwe, M. G., Pitman, A. J., Best, M. J., Abramowitz, G., Haverd, V., Decker, M., and Haughton, N.: Land surface models systematically overestimate the intensity, duration and magnitude of seasonal-scale evaporative droughts, Environ. Res. Lett., 11, 104012,, 2016a. 

Ukkola, A. M., Pitman, A. J., Decker, M., De Kauwe, M. G., Abramowitz, G., Kala, J., and Wang, Y.-P.: Modelling evapotranspiration during precipitation deficits: identifying critical processes in a land surface model, Hydrol. Earth Syst. Sci., 20, 2403–2419,, 2016b. 

Ukkola, A. M., De Kauwe, M. G., Roderick, M. L., Abramowitz, G. and Pitman, A. J.: Robust future changes in meteorological drought in CMIP6 projections despite uncertainty in precipitation, Geophys. Res. Lett., 47, e2020GL087820,, 2020. 

van Dijk, A. I. J. M. J. M., Beck, H. E., Crosbie, R. S., de Jeu, R. A. M., Liu, Y. Y., Podger, G. M., Timbal, B., and Viney, N. R.: The Millennium Drought in southeast Australia (2001–2009): Natural and human causes and implications for water resources, ecosystems, economy, and society, Water Resour. Res., 49, 1040–1057,, 2013. 

Verdon-Kidd, D. C. and Kiem, A. S.: Nature and causes of protracted droughts in southeast Australia: Comparison between the Federation, WWII, and Big Dry droughts, Geophys. Res. Lett., 36, L22707,, 2009. 

Wada, Y.: Impacts of groundwater pumping on regional and global water resources, in: Terrestrial Water Cycle and Climate Change, American Geophysical Union, edited by: Tang, Q. and Oki, T., Hoboken, United States, 71–101,, 2016. 

Wan, Z., Hook, S., and Hulley, G.: MOD11A1 MODIS/Terra Land Surface Temperature/Emissivity Daily L3 Global 1km SIN Grid V006, NASA EOSDIS Land Processes DAAC [data set],, 2015a 

Wan, Z., Hook, S., and Hulley, G.: MYD11A1 MODIS/Aqua Land Surface Temperature/Emissivity Daily L3 Global 1km SIN Grid V006, NASA EOSDIS Land Processes DAAC [data set],, 2015b. 

Wan, Z. and Li, Z.-L.: A physics-based algorithm for retrieving land-surface emissivity and temperature from EOS/MODIS data, IEEE T. Geosci. Remote, 35, 980–996,, 1997. 

Wang, K. and Dickinson, R. E.: A review of global terrestrial evapotranspiration: Observation, modeling, climatology, and climatic variability, Rev. Geophys., 50, RG2005,, 2012. 

Wang, P., Niu, G., Fang, Y., Wu, R., Yu, J., Yuan, G., Pozdniakov, S. P., and Scott, R. L.: Implementing dynamic root optimization in Noah-MP for simulating phreatophytic root water uptake, Water Resour. Res., 54, 1560–1575,, 2018. 

Wang, Y.-P. and Leuning, R.: A two-leaf model for canopy conductance, photosynthesis and partitioning of available energy I: Model description and comparison with a multi-layered model, Agr. Forest Meteorol., 91, 89–111,, 1998. 

Wang, Y. P., Kowalczyk, E., Leuning, R., Abramowitz, G., Raupach, M. R., Pak, B., van Gorsel, E., and Luhar, A.: Diagnosing errors in a land surface model (CABLE) in the time and frequency domains, J. Geophys. Res., 116, G01034,, 2011. 

Warren, J. M., Hanson, P. J., Iversen, C. M., Kumar, J., Walker, A. P., and Wullschleger, S. D.: Root structural and functional dynamics in terrestrial biosphere models – evaluation and recommendations, New Phytol., 205, 59–78,, 2015. 

Watkins, M. M., Wiese, D. N., Yuan, D.-N., Boening, C., and Landerer, F. W.: Improved methods for observing Earth's time variable mass distribution with GRACE using spherical cap mascons, J. Geophys. Res.-Sol. Ea., 120, 2648–2671,, 2015. 

Wiese, D. N., Landerer, F. W., and Watkins, M. M.: Quantifying and reducing leakage errors in the JPL RL05M GRACE mascon solution, Water Resour. Res., 52, 7490–7502,, 2016. 

Wiese, D. N., Yuan, D.-N., Boening, C., Landerer, F. W., and Watkins, M. M.: JPL GRACE Mascon Ocean, Ice, and Hydrology Equivalent Water Height Release 06 Coastal Resolution Improvement (CRI) Filtered Version 1.0. Ver. 1.0., PO.DAAC, CA, USA,, 2018. 

Wiese, D. N., Yuan, D.-N., Boening, C., Landerer, F. W., and Watkins, M. M.: JPL GRACE Mascon Ocean, Ice, and Hydrology Equivalent Water Height RL06 CRI Filtered Version 02, PO.DAAC [data set],, 2019. 

Yang, J., Duursma, R. A., De Kauwe, M. G., Kumarathunge, D., Jiang, M., Mahmud, K., Gimeno, T. E., Crous, K. Y., Ellsworth, D. S., Peters, J., Choat, B., Eamus, D., and Medlyn, B. E.: Incorporating non-stomatal limitation improves the performance of leaf and canopy models at high vapour pressure deficit, Tree Physiol., 39, 1961–1974,, 2019.  

Zencich, S. J., Froend, R. H., Turner, J. V., and Gailitis, V.: Influence of groundwater depth on the seasonal sources of water accessed by Banksia tree species on a shallow, sandy coastal aquifer, Oecologia, 131, 8–19,, 2002. 

Zeng, X. and Decker, M.: Improving the numerical solution of soil moisture–based Richards equation for land models with a deep or shallow water table, J. Hydrometeorol., 10, 308–319,, 2009. 

Zhang, H., Pak, B., Wang, Y. P., Zhou, X., Zhang, Y., and Zhang, L.: Evaluating surface water cycle simulated by the Australian Community Land Surface Model (CABLE) across Different Spatial and Temporal Domains, J. Hydrometeorol., 14, 1119–1138,, 2013. 

Zhou, S., Williams, A. P., Berg, A. M., Cook, B. I., Zhang, Y., Hagemann, S., Lorenz, R., Seneviratne, S. I., and Gentine, P.: Land–atmosphere feedbacks exacerbate concurrent soil drought and atmospheric aridity, P. Natl. Acad. Sci. USA, 116, 18848–18853,, 2019. 

Zipper, S. C., Keune, J., and Kollet, S. J.: Land use change impacts on European heat and drought: remote land-atmosphere feedbacks mitigated locally by shallow groundwater, Environ. Res. Lett., 14, 044012,, 2019. 

Short summary
Groundwater can buffer the impacts of drought and heatwaves on ecosystems, which is often neglected in model studies. Using a land surface model with groundwater, we explained how groundwater sustains transpiration and eases heat pressure on plants in heatwaves during multi-year droughts. Our results showed the groundwater’s influences diminish as drought extends and are regulated by plant physiology. We suggest neglecting groundwater in models may overstate projected future heatwave intensity.
Final-revised paper