Analysis and improvement of the nonlinear iterative techniques for groundwater flow modelling utilising MODFLOW

As groundwater models are being used increasingly in the area of resource allocation, there has been an increase in the level of complexity in an attempt to capture heterogeneity, complex geometries and detail in interaction between the model domain and the outside hydraulic influences. As models strive to represent the real world in ever increasing detail, there is a strong likelihood that the boundary conditions will become nonlinear. Nonlinearities exist in the groundwater flow equation even in simple models when watertable (unconfined) conditions are simulated. This thesis is concerned with how these nonlinearities are treated numerically, with particular focus on the MODFLOW groundwater flow software and the nonlinear nature of the unconfined condition simulation. One of the limitations of MODFLOW is that it employs a first order fixed point iterative scheme to linearise the nonlinear system that arises as a result of the finite difference discretisation process, which is well known to offer slow convergence rates for highly nonlinear problems. However, Newton's method can achieve quadratic convergence and is more effective at dealing with higher levels of nonlinearity. Consequently, the main objective of this research is to investigate the inclusion of Newton's method to the suite of computational tools in MODFLOW to enhance its flexibility in dealing with the increasing complexity of real world problems, as well as providing a more competitive and efficient solution methodology. Furthermore, the underpinning linear iterative solvers that MODFLOW currently utilises are targeted at symmetric systems and a consequence of using Newton's method would be the requirement to solve non-symmetric Jacobian systems. Therefore, another important aspect of this work is to investigate linear iterative solution techniques that handle such systems, including the newer Krylov style solvers GMRES and BiCGSTAB. To achieve these objectives a number of simple benchmark problems involving nonlinearities through the simulation of unconfined conditions were established to compare the computational performance of the existing MODFLOW solvers to the new solution strategies investigated here. One of the highlights of these comparisons was that Newton's method when combined with an appropriately preconditioned Krylov solver was on average greater than 40% more CPU time efficient than the Picard based solution techniques. Furthermore, a significant amount of this time saving came from the reduction in the number of nonlinear iterations due to the quadratic nature of Newton's method. It was also found that Newton's method benefited more from improved initial conditions than Picard's method. Of all the linear iterative solvers tested, GMRES required the least amount of computational effort. While the Newton method involves more complexity in its implementation, this should not be interpreted as prohibitive in its application. The results here show that the extra work does result in performance increase, and thus the effort is certainly worth it