Integrated hydrologic models can simulate coupled surface and subsurface processes but are computationally expensive to run at high resolutions over large domains. Here we develop a novel deep learning model to emulate continental-scale subsurface flows simulated by the integrated ParFlow-CLM model. We compare convolutional neural networks like ResNet and UNet run autoregressively against our novel architecture called the Forced SpatioTemporal RNN (FSTR). The FSTR model incorporates separate encoding of initial conditions, static parameters, and meteorological forcings, which are fused in a recurrent loop to produce spatiotemporal predictions of groundwater. We evaluate the model architectures on their ability to reproduce 4D pressure heads, water table depths, and surface soil moisture over the contiguous US at 1km resolution and daily time steps over the course of a full water year. The FSTR model shows superior performance to the baseline models, producing stable simulations that capture both seasonal and event-scale dynamics across a wide array of hydroclimatic regimes. The emulators provide over 1000x speedup compared to the original physical model, which will enable new capabilities like uncertainty quantification and data assimilation for integrated hydrologic modeling that were not previously possible. Our results demonstrate the promise of using specialized deep learning architectures like FSTR for emulating complex process-based models without sacrificing fidelity.

Water age is a fundamental descriptor of source, storage, and mixing of water parcels in a watershed. The Lagrangian, particle tracking, approach is a powerful tool for physically-based modeling of water age distributions, but its application has been hampered since it is computationally demanding. In this study, we present a parallel approach for particle tracking simulations. This approach uses multi-GPU with MPI parallelism based on domain decomposition. An inherent challenge of distributed parallelization of Lagrangian approaches is the disparity in computational work or load imbalance (LIB) among different processing elements (PEs). Here, load balancing (LB) schemes were proposed to dynamically balance the distribution of particles across PEs during runtime. In the followed hillslope simulations, LIB was observed in all LB-disabled runs, e.g., with a load ratio of 423.62% by using 2-GPU in LW_Shrub case. LB schemes then accurately balanced the load distribution and improved the parallel scaling. Additionally, the parallel approach showed excellent overall speedup: a 60-fold improvement using 4-GPU relative to the serial run. A regional scale application further demonstrated the LB performance. The parallel time used by 8-GPU without LB was 31.33% reduced after LB was activated. When increasing 8-GPU with LB to 16-GPU with LB, it showed parallel scalability by reducing the parallel time of ~50%. This work shows how massively parallel computing can be applied to particle tracking in water age simulations. It also demonstrates the practical importance of load balancing in this context, which enables the large-scale simulations with an increased complexity of flow paths.

Stable isotopes of water are important tracers in hydrologic research for understanding water partitioning between vegetation, groundwater, and runoff, but are rarely applied to large watersheds with persistent snowpack and complex topography. We combined an extensive isotope dataset with a coupled hydrologic and snow isotope fractionation model to assess mechanisms of isotopic inputs into the soil zone and implications on recharge dynamics within a large, snow-dominated watershed of the Upper Colorado River Basin. Results indicate seasonal isotopic variability and isotope lapse rates of net precipitation are the dominant control on isotopic inputs to the basin. Snowpack fractionation processes account for <5% annual isotope influx variability. Isotopic fractionation processes are most important in the shrub-dominated upper montane. Effects of isotopic fractionation are less important in the low-density conifer forests of the upper subalpine due to vegetative shading, low aridity, and a deep, persistent snowpack that buffers small sublimation losses. Melt fractionation can have sub-seasonal effects on snowmelt isotope ratios with initial snowmelt depleted but later snowmelt relatively enriched in heavy isotopes through the isotopic mass balance of the remaining snowpack, with the efficiency of isotopic exchange between ice and liquid water declining as snow ablation progresses. Hydrologic analysis indicates maximum recharge in the upper subalpine with wet years producing more isotopically depleted snowmelt (1-2‰ reduction in d18O) through reduced aridity when energy-limited. The five-year volume-weighted d18O in this zone (18.2±0.4‰) matches groundwater observations from multiple deep wells, providing evidence that the upper subalpine is a preferential recharge zone in mountain systems.

