Numerical behavior of a coupled surface/subsurface, flow/transport hydrological model

Integrated surface/subsurface hydrological models (ISSHMs) are by now widely used in research and applied hydrology. While most studies have so far focused on water flow alone, ISSHMs that include also solute transport are beginning to get attention (e.g., Scudeler et al., 2016, doi:10.5194/hess-20-4061-2016; Gatel et al., 2019, doi:10.1016/j.envsoft.2018.12.006). Numerous numerical challenges are associated with these "doubly coupled" systems: correct treatment of surface boundary conditions and other mass and flux exchange terms; appropriate time stepping schemes across subsystems that are characterized by different dynamic time scales and often also widely different numerical discretization approaches; performance assessments that can be highly sensitive to the response variables of interest; and so on. We will illustrate some of these challenges via test case simulations of an experimental hillslope using the CATHY (CATchment HYdrology) model (Camporese et al., 2010, doi:10.1029/2008WR007536; Weill et al., 2011, doi:10.1016/j.advwatres.2010.10.001). The boundary condition-based coupling strategy used in this model (Putti and Paniconi, 2004, doi:10.1016/S0167-5648(04)80152-7) has been shown to be mathematically rigorous and mass-conservative for the flow model (Sochala et al., 2009, doi:10.1016/j.cma.2009.02.024). The convergence-based time step adaptation strategy used for the nonlinear flow equation (Paniconi and Putti, 1994, doi:10.1029/94WR02046) is likewise thoroughly tested (e.g., D'Haese et al., 2007, doi:10.1002/fld.1369) and widely used. Nonetheless, these schemes, and analogous approaches used in other ISSHMs, need to be adapted and thoroughly tested for coupled systems that include solute transport.