Lesson 3: Pump and Treat at the St. Joseph Superfund Site, Michigan, USA

Sections in this lesson:

• 3.1. Lesson Introduction
• 3.2. Sensitivity Analysis
  • Justification of Final Baseline Values
  • Site Variability
• 3.3. Pump-and-Treat Remediation
  • Ten Well Solution
  • Two Well Solution
• 3.4. Comparison of IGW to MODFLOW Pump-and-Treat Simulations
• 3.5. Summary
• 3.6. Step-by-Step Tutorial

Acknowledgement: This lesson and the accompanying tutorials were originally written by Carmen Nale in 2004 as part of her Master's Project in Environmental Engineering.

3.1 Lesson Introduction

This lesson uses IGW to recreate the modeling exercise conducted by Claire Tiedeman and Steven M. Gorelick (1993) in their Water Resources Research paper, “Analysis of Uncertainty in Optimal Groundwater Contaminant Capture Design”. Tiedeman and Gorelick modeled a Superfund site in Saint Joseph, Michigan with the three-dimensional, U.S. Geological Survey model called MODFLOW. The site is contaminated with TCE and its anaerobic transformation products. Tiedeman and Gorelick modeled plume response to different combinations and placements of extraction wells. Their goal was to design an optimal plume containment strategy and determine the reliability of containment.

Tiedeman and Gorelick's approach included the following steps. First, they compiled field data for aquifer properties such as porosity, conductivities, recharge, heads, and contaminant concentration. Then they used this information to determine the unknown site parameters (aquifer and creek conductivities and recharge) using STARPAC, a nonlinear weighted least squares regression routine. They then used the parameters values and their uncertainties in MODFLOW simulations to compare model heads and contaminant transport to measured values. Finally, they used MODFLOW with a nonlinear optimization solver, NPSOL, to design and evaluate the effects of different combinations and placements of extraction wells.

This lesson uses IGW to model two of Tiedeman and Gorelick's optimal plume containment designs. It begins with a site description and information on how to set up IGW to model the St. Joseph's site. In Section 3.2, it then describes a sensitivity analysis that was conducted to determine baseline values for use in the exercise. This section also compares the head contours and particle pathways modeled in IGW with those described by Tiedeman and Gorelick and explores the influence of site variability on the model setup. Section 3.3 is the heart of the lesson -- it recreates Tiedeman and Gorelick's 10 and 2 extraction well solutions within IGW. Finally, Section 3.4 and 3.5 compares the IGW results to those Tiedeman and Gorelick obtained with MODFLOW and summarizes the lesson. Finally, you can setup and run the simulation by following the step-by-step instructions that are in Section 3.6, the tutorial. You can also download the igw file that will run the simulation.

Site Description

In this lesson the IGW model is used to model a site in Saint Joseph, Michigan, which is contaminated with TCE and its anaerobic transformation products. The location of the site on the U.S.G.S. quadrangle map, Stevensville, MI, is shown in Figure 3.1 (Michigan, 1970).

Figure 3.1. The site at Saint Joseph, Michigan. The Saint Joseph Tutorial contains detailed instructions on how to include the topographic map into the IGW model.

The site has a shallow, unconfined sandy aquifer with two groundwater sinks, Lake Michigan to the west and Hickory Creek to the east. These two sinks create a groundwater divide. The main source of recharge is precipitation. The area is generally flat with a steep cliff by Lake Michigan and an incised channel at Hickory Creek. The source of the contamination is an automotive brake plant that disposed of trichloroethylene (TCE) wastewater into unlined lagoons. The flow amount and frequency of disposal are unknown so the recharge into the lagoons was not modeled by Tiedeman and Gorelick (1993) or in this lesson. The lagoon was located near the groundwater divide. This resulted in two plumes of contamination, one migrating towards Lake Michigan and the other towards Hickory Creek. Through intrinsic anaerobic transformations, some TCE transformed to vinyl chloride (McCarty et. al., 1990, and Semprini et. al, 1995).

Model Setup

