Unsaturated subsurface flow with surface water and nonlinear in- and outflow conditions

We analytically and numerically analyze groundwater flow in a homogeneous soil described by the Richards equation, coupled to surface water represented by a set of ordinary differential equations (ODE's) on parts of the domain boundary, and with nonlinear outflow conditions of Signorini's type. The coupling of the partial differential equation (PDE) and the ODE's is given by nonlinear Robin boundary conditions. This article provides two major new contributions regarding these infiltration conditions. First, an existence result for the continuous coupled problem is established with the help of a regularization technique. Second, we analyze and validate a solver-friendly discretization of the coupled problem based on an implicit-explicit time discretization and on finite elements in space. The discretized PDE leads to convex spatial minimization problems which can be solved efficiently by monotone multigrid. Numerical experiments are provided using the DUNE numerics framework.Comment: 34 pages, 5 figure

Obstacles, Interfacial Forms, and Turbulence: A Numerical Analysis of Soil–Water Evaporation Across Different Interfaces

AbstractExchange processes between a turbulent free flow and a porous media flow are sensitive to the flow dynamics in both flow regimes, as well as to the interface that separates them. Resolving these complex exchange processes across irregular interfaces is key in understanding many natural and engineered systems. With soil–water evaporation as the natural application of interest, the coupled behavior and exchange between flow regimes are investigated numerically, considering a turbulent free flow as well as interfacial forms and obstacles. Interfacial forms and obstacles will alter the flow conditions at the interface, creating flow structures that either enhance or reduce exchange rates based on their velocity conditions and their mixing with the main flow. To evaluate how these interfacial forms change the exchange rates, interfacial conditions are isolated and investigated numerically. First, different flow speeds are compared for a flat surface. Second, a porous obstacle of varied height is introduced at the interface, and the effects the flow structures that develop have on the interface are analyzed. The flow parameters of this obstacle are then varied and the interfacial exchange rates investigated. Next, to evaluate the interaction of flow structures between obstacles, a second obstacle is introduced, separated by a varied distance. Finally, the shape of these obstacles is modified to create different wave forms. Each of these interfacial forms and obstacles is shown to create different flow structures adjacent to the surface which alter the mass, momentum, and energy conditions at the interface. These changes will enhance the exchange rate in locations where higher velocity gradients and more mixing with the main flow develop, but will reduce the exchange rate in locations where low velocity gradients and limited mixing with the main flow occur

Numerical modeling of flow and solute transport phenomena in subsurface and coupled surface-subsurface hydrology

The overall aim of the work described in this thesis is to bring a number of contributions to hydrology and hydrological modeling in the framework of a specific physically-based numerical model for integrated surface subsurface and flow-transport processes, the CATchment-HYdrology Flow-Transport (CATHY_FT) model. These contributions revolve around three main themes: the enhancement of the numerical performance of hydrological models for flow and transport phenomena, the improvement of our current understanding of complex boundary conditions in order to reduce the errors associated with their modeling, and the testing and benchmarking of distributed physically-based models for groundwater flow and transport processes. The work to achieve the general objective is elaborated into four stages. First, the Larson-Niklasson post-processing algorithm is implemented in CATHY_FT to reconstruct mass-conservative velocities from a linear, or P1, Galerkin solution of Richards' equation. This is done to improve the accuracy and mass balance properties of the companion advective transport model (finite volume-based), which rely on accurate velocity fields as input. Through a comparison between the results from the reconstructed velocities and the P1 Galerkin velocities, it is shown that a locally mass-conservative velocity field is necessary to obtain accurate transport results. Second, a detailed and novel analysis of the behavior of seepage face boundaries is performed with the flow model of CATHY_FT. The numerical simulations examine the model's performance under complex conditions such as heterogeneity and coupled surface/subsurface flow. It is shown that the overall numerical solution can be greatly affected by the way seepage face boundaries are handled in hydrological models and that careful considerations are required when using simple approximations, in the presence of heterogeneous slopes, and for seepage faces forming on a portion of the land surface. Third, CATHY_FT is implemented and run at the Landscape Evolution Observatory of the Biosphere 2 facility, Arizona. A detailed modeling analysis is performed of the experimental data collected during an isotope tracer experiment and from an intensively-measured hillslope, including quantity and quality of groundwater discharge and point-scale flow and transport data. This flow and tracer data is used to incrementally explore complex phenomena and associated hypotheses (e.g., heterogeneity, fractionation, and dispersion), progressing from flow to transport and from integrated to point-scale response analysis. This incremental approach highlights the challenges in testing and validating the new generation of integrated hydrological models when considering many types and levels of observation data. Finally, a concluding analysis is performed that relates to all three themes of the thesis, describing some of the features of the CATHY_FT model, discussing key issues associated to its further development, and testing its physical and numerical behavior for both real and synthetic scenarios. This final stage of the thesis addresses the myriad challenges faced in accurately and efficiently resolving the difficult behavior of the advection-dispersion equation for subsurface solute transport, in properly handling the complex boundary conditions for solute interactions across the land surface, and generally in capturing process interactions and feedbacks between flow and transport phenomena in surface and subsurface environments

Structural joint inversion of electrical and seismic tomography data

This research project has been focused on the achievement of the structural joint inversion of two geophysical methods. The final target is to obtain a high resolution characterization of the shallow subsurface. The aim of determining petrophysical properties, structural boundaries, etc, can be obtained through the integration of different information that derives from various geophysical methods. In fact, since each method is sensitive to a specific physical property, their integration can lead to an accurate final model. However, if such integration is conducted individually inverting the data sets, the final model will be affected by the resolution limitations of each method. For this reason, an important tool has been developed in geophysical applications: the joint inversion. Two different approaches can be used to carry out the joint inversion: the petrophysical one, in which a petrophysical relationship is used, and the structural one, in which a structural similarity between models is imposed (Gallardo and Meju, 2004). Specifically, I decided to implement the algorithm for the structural joint inversion and specifically the structural approach developed by Gallardo and Meju (2003, 2004), since from literature it results to be the most robust method in the joint inversion (Moorkamp, 2017). In this process, an objective function that includes the objective function of each geophysical method is build and simultaneously minimized. In conclusion, the joint inversion may improve the resolution of each geophysical model and bring to models that are more accurate and easier to interpret. Specifically, in this thesis, the electrical resistivity tomography (ERT) and the seismic refraction tomography (SRT) have been used to carry out the joint inversion. Both these high-resolution methods can be crucial in environmental and engineering applications, as for the geotechnical characterization of a site or for the detection of hydrological resources. Since the resistivity range overlaps for the different materials, resistivity measurements cannot be related to a specific soil or rock. Because of that, it would be better to obtain other information, for example from the seismic tomography. In fact, this method allows not only the reconstruction of the seismic wave velocities with depth, but also to obtain a good lateral resolution. Instead, the Ground Penetrating Radar (GPR) has not been considered since it presents some limits in the investigation depth, due to the high attenuation of electromagnetic energy in porous conductive media. In addition to the integrated inversion, another goal has been obtained in this thesis: the implementation of the forward modeling for the seismic method and specifically, the Multistencils Fast Marching Method (MSFMM). This method can be seen as an extension of the FMM, that is considered from literature the fastest and the most efficient method for the solution of the eikonal equation and accordingly for the computation of the first arrivals traveltimes. In particular, the MSFMM improves the accuracy and the efficiency of the FMM, since it considers also the information that derives from the diagonal directions. Both the algorithms, the one of the joint inversion and the one of the forward modeling for the seismic method, have been implemented in Python language and integrated in the open-source software pyGIMLi

A Case Study for Assessing the Hydrologic Impacts of Climate Change at the Watershed Scale

Since the advent of the industrial era atmospheric concentrations of greenhouse gases have been on the rise leading to increasing global mean temperatures. Through increasing temperatures and changes to distributions of precipitation, climate change will intensify the hydrologic cycle which will directly impact surface water sources while the impacts to groundwater are reflected through changes in recharge to the water table. The IPCC (2001) reports that limited investigations have been conducted regarding the impacts of climate change to groundwater resources. The complexity of evaluating the hydrologic impacts of climate change requires the use of a numerical model. This thesis investigates the state of the science of conjunctive surface-subsurface water modeling with the aim of determining a suitable approach for conducting long-term transient simulations at the watershed scale. As a result of this investigation, a coupled modeling approach is adopted using HELP3 to simulate surface and vadose zone processes and HydroSphere to simulate saturated flow of groundwater. This approach is applied to the Alder Creek Watershed, which is a subwatershed of the Grand River Watershed and located near Kitchener-Waterloo, Ontario. The Alder Creek Watershed is a suitable case study for the evaluation of climate change scenarios as it has been well characterized from previous studies and it is relatively small in size. Two contrasting scenarios of climate change (i.e., drier and wetter futures) are evaluated relative to a reference scenario that is based on the historical climatic record of the region. The simulation results show a strong impact upon the timing of hydrologic processes, shifting the spring snow melt to earlier in the year leading to an overall decrease in runoff and increase in infiltration for both drier and wetter future climate scenarios. Both climate change scenarios showed a marked increase to overall evapotranspiration which is most pronounced in the summer months. The impacts to groundwater are more subdued relative to surface water. This is attributed to the climate forcing perturbations being attenuated by the shift of the spring snow melt and the transient storage effects of the vadose zone, which can be significant given the hummocky terrain of the region. The simulation results show a small overall rise of groundwater elevations resulting from the simulated increase in infiltration for both climate change scenarios

On the use of surrogate-based modeling for the numerical analysis of Low Impact Development techniques

Mechanistic models have proven to be accurate tools for the numerical analysis of the hydraulic behavior of Low Impact Development (LIDs) techniques. However, their widespread adoption has been limited by their computational cost. In this view, surrogate modeling is focused on developing and using a computationally inexpensive surrogate of the original model. While having been previously applied to various water-related and environmental modeling problems, no studies have used surrogate models for the analysis of LIDs. The aim of this research thus was to investigate the benefit of surrogate-based modeling in the numerical analysis of LIDs. The kriging technique was used to approximate the deterministic response of the widely used mechanistic model HYDRUS-2D, which was employed to simulate the variably-saturated hydraulic behavior of a contained stormwater filter. The Nash-Sutcliffe efficiency (NSE) index was used to compare the simulated and measured outflows and as the variable of interest for the construction of the response surface. The validated kriging model was first used to carry out a Global Sensitivity Analysis of the unknown soil hydraulic parameters of the filter layer, revealing that only the shape parameter α and the saturated hydraulic conductivity Ks significantly affected the model response. Next, the Particle Swarm Optimization algorithm was used to estimate their values. The NSE value of 0.85 indicated a good accuracy of estimated parameters. Finally, the calibrated model was validated against an independent set of measured outflows with a NSE value of 0.8, which again corroborated the reliability of the surrogate-based optimized parameters