A vector‐river network explicitly uses realistic geometries of river reaches and catchments for spatial discretization in a river model. This enables improving the accuracy of the physical properties of the modeled river system, compared to a gridded river network that has been used in Earth System Models. With a finer‐scale river network, resolving smaller‐scale river reaches, there is a need for efficient methods to route streamflow and its constituents throughout the river network. The purpose of this study is twofold: (1) develop a new method to decompose river networks into hydrologically independent tributary domains, where routing computations can be performed in parallel; and (2) perform global river routing simulations with two global river networks, with different scales, to examine the computational efficiency and the differences in discharge simulations at various temporal scales. The new parallelization method uses a hierarchical decomposition strategy, where each decomposed tributary is further decomposed into many sub‐tributary domains, enabling hybrid parallel computing. This parallelization scheme has excellent computational scaling for the global domain where it is straightforward to distribute computations across many independent river basins. However, parallel computing for a single large basin remains challenging. The global routing experiments show that the scale of the vector‐river network has less impact on the discharge simulations than the runoff input that is generated by the combination of land surface model and meteorological forcing. The scale of vector‐river networks needs to consider the scale of local hydrologic features such as lakes that are to be resolved in the network.

Earth System Models (ESMs) are essential tools for understanding and predicting global change, but they cannot explicitly resolve hillslope‐scale terrain structures that fundamentally organize water, energy, and biogeochemical stores and fluxes at subgrid scales. Here we bring together hydrologists, Critical Zone scientists, and ESM developers, to explore how hillslope structures may modulate ESM grid‐level water, energy, and biogeochemical fluxes. In contrast to the one‐dimensional (1‐D), 2‐ to 3‐m deep, and free‐draining soil hydrology in most ESM land models, we hypothesize that 3‐D, lateral ridge‐to‐valley flow through shallow and deep paths and insolation contrasts between sunny and shady slopes are the top two globally quantifiable organizers of water and energy (and vegetation) within an ESM grid cell. We hypothesize that these two processes are likely to impact ESM predictions where (and when) water and/or energy are limiting. We further hypothesize that, if implemented in ESM land models, these processes will increase simulated continental water storage and residence time, buffering terrestrial ecosystems against seasonal and interannual droughts. We explore efficient ways to capture these mechanisms in ESMs and identify critical knowledge gaps preventing us from scaling up hillslope to global processes. One such gap is our extremely limited knowledge of the subsurface, where water is stored (supporting vegetation) and released to stream baseflow (supporting aquatic ecosystems). We conclude with a set of organizing hypotheses and a call for global syntheses activities and model experiments to assess the impact of hillslope hydrology on global change predictions.