The basic IGW model consists of four zones and a polyline (i.e. a small zone used to simulate small rivers, streams, and creeks) -- these are displayed on the topographical map shown in Figure 3.2. The constant head contours (red lines) and the groundwater velocity vectors (blue arrows) were generated using the IGW model. The K1 Zone (red) represents the basic characteristics of the water table aquifer system (i.e., conductivity, recharge, porosity, specific yield, surface elevation, bottom elevation of aquifer, etc.). The area surrounding Hickory Creek (K2 Creek Zone, yellow) consists of fine-grained floodplain deposits. The soil within the lining of the lagoons (K2 Lagoons, purple) is less permeable relative to the rest of the site (Tiedeman and Gorelick, 1993). The Lake Michigan Zone (blue) is modeled as having a constant head. The Hickory Creek Zone (light blue) is represented as a polyline due to its small area and given a value of constant head.

Figure 3.2. Baseline map of the Saint Joseph used for the sensitivity analysis.

Figure 3.2 illustrates the groundwater divide phenomena with the lagoons directly on the divide with one groundwater receptor on each side, Lake Michigan to the west, and Hickory Creek to the east. The shape of the head contours is created by recharge infiltrating into the aquifer then flowing from the groundwater divide (higher head) to the receptors (lower head). Figure 3.2 was the baseline figure used in the following simulations with variations in certain model parameters. The baseline figure is converted into four copies using the ‘Create Multiple Models’ button. It allows the user to change one or more parameters in any of the four models to visually compare the results at the same time. Be aware that the models may interact with each other and produce an error in output if the spacing between the models is too small and/or the grid size is too coarse.

3.2 Sensitivity Analysis

This section describes the sensitivity analysis that was conducted set up conditions within IGW that will match the head and plume conditions Tiedeman and Gorelick (1993) reported for the site. The sensitivity analysis varies the following parameters: K1 (aquifer hydraulic conductivity), K2 (creek sediment and lagoon hydraulic conductivity), and Recharge (constant annual recharge). The goal of the sensitivity analysis was to find parameter values that will match the measured head contours and create particle pathways similar to the direction of travel and location of the plumes reported by Tiedeman and Gorelick (1993) (Figure 3.3).

Figure 3.3. Contaminant plumes and observed heads as reported by Tiedeman and Gorelick (1993).

The sensitivity analysis makes use of the statistics given for the parameter estimates in the Tiedeman and Gorelick paper (1993), shown in Table 3.1. The values they provided for the mean, standard deviation, and 95% confidence interval for the K1, K2, recharge multiplier, and annual recharge parameters were used throughout the sensitivity analysis.

The mean values for each parameter were used as a reference point and are shown in the middle of each simulation graphic. Table 3.1 shows these values along with other site characteristics used in the model. The simulation graphics also show particle pathways. These were developed using IGW's particle tracking feature. The particles are added to the lagoon area and then selected as continuous particles, in order to keep track of their pathway. Detailed instructions for all Saint Joseph site simulations are provided in the tutorial.

Variation of K1

Figure 3.4 displays the model's sensitivity on the head field and the resulting particle pathways with the range of K1 values given in Table 3.1. The figure shows the result of an increasing K1 value. The values increase clockwise from the lowest value of 4.92 m/day in the bottom left-hand quadrant to the highest value of 10.09 m/day in the bottom right-hand quadrant, the mean K1 value from Table 3.1 is shown in the middle of the graphic.

Figure 3.4. Sensitivity analysis on K1.

Notice that the head value decreases as KI increases since there is less resistance to fluid flow. Also as K1 increases, the groundwater divide shifts to the right causing about half the particles to travel toward Lake Michigan when K1 equals 10.09 m/day instead of Hickory Creek when K1 equals 4.92 m/day. The low K1 value appears to be the best fit for the head distribution shown in Figure 3.4 since the value of the measured head near the lagoon is around 183 to 184 m. However, the high K1 value shows a better distribution of the particles path since there is an equal distribution of particles traveling to Lake Michigan and Hickory Creek. This suggests that further analysis is necessary in order to determine the best K1 value.

Variation in Recharge

Figure 3.5 illustrates the effects of the variation of constant recharge on the other baseline parameters. The values of constant recharge increases from the bottom left-hand corner value of 22.7 cm/year then clockwise to the highest value of 57.1 cm/year. The mean value from Table 3.1 is shown in the middle.

Figure 3.5. Sensitivity analysis on constant recharge.

Notice that as recharge increases the head also increases -- this is expected since more water percolating into the aquifer would increase the overall head values. In terms of a best fit, the high recharge value of 57.1 cm/year appears to create the desired head profile but the placement of the divide is too far to the left and causes the particles to go only into Hickory Creek. The desired particle pathway is produced by the simulation with the lowest recharge value of 22.7 cm/year, but the head profile does not match the measured values as well. Further analysis is required in order to determine the best recharge value that fits both the head profile and the particle pathway within reasonable limits.

