Numerical modeling of groundwater flow

Numerical modeling of groundwater flow, transport and remediation

Large scale modeling of groundwater flow, transport and remediation is one my major research focus areas. A suite of parallel groundwater flow, transport and remediation codes, PGREM3D, have been developed for this purpose. These codes have been used for a variety of applications including field simulations. These codes can be obtained free of charge for research purposes. Please send e-mail to gmkumar@eos.ncsu.edu if you are interested in obtaining these codes. The entire suite contains the following:

Parallel groundwater flow code: three-dimensional steady state saturated flow, hexahedral finite-elements, conjugate gradient and multi-grid solvers, capable of handling horizontal and vertical wells in heterogeneous media.
Parallel groundwater transport code: three-dimensional multi-component transport, hexahedral finite-elements, mobile and immobile components, sorption, bioremediation, and reversible kinetic reactions, BiCGSTAB and GMRES solvers, handles flow field from flow code.
Other codes that may be available on request (not automatically available in the suite):
- A parallel turning bands code: generates three-dimensional random hydraulic conductivity fields for the flow code using the turning bands method (a slightly modified version of the original serial code from Andy Thompson of LLNL).
- A parallel random walk particle tracking code: three-dimensional, single component, handles dispersion and sorption, handles flow field from flow code.
- A hybrid genetic algorithm - local search optimizer: currently available only as a coupled module with the transport code, specifically tuned for parameter estimation problems.

Some Groundwater Contamination and Remediation Animations

The animations show contaminant plume evolution and remediation using different types of extraction wells in a highly heterogeneous subsurface environment. In these simulations we have tried to capture and represent the effects of natural formation heterogeneity on transport and recovery processes which require a high degree of spatial resolution.

The natural formation heterogeneity is generated for a hypothetical aquifer using a non-gaussian fractal geostatical simulator based on fractional Levy motion. Non-Gaussian fields represent true heterogeneous formations in natural subsurface environments that cannot be easily obtained by traditional guassian simulators.

The simulations were performed on a 200 x 100 x 40 rectangular grid using 72 processors of the Intel Paragon XPS/150. The flow calculations were performed for the steady state saturated case. It required less than a minute to solve the flow problem and approximately 2 hours to perform 5000 timesteps of the transport simulations. The flow code used a multigrid preconditioned conjugate gradient solver and the transport code a diagonally preconditioned BiCGSTAB solver. All the wells used in these simulations used a prescribed flow rate condition. The flow rates for individual well nodes are distributed using an iterative procedure requiring several matrix solves. One contamination scenario and three remediation scenarios using both horizontal and vertical wells are shown in the animations.

Contamination:
Contamination from a single rectangular source located at the top of the aquifer over 2000 days. The source concentration is maintained constant at 100 and isosurfaces are shown for concentrations of 1 (blue), 4 (yellow), and 25 (red). The heterogeneous permeability field is shown in the background.

Remediations:

Single horizontal extraction well: One extraction horizontal well centered approximately around the concentration weighted centroid of the plume. The remediation is stopped when approximately 90% of the contaminants are removed. Note that most of the contaminant that is present are at a concentration level of 1 even though the average initial plume concentration is around 100. The following two animations are identical except for the transparency property of the isosurfaces. The left animation has opaque isosurfaces, and the right animation's isosurfaces have varying degrees of transparency.

Two horizontal extraction wells: Two extraction horizontal wells, with the centroid of these wells centered approximately around the concentration weighted centroid of the plume. The combined flowrate from these two wells is the same as the single well flow rate in Remediation 1. The remediation is stopped when approximately 90% of the contaminants are removed.

Two vertical extraction wells: Two extraction vertical wells are used for plume extraction, with the centroid of these wells centered approximately around the concentration weighted centroid of the plume. The combined flowrate from these two wells is the same as the single well flow rate in Remediations 1 and 2. The remediation is stopped when approximately 90% of the contaminants are removed.

From the above simulations we can see that heterogeneity has a profound impact on remediation. Using multiple wells with the same combined flow rate does not seem to provide improved removal efficiency for the cases tested here. However, additional simulations that were performed later (but not shown here) show that they provide better removal efficiency if they are farther apart than the cases shown here.