A global optimization (metaheuristic) method, tabu search, is integrated with linear programming to solve remediation design problems. This integrated approach takes advantage of the fact that the global optimization approach is most effective for optimizing discrete well location variables, while linear programming is much more efficient for optimizing continuous pumping rate variables. In addition, an efficient forward solution updating procedure is used to lessen the computational burden of the global optimization approach. With this procedure the new solution to a linear flow model perturbed by pumping is obtained as the sum of a nonperturbed base solution and the solution to the perturbed portion of the flow system, which can be derived directly without running the flow model. Numerical results, based on a two-dimensional capture zone design problem, show that the computation time can be reduced to a small fraction of that required by the conventional approach, in which a forward simulation model is run each time the objective function needs to be evaluated. It is also demonstrated that the maximum number of wells allowed in a given design has a significant effect on the total remediation costs. (The total remediation costs are nearly doubled when only one well is allowed instead of the optimal number of six for the test problem.) A Monte Carlo analysis, based on 200 realizations of a lognormally distributed random hydraulic conductivity field (the variance of ln T=1.0), further reveals that the total remediation costs determined for the heterogeneous aquifer have a large uncertainty (the ratio of standard derivation over mean is 0.4). The total remediation costs and associated uncertainty are also shown to increase with the uncertainty of the hydraulic conductivity field. ¿ 1999 American Geophysical Union