This study explores the impact of deep ( >5 m) critical zone (CZ) architecture on vegetation distribution in a semi-arid snow-dominated climate. Utilizing seismic refraction surveys, we identified a significant correlation between saprolite thickness and LiDAR-derived canopy heights (R²=0.47). We argue that CZ structure, specifically shallow fractured bedrock under valley bottoms, redirects groundwater to locations where trees are established—suggesting they are located in specific locations with access to nutrients and water. This work provides a unique spatially exhaustive perspective and adds to growing evidence that in addition to other factors such as slope, aspect, and climate, deep CZ structure plays a vital role in ecosystem development and resilience.

Topography is a key control on runoff generation, as topographic slope affects hydraulic gradients and curvature affects water flow paths. At the same time, runoff generation shapes topography through erosion, which affects landscape morphology over long timescales. Previous modeling efforts suggest that subsurface hydrological properties, relative to climate, are key mediators of this relationship. Specifically, when subsurface transmissivity and water storage capacity are low, (1) saturated areas and storm runoff should be larger and more variable, and (2) hillslopes shorter and with less relief, assuming other geomorphic factors are held constant. While these patterns appear in simulations, it remains uncertain whether subsurface properties can exert such a strong control on emergent properties in the field. We compared emergent hydrological function and topography in two watersheds that have very similar climatic and geologic history, but very different subsurface properties due to contrasting bedrock lithology. We found that hillslopes were systematically shorter and saturated areas more dynamic at the site with lower transmissivity. To confirm that these differences were due to subsurface hydrology rather than differences in geomorphic process rates, we estimated all parameters of a coupled groundwater-landscape evolution model without calibration. We showed that the difference in subsurface properties has a profound effect on topography and hydrological function that cannot be explained by differences in geomorphic process rates alone. The comparison to field data also exposed model limitations, which we discuss in the context of future efforts to understand the role of hydrology in the long-term evolution of Earth’s critical zone.

Features of landscape morphology—including slope, curvature, and drainage dissection—are important controls on runoff generation in upland landscapes. Over long timescales, runoff plays an essential role in shaping these same features through surface erosion. This feedback between erosion and runoff generation suggests that modeling long-term landscape evolution together with dynamic runoff generation could provide insight into hydrological function. Here we examine the emergence of variable source area runoff generation in a new coupled hydro-geomorphic model that accounts for water balance partitioning between surface flow, subsurface flow, and evapotranspiration as landscapes evolve over millions of years. We derive a minimal set of dimensionless numbers that provide insight into how hydrologic and geomorphic parameters together affect landscapes. We find an inverse relationship between the dimensionless local relief and the fraction of the landscape that produces saturation excess overland flow, in agreement with the synthesis described in the “Dunne Diagram.’ Furthermore, we find an inverse, nonlinear relationship between the Hillslope number, which describes topographic relief relative to aquifer thickness, and the proportion of the landscape that variably saturated. Certain parameter combinations produce features with wide valley bottom wetlands and nondendritic, diamond-shaped drainage networks, which cannot be produced by simple landscape evolution models alone. With these results, we demonstrate the power of coupled hydrogeomorphic models for generating new insights into hydrological processes, and also suggest that subsurface hydrology may be integral for modeling aspects of long-term landscape evolution.

Unsteady transit time distribution (TTD) theory is a promising new approach for merging hydrologic and water quality models at the catchment scale. A major obstacle to widespread adoption of the theory, however, has been the specification of the StorAge Selection (SAS) function, which describes how the selection of water for outflow is biased by age. In this paper we hypothesize that some unsteady hydrologic systems of practical interest can be described, to first-order, by a “shifted-uniform” SAS that falls along a continuum between plug flow sampling (for which only the oldest water in storage is sampled for outflow) and uniform sampling (for which water in storage is sampled randomly for outflow). For this choice of SAS function, explicit formulae are derived for the evolving: (1) age distribution of water in storage; (2) age distribution of water in outflow; and (3) breakthrough concentration of a conservative solute under either continuous or impulsive addition. Model predictions conform closely to chloride and deuterium breakthrough curves measured previously in a sloping lysimeter subject to periodic wetting, although refinements of the model are needed to account for the reconfiguration of flow paths at high storage levels (the so-called inverse storage effect). The analytical results derived in this paper should lower the barrier to applying TTD theory in practice, ease the computational demands associated with simulating solute transport through complex hydrologic systems, open up new opportunities for real-time control, and provide physical insights that might not be apparent from traditional numerical solutions of the governing equations.

The hydrologic dynamics and geomorphic evolution of watersheds are intimately coupled – runoff generation and water storage are controlled by topography and properties of the surface and subsurface, while also affecting the evolution of those properties over geologic time. However, the large disparity between their timescales has made it difficult to examine interdependent controls on emergent hydro-geomorphic properties, such as hillslope length, drainage density, extent of surface saturation. In this study, we develop a new model coupling hydrology and landscape evolution to explore how runoff generation affects long-term catchment evolution, and analyze numerical results using a nondimensional scaling framework. We focus on hydrologic processes dominating in humid climates where storm runoff primarily arises from shallow subsurface flow and from precipitation on saturated areas. The model solves hydraulic groundwater equations to predict the water table location given prescribed, constant groundwater recharge. Water in excess of the subsurface capacity for transport becomes overland flow, which generates shear stress on the surface and may detach and transport sediment. This affects the landscape form that in turn affects runoff generation. We show that (1) three dimensionless parameters describe the possible steady state landscapes that coevolve under steady recharge; (2) hillslope length increases with increasing transmissivity relative to the recharge rate; (3) three topographic metrics—steepness index, Laplacian curvature, and topographic index—provide a basis to recover key model parameters from topography (including subsurface transmissivity). These results open possibilities for topographic analysis of humid upland landscapes that could inform quantitative understanding of hydrological processes at the landscape scale.