A coupled hydrologic and snowpack stable water isotope model assesses controls on isotopic inputs across a large, mountainous basin. The most depleted isotope conditions occur in the upper subalpine where snow accumulation is high and rainfall is low. Snowmelt evolution over time indicates isotopic enrichment is not dictated by melt fractionation but is determined by elevation which controls the amount, phase and isotopic mass of spring precipitation coincident with the ablation period. With respect to snowpack kinetic fractionation, its effect on snowmelt is a balance between energy and snow-availability. It is highest above treeline and in the shrub-dominated upper montane where vegetation shading is low, while deep snowpack and conifer forests limit the influence of kinetic fractionation in the subalpine. Wet years reduce the effects of snowpack fraction on snowmelt across the basin, except in the lower montane where added snowfall bolsters water-limited conditions.

Despite the importance of headwater catchments for western United States’ water supply, these regions are often poorly understood, particularly with respect to quantitative understanding of evapotranspiration (ET) fluxes. Heterogeneity of land cover, physiography, and atmospheric patterns in these high-elevation regions lead to difficulty in developing spatially-distributed characterization of ET. As the largest terrestrial water flux behind precipitation, ET represents a significant fraction of the water budget for any watershed. Likewise, groundwater is the largest available freshwater store and has been shown to play a large role in the water balance, even in headwater systems. Using an eddy covariance tower in the East River Catchment, a Colorado River headwaters basin, this study estimates water and energy fluxes in high-elevation, complex systems to better constrain ET estimates and calculate overall water and energy budgets, including losses from groundwater. The eddy covariance method is used to estimate ET from years 2017 through 2019 at a saturated, riparian end-member site. Owing to complexities in near surface atmospheric structure such as stable boundary layers over snowpack and shallow terrain driven flow from surrounding landscape features, energy flux and ET estimates were limited to the warm season when energy closure residuals from the eddy-covariance system were reliably less than 30 %, a threshold commonly used in eddy covariance energy flux estimation. The resulting ET estimations are useful for constraining water budget estimates at this energy-limited site, which uses groundwater for up to 84 % of ET in the summer months. We also compared East River ET magnitudes and seasonality to two other flux towers (Niwot Ridge, CO and Valles Caldera, NM), located in the Rocky Mountains. This data is useful for constraining ET estimates in similar end-member locations across the East River Catchment. Our results show that groundwater-fed ET is a significant component of the water balance and groundwater may supply riparian ET even during low-snow years.

Despite the importance of headwater basins for western United States’ water supply, these regions are often poorly understood, particularly with respect to quantitative understanding of evapotranspiration (ET) fluxes. Heterogeneity of land cover, topography, and atmospheric patterns in these high-elevation regions lead to difficulty in developing spatially distributed characterization of ET. As a significant fraction of the water budget, ET contributes to overall water and energy availability in the basin. Using an eddy covariance tower in the East River Basin, a Colorado River headwaters basin, this study improves the quantification of water and energy fluxes in high-elevation, complex systems to better constrain ET estimates and calculate overall water and energy budgets. The eddy covariance method estimates ET from years 2017 through 2019 at a saturated, riparian end-member site. During the late spring, summer, and early fall months, due to strong variations in lower atmospheric stability and evidenced by a less than 30% energy balance closure error in these months (within the range of closure error reported at other riparian locations) we conclude that the eddy covariance method is useful in high-elevation, complex areas such as the East River Basin and helps bound regional ET estimates. We also compared East River ET magnitudes and seasonality to two other eddy covariance towers (Niwot Ridge, CO and Valles Caldera, NM), with similar site characteristics, located in the Rocky Mountains. East River ET estimations are useful for constraining water budget estimates at this energy-limited site, which uses groundwater for up to 76% of ET in the summer months. This data is useful for constraining ET estimates in similar end-member locations; however, to better constrain ET estimates across the entire East River basin, additional sampling is needed. This study helps constrain both the energy and water budgets in locations that are underrepresented by observations and where indirect estimates of ET may perform poorly.