Variation of K2

Figure 3.6 shows the results of varying K2 to the values given in Table 3.1. In the figure, the value of K2 increases going clockwise, starting at the lowest value of 1.35 m/day in the bottom left-hand corner and then increases to the highest value of 3.13 m/day. The mean value of 2.24 m/day is shown in the middle.

Figure 3.6. Sensitivity analysis on K2.

Notice that as the value of K2 decreases, fewer particles tend to flow towards the creek because of the increase in resistance to flow. As K2 decreases, the head values increase and the contour lines are closer together indicating a steep head gradient. Also, the relatively low permeability of the soil underlying Hickory Creek and the lagoon compared to K1 causes an increase in resistance for water to flow towards the creek and the particles to travel out of the lagoon. The results show that using the lowest K2 value best fits both the head distribution and the particle pathways. Therefore, the K2 value of 1.35 m/day will be used.

Variation of K1 with K2 = 1.35 m/day

With a value for K2 selected, the next step is to redo the sensitivity analysis for K1 with K2 equal to 1.35 m/day. The following simulations use the baseline parameters shown in Table 3.1 with an exception of K2 equaling 1.35 m/day. Figure 3.7 shows the results of this sensitivity analysis.

Figure 3.7. Sensitivity analysis on K1 with K2 = 1.35 m/day.

Notice that as K1 increases, the value of the head contours tend to decrease due to less resistance to flow. However, more of the particles tend to flow towards Lake Michigan as K1 increases. The baseline K1 value of 7.51 m/day gives the best results when considering the particle pathways, but the head values are still low near the lagoon.

Variation in Recharge with K2 = 1.35 m/day

Figure 3.8 shows a repeat of the sensitivity analysis on recharge with the values given in Table 3.1 and the modification of the K2 value to 1.35 m/day. As recharge increases from 22.7 cm/year to 57.1 cm/year, the value of the head profile increases and shifts to the left since more water is available within the aquifer. When comparing Figure 3.8 with Figure 3.5, notice that recharge has a profound effect on the heads but very little effect on the particle pathways. A recharge value of 48.5 cm/year gives the best head distribution.

Figure 3.8. Sensitivity analysis on recharge with K2 = 1.35 m/day.

In order to determine if this value is reasonable, a water balance was conducted on the area using the following precipitation, evaporation, and pan coefficient values needed for the recharge equation (Bedient and Huber, 1992) given that:

Average Precipitation (NOAA, Benton Harbor, 1951-1980) = 36.41 inches/yr
Average Evaporation (NCDC, South Haven, 1952-1978) = 35.41 inches/yr
Pan Coefficient (Huber, 1992) = 0.7 to 0.8
Recharge = Precipitation – Pan Coefficient*Evaporation = 30 to 21 cm/yr, respectively

The range of recharge values is from average to conservative in respect to the pan coefficient range given. The recharge range calculated agrees with the value determined by Tiedeman and Gorelick (1993) of 29 cm/year. For the MODFLOW model simulations, they applied a recharge multiplier of 1.4 to get 40 cm/year compared to the IGW value of 48.5 cm/year with a recharge multiplier of 1.7. This seems to be a reasonable value due to many unknown factors including land use, the pan coefficient value, and detailed aquifer properties and water sources/uses. The value of constant recharge used for the following simulation is 48.5 cm/year.

Variation of K1 with K2 = 1.35 m/day and Recharge = 48.5 cm/year

With a value seclected for constant recharge, the next step is to repeat the sensitivity analysis for K1 and select a hydraulic conductivity of the aquifer. Figure 3.9 shows this analysis with the range of K1 values from Table 3.1 and the chosen values for K2 and recharge.

Figure 3.9. Sensitivity analysis on K1 with K2 = 1.35 m/day, Recharge = 48.5 cm/year.

As expected as K1 increases the value of the head contours decrease and the particles move away from Hickory Creek toward Lake Michigan. The best range of K1 values is 7.51 to 8.80 m/day -- these most closely match both the value of the head contours and the pathway of the particles for the known contaminant profile at the site. The results from the more refined range are shown in Figure 3.10.

Figure 3.10. Sensitivity analysis on K1, K2 = 1.35 m/day, Recharge = 48.5 cm/year.

