WaterSubject

Distributed parameter models are used to simulate groundwater flow, but the reliability of results is strongly dependent on whether the model parameters are properly identified. However, the real geological structure of an aquifer is generally very complex and mostly unknown. More data are required for identifying a more complex parameter structure, and there is always the problem of over-parameterization. Additionally, the inverse problem of parameter structure identification is highly non-linear and non-convex; hence, gradient-based algorithms often get trapped in local optima. On the other hand, heuristic search algorithms usually require a very large number of simulation runs. There remains a need to develop efficient inverse solution algorithms.

Although the inverse problem in groundwater modeling has been studied for over four decades, identifying the spatial distribution of a heterogeneous aquifer remains formidable because of the limitation in both quantity and quality of data. The inverse problem is inherently ill-posed in that the solution is usually non-unique or unstable. As pointed out by Yeh [1], parameterization is always necessary in order to reduce the parameter dimension of a distributed parameter to a finite dimensional form so that a stable and unique solution of the inverse problem can be sought.

Parameterization as well as the efficiency and effectiveness of the inverse algorithm are the two most important components in parameter structure identification. In this study, Tabu search is allied with the adjoint state method to improve the search efficiency and parameterization is achieved by Voronoi tessellation (VT).

Sykes et al. [2] developed adjoint sensitivity theory for both the continuous and discrete (numerical) equations of two-dimensional steady state flow in a confined aquifer. They derived a performance measure with respect to the system parameters, such as the piezometric heads, velocities, travel time, and mass discharge with respect to recharge–discharge rates, prescribed boundary heads or fluxes, thicknesses, and hydraulic conductivities. As Yeh [1] pointed out, there are basically three methods, the influence coefficient method, the sensitivity equation method, and the variational or adjoint state method (ASM) for calculating the Jacobian matrix (parameter sensitivity matrix). The ASM is particularly efficient when the sensitivities of a large number of parameters are required. Sun and Yeh [3] proposed a general procedure for deriving the adjoint state equations and described their associated conditions for solving coupled inverse problems in groundwater modeling.

In recent years, heuristic algorithms have been used to identify model parameters in groundwater modeling. Zheng and Wang [4] used Tabu search (TS) and simulated annealing (SA) to search for an optimal hydraulic conductivity zonation structure in one-dimensional problems. Tsai et al. [5] developed a sequential global–local optimization method, which consists of a genetic algorithm (GA), a quasi-Newton method, and local search to solve the combinatorial problem. More recently, Tung and Chou [6] and Tung and Tan [7] applied heuristic optimization algorithms, such as TS and SA, with different zonation methods to identify the spatial distribution of hydraulic conductivity or transmissivity. This study proposes to use the ASM to improve computational efficiency of TS. We note that TS requires running the simulation model N times to evaluate N neighbor solutions. More neighbor solutions require more simulation runs. When TS is allied with the ASM, it only requires running the simulation model twice regardless of how many neighbor solutions there are. Details will be explained in a later section.

To reduce the parameter dimension of a distributed parameter to a finite dimensional form, either zonation or interpolation can be used [1]. Doherty [8] used the pilot point method as an interpolation method to identify the spatial distribution of transmissivity. Tsai et al. [5] and Tung and Tan [7] applied Voronoi tessellation (VT) to automatically parameterize the zonation pattern of hydraulic conductivity or transmissivity. Tung and Chou [6] proposed a pattern zonation based on the pattern classification to identify the spatial distribution of pumping.

In this paper, we first use VT to automatically parameterize the distributed hydraulic conductivity for a given parameter dimension, represented by the number of basis points. Second, the ASM is used to improve the efficiency of the TS search. To further improve computational efficiency, we apply a coarse-fine grid search over the modified TS. For each parameter dimension, we calculate the fitting residual error, the parameter uncertainty error and the modified Akaike Information Criterion. We use these three indices to determine the optimal level of parameter complexity. Numerical experiments are carried out to demonstrate the validity and computational efficiency of the proposed approach.

The rest of this paper is organized as follows. In Section 2, the inverse problem and its formulation for parameter structure identification in groundwater modeling are introduced. The stopping criteria used to determine the optimal level of parameter complexity are also described. The proposed integrated optimization algorithm, which consists of VT, TS and the ASM are illustrated in Section 3. Section 4 outlines the numerical experiments with hypothetical study fields and the proposed integrated optimization algorithm is applied with the coarse-fine grid search technique. We then summarize the results with conclusions.