Quantifying solute transport numerical dispersion in integrated surface-subsurface hydrological modeling

Numerical dispersion is a well-known problem that affects solute transport in groundwater simulations and can lead to wrong results, in terms of plume path overestimation and overprediction of contaminant dispersion. Numerical dispersion is generally introduced through stabilization techniques aimed at preventing oscillations, with the side effect of increasing mass spreading. Even though this issue has long been investigated in subsurface hydrology, little is known about its possible impacts on integrated surface–subsurface hydrological models (ISSHMs). In this study, we analyze numerical dispersion in the CATchment HYdrology (CATHY) model. In CATHY, a robust and computationally efficient time-splitting technique is implemented for the solution of the subsurface transport equation, whereby the advective part is solved on elements with an explicit finite volume scheme and the dispersive part is solved on nodes with an implicit finite element scheme. Taken alone, the advection and dispersion solvers provide accurate results. However, when coupled, the continuous transfer of concentration from elements to nodes, and vice versa, gives rise to a particular form of numerical dispersion. We assess the nature and impact of this artificial spreading through two sets of synthetic experiments. In the first set, the subsurface transport of a nonreactive tracer in two soil column test cases is simulated and compared with known analytical solutions. Different input dispersion coefficients and mesh discretizations are tested, in order to quantify the numerical error and define a criterion for its containment. In the second set of experiments, fully coupled surface–subsurface processes are simulated using two idealized hillslopes, one concave and one convex, and we examine how the additional subsurface dispersion affects the representation of pre-event water contribution to the streamflow hydrograph. Overall, we show that the numerical dispersion in CATHY that is caused by the transfer of information between elements and nodes can be kept under control if the grid Péclet number is less than 1. It is also suggested that the test cases used in this study can be useful benchmarks for integrated surface–subsurface hydrological models, for which thus far only flow benchmarks have been proposed.