Exploring Variably Saturated Flow Formulations for MODFLOW: Implementing a Fully 3D Richards' Equation in MODFLOW 6

This thesis explores the implementation of a fully 3D Richards’ equation solution into MODFLOW 6 to solve for saturated/unsaturated flow. Currently, MODFLOW 6 has limited capability for solving unsaturated zone flow problems by using the Unsaturated Zone Flow (UZF) package, which neglects capillarity and assumes vertical 1-dimensional flow. Various formulations of the Richards’ equation were first evaluated for their strengths and weaknesses. Various numerical solution approaches to Richards’ equation were also studied to understand the robustness and efficiency of different approaches and to select one that would best fit into the MODFLOW 6 framework. The MODFLOW 6 framework was also studied to understand how the Unsaturated Zone Richards (UZR) formulation would fit into it. A UZR (Unsaturated Zone using Richards’ equation) package was then developed to provide MODFLOW 6 with the capability of solving three-dimensional saturated/unsaturated flow. The UZR package solves for three-dimensional groundwater and unsaturated zone flow using a mixed form of Richards’ equation. Implementation into MODFLOW 6 involved replacing the relative permeability and saturation terms in place of the upstream transmissivity and grid-block saturation fraction of unconfined flow. This was done for Picard and Newton-Raphson implementations of the flow equation. The Newton-Raphson formulation provides a computationally robust and efficient approach for solving Richards’ equation. To evaluate the UZR package, a series of test problems were conducted to benchmark its performance against existing software, particularly HYDRUS, MODFLOW VSF, and USG-T. Results indicate that the UZR package successfully provides accurate simulations of variably saturated flow compared with other software, especially USG-T, even under challenging conditions like perched aquifers and transient recharge events. Additionally, runtime comparisons demonstrate that UZR maintains computational efficiency, positioning it as a practical tool for large-scale applications. This study contributes a robust solution to MODFLOW 6, providing a foundation for future advancements in groundwater modeling and offering a valuable resource for hydrogeologists, environmental engineers, and researchers seeking to enhance their modeling of unsaturated flow dynamics. Advisors: Tirthankar Roy and Sorab Panday