Spatially-integrated water transport dynamics at the hillslope scale have rarely been observed directly, and underlying physical mechanisms of those dynamics are poorly understood. We present time-variable transit time distributions (TTDs) and StorAge Selection (SAS) functions observed during a 28 days tracer experiment conducted at the Landscape Evolution Observatory (LEO), Biosphere 2, University of Arizona, AZ, USA. The observed form of the SAS functions is concave, meaning that older water in the hillslope was discharged more than younger water. The concavity is, in part, explained by the relative importance of advective and diffusive water dynamics and by the geomorphologic structure of the hillslopes. A simple numerical examination illustrates that, for straight plan shape hillslopes, the saturated zone SAS function is concave when hillslope Peclet (Pe) number is large. We also investigated the effect of hillslope planform geometry on the SAS function: The more convergent the plan shape is, the more concave the SAS function is. A numerical examination also indicates that the unsaturated zone SAS function is concave for straight and convergent hillslopes, when the soil thickness is constant. The concavity of those sub-component SAS functions signifies that the hillslope scale SAS function is concave for straight or convergent plan shape hillslopes when the hillslope Pe number is high.

Uncovering the hillslope scale flow and transport dynamics in an experimental hydrologic systemMinseok Kim1, Till H. M. Volkmann1,2, Aaron Bugaj1, Yadi Wang3, Antônio A. Meira Neto4, Katarena Matos4, Ciaran J. Harman5,6, Peter A. Troch1,41Biosphere 2, University of Arizona, Tucson, AZ, USA,2Applied Intelligence, Accenture, Kronberg im Taunus, Germany, 3Department of Environmental Science, University of Arizona, Tucson, AZ, USA, 4Department of Hydrology and Atmospheric Sciences, University of Arizona, Tucson, AZ, USA, 5Department of Environmental Health and Engineering, Johns Hopkins University, Baltimore, MD, USA,6Department of Earth and Planetary Sciences, Johns Hopkins University, Baltimore, MD, USAHillslope scale water flow and transport dynamics have been extensively studied (Burt & McDonnell, 2015; Hewlett & Hibbert, 1963), but observing those internal dynamics in high spatial and temporal resolutions remains challenging. In this study, we uncover internal water flow and transport dynamics in an artificial hillslope in the Landscape Evolution Observatory (LEO), Biosphere 2, University of Arizona, Tucson, USA, using the experimental dataset collected in December 2016. Complete information about the hillslope and experiment can be found elsewhere (Pangle et al., 2015; Till H. M. Volkmann et al., 2018); Here, we only summarize some relevant information.The first part of the animation describes the experimental system and setup (time 00:12 – 04:14 in Animation S1). The LEO hillslope is 330 m3 (30 m long, 11 m wide, and 1 m deep) sloping soil lysimeter. The hillslope is primarily made up of loamy-sand textured basaltic tephra, and the most downslope 5.5 m3 is filled with gravel-textured basaltic tephra. A custom irrigation system supplies reverse osmosis filtered water onto the LEO surface. The downslope boundary is exposed to atmospheric pressure, creating the seepage face boundary condition. The sensor networks (including pressure transducers and volumetric water content sensors) and the water isotope sampling locations and intervals (7 hrs to 101 hrs) are illustrated in Animation S1 (time 02:09 – 03:01). The isotope composition of subsurface water is obtained from laser-based online measurements of vapor that is extracted via custom gas probes through equilibrium calculation (T. H.M. Volkmann & Weiler, 2014). The irrigation sequence of this experiment was designed to generate a periodic steady state, which allows the application of the PERidoic Tracer Hierarchy method (Harman & Kim, 2014) for the observation of the time-variable transit time distributions and the StorAge Selection functions. Deuterium-labeled water was irrigated during the first two irrigation events.The second part of the animation shows the dynamics of the perched water table and soil water content (time 04:15 – 06:53). The extent of the saturated zone was estimated using the pressure transducer data and Delaunay triangulation (Delaunay, 1934). The experimental data show the saturation from below mechanisms—wetting up from the bedrock surface into the soil profile (McDonnell, 1997)—and the saturation from downslope to upslope. The water table profile forms a wedge-like shape, which is a characteristic of hillslope with a high hillslope (Peclet) number (Berne et al., 2005; Brutsaert, 1994). The hillslope Peclet number of the LEO hillslope during the experiment is high (> 10) (Kim et al., 2020). Significant time delays in the water table dynamics are observed at some upslope locations (e.g., at 13 m upslope), which is mostly due to the delayed water supply from the convergent upslope area. The water content data indicates that the convergent upslope water content began to decrease around the timing of the water table peak at 13 m upslope.The third part of the animation shows the tracer dynamics (from time 06:43). The animated experimental data reveal two notable water transport dynamics. First, the vertical tracer movement is faster at the upslope. This faster movement at the upslope is, in a sense, counter-intuitive because the upslope region is drier than the downslope. This is due to the lateral flow in the saturated zone and the tension saturated zone, that are thicker at the downslope. While water velocity is higher at the downslope, the direction of velocity is not vertical but rotated towards the downslope in those zones.Second, the animated data illustrate that old water is present only at the downslope. This observation is a characteristic of hillslope with a high hillslope number, in which old water is preferentially discharged (Kim et al., 2020). Indeed, the observed SAS function in this hillslope is concave (Kim et al., 2020), indicating that the hillslope preferentially discharges old water that is stored at the downslope.