From the head contour values and the pathway of the particles in these simulations, it appears that the best K1 value is found to be between 7.51 and 7.94 m/day. As a check, the particle tracking shown in Figure 3.10 is replaced with a continuous contaminant source with an arbitrary concentration. This simulation is illustrated in Figure 3.11 (A) and (B) for the range of K1 values shown in Figure 3.10. In all simulations the plume's shape is very similar only the plume's extent (time traveled) varies. For our simulation, travel time of the plume is unknown. Given the particle pathway and head contour values shown in Figure 3.10, and the unknown time of plume transport shown in Figure 3.11, we chose to use a K1 value of 7.7 m/day.

Figure 3.11. Continuous plume sensitivity analysis on K1, K2 = 1.35 m/day, Recharge = 48.5 cm/year. (A) shows a picture time of 18 years while (B) shows a picture time of 28 years.

Table 3.3 shows the final baseline values (italics) for K1, K2 and Constant Recharge determined from the sensitivity analysis.

Figure 3.12. Illustration of final baseline values.

Justification of Final Baseline Values

Comparing IGW to Observed Head Values

To determine the accuracy of the IGW model, the final baseline values are used to compare the observed (black contours) to the model heads (filled-in colored contours) shown in Figure 3.13.

Figure 3.13. Comparing measured (black) to IGW (filled-in color) head contours.

In general, the measured head contours are abrupt compared to IGW’s smooth and rounded head contours. The abrupt head contours close to the lagoon could be replicated in the model if more head information was available to calibrate the model. However, the measured and IGW head values match closely and increase in accuracy moving from the lagoon toward Lake Michigan. Figure 3.13 shows that Tiedeman and Gorelick’s (1993) illustration of the measured head contours was not scaled appropriately in respect to the U.S.G.S. topographical map used in the IGW simulation (i.e., the location of Hickory Creek does not match) this may explain the difference in head contours from the lagoon to Hickory Creek.

Figure 3.14. Monitoring well number and location used to compare measured (black) to IGW head contours.

To determine the accuracy of the model's head distribution, we compared the head value at the location of each monitoring well in the Tiedeman and Gorelick paper (shown in Figure 3.14). Table 3.4 displays the results of the measured heads, model heads, head difference, and percent difference in reference to the measured head values. As expected, the head variation increases close to the lagoon location due to the inability of the model to simulate those heads accurately. The overall head field could be improved by calibrating the model with more head data for the site. However, the model does fairly well overall at modeling the heads -- the maximum difference between IGW and the Tiedeman and Gorelick values is 4.5%.

Comparing IGW Particle Pathway to Plume Location on Site

The particle pathways are compared to the plume location shown in Figure 3.15 (Tiedeman and Gorelick, 1993).

Figure 3.15. Comparison of observed head and plume location to final baseline parameter values.

The location of the plume traveling to Hickory Creek matches the pathway of the particles modeled in IGW. The plume heading towards Lake Michigan is just south of the particle pathways, this is most likely due to the variation between model and measured heads around the lagoons. Since the lagoon is on the groundwater divide, the head field around this area has a profound effect on the particle pathways. Again, to correct this variation, the model would need to be calibrated in more detail.

Comparing IGW to MODFLOW's Profile of a Particle Pathway

A cross-section from A to A' in Figure 3.15 helps show the particle pathway from the lagoon to Lake Michigan. This cross-section was compared to one given in Tiedeman and Gorelick’s paper (Figure 3.16). The main difference between the two models is the ability to accurately characterize the bottom of the aquifer. Tiedeman and Gorelick (1993) used the three-dimensional model MODFLOW that incorporates vertical characteristics of the aquifer. The IGW model is a two-dimensional model that only has the ability to model the aquifer bottom by adding scatter points. In this lesson an average elevation was used since that was the only information available. Despite this difference, the resulting head distribution from the lagoon to Lake Michigan is similar for both models. When looking at Tiedeman and Gorelick’s (1993) profile (Figure 2.16(B)), notice that the head values for the lake and lagoon are 1 meter too high, they gave an average elevation value for Lake Michigan as 176.5 m and the lagoon as 184 m.

Figure 3.16. (A) IGW and (B) paper profile of particle pathway from lagoon to Lake Michigan.

When comparing the correct head values, both profiles match well. Even the pathway of the particle is comparable when considering that the IGW model has no detailed information of the aquifer’s bottom. Notice that for both models recharge is causing a strong vertical migration of the particle towards the bottom of the aquifer at the lagoon location, and then travels relatively along the bottom of the aquifer to Lake Michigan.

