Fast iterative implementation of large-scale nonlinear geostatistical inverse modeling

In nonlinear geostatistical inverse problems, it often takes a significant amount of computational cost to form linear geostatistical inversion systems by linearizing the forward model. More specifically, the storage cost associated with the sensitivity matrix H (m x n, where m and n are the numbers of measurements and unknowns, respectively) is high, especially when both m and n are large in for instance, 3-D tomography problems. In this research, instead of explicitly forming and directly solving the linear geostatistical inversion system, we use MINRES, a Krylov subspace method, to solve it iteratively. During each iteration in MINRES, we only compute the products Hx and H(T)x for any appropriately sized vectors x, for which we solve the forward problem twice. As a result, we reduce the memory requirement from O(mn) to O(m)+O(n). This iterative methodology is combined with the Bayesian inverse method in Kitanidis (1996) to solve large-scale inversion problems. The computational advantages of our methodology are demonstrated using a large-scale 3-D numerical hydraulic tomography problem with transient pressure measurements (250,000 unknowns and approximate to 100,000 measurements). In this case, approximate to 200 GB of memory would otherwise be required to fully compute and store the sensitivity matrix H at each Newton step during optimization. The CPU cost can also be significantly reduced in terms of the total number of forward simulations. In the end, we discuss potential extension of the methodology to other geostatistical methods such as the Successive Linear Estimator.