Site Variability

Lake Michigan: Steady State versus Transient State

The final baseline parameters shown in Table 3.3 are used to develop Figure 3.17. This will be the reference figure for the following simulations. One major concern is the influence of Lake Michigan’s fluctuating head on the pathway of the particle.

Figure 3.17. Final baseline figure used for reference.	Figure 3.18. Transient head for Lake Michigan.

From 1970 to 2002, the minimum, average, and maximum head values of Lake Michigan at Holland, MI were determined to be 175.7 m and, 176.7 m, 177.6 m, respectively. To determine the impact of a fluctuating lake head, a constant head value of 176.5 m is replaced with a transient head trend shown in Figure 3.18. The input of the actual measured lake levels over time can not be put into the model simulations. Therefore, the random fluctuation box is activated and then the trend data is set with 176.5 m as the head value at 0 and 360 days with an amplitude value of 2 m, about the difference between the maximum and minimum head values. When comparing steady state (Figure 3.17) to transient state (Figure 3.18), more particles tend to go to the creek due to the head fluctuation of Lake Michigan. This is good to know, since the actual conditions could be at either extreme. The high and low head values are most likely to occur during periods of floods and droughts, respectively. Therefore, the steady state solution will be used for the remediation section for comparison reasons and since the transient state does not significantly change the course of the particle pathways.

Uniform versus Random K1 and K2 Fields

Influence on Particle Pathway

When characterizing a site, it is difficult to determine accurate aquifer parameters using limited information collected onsite. Therefore, the site's conductivities (K1 and K2) are changed from a uniform to random conductivity field to determine the sensitivity of the particle pathways and head distribution as a result of a heterogeneous distribution. Results are shown in Figure 3.19 using (A) discrete and (B) continuous particles.

Figure 3.19. (A) Influence of a random K1 and K2 field on discrete particles and (B) influence of a random K1 and K2 field on continuous particles.

The random conductivity field causes the head contours to be wavy resulting in different particle pathways compared to the uniform field. When comparing the discrete particles with a random K1 and K2 field to the uniform field shown in Figure 3.19 (A), notice that the particles motion is more sporadic compared to the evenly spaced particles, respectively. The main difference between the two fields is the general pathway of the particles. With a uniform field they travel to both the lake and creek, but the random field causes the particles to go only to Hickory Creek shown in Figure 3.19 (B). This difference could be caused by the slight shift of the groundwater divide to the left in the random field simulations causing the particles to travel to the creek only. The random field simulation seems to have a greater effect on the head distribution, causing the head values to slightly increase. This could be a result of having a mixture of low and high conductivities creating an increase in fluid resistance. More simulations are needed to show the effects of different random solutions, which is beyond the scope of this lesson. The influence of the random field on a contaminant concentration from an instantaneous source is discussed next.

Influence on an Instantaneous Plume

Figure 3.20 shows the influencing of a uniform versus random K1 and K2 field at time 12, 20, and 30 years on a plume created by an instantaneous source, respectively. When comparing the head contours between Figure 3.20 to 3.19, notice that the head contour values are generally the same with more variability when the field is random. In Figure 3.20, the uniform field has very smooth, even concentration contours compared to the random field that has sporadic contour lines. Even with sporadic contour lines, the general shape and size of the plume is similar for both fields. After 12 years the plume has reached Hickory Creek for both the uniform and random field, while it takes about 30 years for the plume to reach Lake Michigan. Notice that the highest contaminant concentration contour (red) only goes to the creek for the random field, following a similar pathway as the particles shown in Figure 3.19. While the lower contaminant concentration contours (green and blue) go to both the creek and lake. This suggests that for a random conductivity field, a higher concentration of the contaminant goes towards the creek, while a lower concentration goes towards the lake. For the uniform conductivity field, most of the high contaminant concentration goes to the creek but some also goes to the lake. This indicates that the particle pathways represented in Figure 3.19 show the effects of advection while the plume represented in Figure 3.20 also includes the effects of diffusion and dispersion.

Figure 3.20. Random K1 field, both K1 and K2 fields influence on an instant plume at time various times.

Generally, the average and random conductivity fields show similar transport of the contaminant. Therefore, the uniform K1 and K2 field are used for the following simulations because using a random field would not necessarily give us more accurate results.

Continue to Section 3.3. Pump and Treat Remediation

Back to top