Capture Zone Analysis Using MODFLOW with Groundwater Vistas
Capture Zone Analysis using MODFLOW with Groundwater Vistas
We will simulate the flow system described below using a computer program called MODFLOW with a pre- and post-processor known as Groundwater Vistas. MODFLOW was developed by the USGS; it is the most popular groundwater flow model in the U.S. It solves the general form of the groundwater flow equation and therefore, can simulate three-dimensional, heterogeneous, transient or steady-state conditions.
Directions for using MODFLOW with Groundwater Vistas are attached below. A description of the problem is given below.
Description of the Problem
Conceptual Model
The unconfined aquifer shown in Fig. 1 consists of outwash filling a buried glacial valley. The outwash has a maximum thickness of 100 ft. The hydraulic conductivity of the aquifer is estimated to be 30 ft/day. Recharge from precipitation is estimated to be 0.006 ft/day (around 24 in/yr). Runoff from upland areas provides an additional recharge of 0.003 ft/day along the sides of the valley, for a total of 0.009 ft/day. A shallow perennial river flows through the center of the valley. There is a two foot thick (b) layer of sediment at the bottom of the river. The vertical hydraulic conductivity (Kz) of the riverbed sediments is estimated to be 2 ft/day, so that the leakance (Kz/b) beneath the river is is 1 day-1. A shallow well located close to the river pumps at a rate of 50,000 ft3/day (around 0.4 MGD).
Modeling Objective
The objectives are to delineate the capture zone of the well and to quantify the effects of pumping on groundwater flow to the river.
Methodology
Use MODFLOW to:
1. calculate the steady-state distribution of heads without pumping;
2. calculate the steady-state distribution of heads with pumping;
3. calculate the drawdown caused by pumping;
4. calculate the fluxes in and out of the river under pumping conditions.
This will require two MODFLOW runs: run 1 is a steady-state simulation without the pumping well and run 2 is a steady-state simulation with the pumping well. Run 2 uses the head distribution computed in run 1 as the starting heads so that drawdown may be calculated.
MODFLOW Design
Grid
We will use a three layer model with 11 columns and 23 rows. The grid is shown in Fig. 2. Note that x = y = 200 ft. Layer 3 is 60 ft thick; layer 2 is 20 ft thick and layer 1 has variable thickness as defined by the water table. Note the location of the river and the well. Also note that the simulated aquifer narrows with depth to approximate the geometry of the buried valley. That is, layer 1 has 11 columns of active nodes, layer 2 has 7 columns and layer 3 has only 5 columns of active nodes.
Boundary Conditions
(1) The aquifer is underlain by impermeable bedrock that forms a no flow boundary. (2) At the top of the problem domain, water is recharged to the aquifer forming a specified flow boundary condition. (3) The head in the river is equal to 100 ft and forms a specified head boundary in layer 1, column 6. (For our purposes we will assume that the gradient of the river is zero within the problem domain.) (4) The boundaries at the east and west edges of the aquifer represents no flow conditions since the surrounding bedrock is relatively impermeable (Fig. 1). (5) The aquifer extends farther to the north and south beyond the problem domain shown in Fig. 2. Therefore, the boundaries at the north and south ends of the model are not physical boundaries. We will use hydraulic boundary conditions here; we assume no flow hydraulic boundaries at the north and south ends of the model so that flow lines will be directed toward the river.
Parameters
Hydraulic Conductivity. We will assume that the aquifer is homogeneous and isotropic and the hydraulic conductivity is therefore a constant and equal to 30 ft/day.
Leakance. Leakance is equal to vertical hydraulic conductivity divided by thickness. Confirm that leakance for layer 1 is 1.5 day-1, except for the river node where leakance is equal to 1.0 day-1. Leakance for layer 2 is 0.75 day-1. By convention in MODFLOW, the bottom layer (our layer 3) has no leakance.
Recharge rate. Recharge is equal to 0.006 ft/day everywhere in layer 1, except for the river nodes in column 6 where recharge is equal to zero and for all the nodes in columns 1 and 11, where recharge is equal to 0.009 ft/day.
Well discharge rate. The well is located in layer 1, column 7, row 12. It has a discharge rate of -50,000 ft3/day.
Layer type. Layer 1 is a type 1 (unconfined) layer. We will simulate layers 2 and 3 as type 3 (variable) layers.
Convergence criterion. We will use a convergence criterion of 0.001 ft.
Figure 1.
Using Groundwater Vistas to solve the Capture Zone Problem
Insert a disk in the a: drive. Create two directories on the a: drive in which you will store the files from the GW Vistas simulations. Call one of them nopump and the other pump. Open GW Vistas.
Run 1. Steady-state Flow without Pumping (nopump).
Design the grid and assign default parameter values. When the initial screen appears, click file>new. The "Initialize Model Grid" screen will appear. Go through the screen, entering the correct default parameter values. (Double click on the parameter value in the box, and then enter a new value to replace the existing value.) We will adjust the z spacing later. So for now, just let that value equal 100. Remember to enter default values for hydraulic conductivity in the x- y- and z- directions (30), for leakance (1.5) and for recharge (0.006). Storage and porosity are not used in this problem and do not need to be entered. Click on the OK box when finished with this screen. (Length units are not specified when using GWVistas; the user should take care to use consistent units. The selection of the time unit is done later.)
The grid will appear on the screen. You can change the font size or the numbering of the rows and columns by selecting grid>options.
Define property zones. Now select Props (for properties). Note that there is a check mark next to hydraulic conductivity. Since the hydraulic conductivity is uniform and we have already set it as a default value we don't need to change it. Instead click on top elevation. Click on Props again and notice that the check mark is now on top elevation. Then choose property values>database. Set the top of zone 1 to be 110, zone 2 to be 80 and zone 3 to be 60. Note the cube pictured on the left side of the screen. Note that the cube shows that we are now on layer 1. Now click on set zone numbers>window. You are now going to set the top elevations for layer 1 using the window option. Drag the window over the grid to cover the entire layer. Then set the zone number for this layer to be zone 1 (i.e., top elevation = 110). Then go to the cube and click on layer number to get layer 2. Repeat the process of setting the top elevation for layer 2 by dragging the window over the grid and then setting the zone equal to zone 2. Repeat this process for layer 3, setting the zone number to 3.
Now go back to Props and click on bottom elevation. Repeat the process used in setting the top elevations, but now set bottom elevations to be 80 for layer 1, 60 for layer 2 and zero for layer 3.
When you are finished, note that as you drag the arrow across the grid, the value of the parameter (i.e., the bottom elevation) for each cell appears in the lower left hand corner of the screen. You can also drag the arrow over the cells shown in the cross section.
Now go through similar steps to set recharge rates and the leakance values, remembering to zone these parameters properly.
Boundary Conditions. Click on BC's. Note that the check mark is on "Constant Head/Conc". We want to set constant head cells in layer 1, column 6 to represent the river. Click on insert>digitize polyline. Position the arrow on the cell in row 1, column 6 of layer one and depress the left mouse button. Drag the arrow down column 6 across all the rows to row 23. You will see a red line appear. Double click on the cell in column 6, row 23. A screen will appear for the head at the beginning of the polyline. Enter the head at the boundary to be 100 and then click on the OK box. Another screen will appear for the end of the polyline. Enter the head to be 100 and again click on the OK box. If you make a mistake and want to start over, click on delete>clear all.
Now click on "No-flow". The screen will disappear. Go to the cube and move to layer 2. Click on BCs>insert. (Note that the check mark is now on no-flow.) You are now ready to enter the inactive cells as no flow cells in layer 2 using the window option. Repeat the process for layer 3. (The default boundary condition is no-flow.)
Save your data. Click on "Model>paths to models". Click on the browse box under working directory and choose the a: drive. Enter the working directory to be nopump. The input and output files generated by GW Vistas will be written to a:\nopump.
Save your results by clicking on file>save as. Choose a file name or root file name to identify this simulation. This name will be used to identify the master file created by GW Vistas, which has the extension .gwv. All of the files created during the simulation will use this root file name. For example if the root name you select is run1, GWVistas will create a master file called run1.gwv and a series of input files with that same root name (e.g., run1.bas, run1.bcf, run1.sip, etc.) and output files as well (e.g., run1.hds, run1.ddn, run1.cbb, run1.out).
Run the Model. You are now ready to run the model. Click on Model>MODFLOW>Packages. Change the root file name on this screen to whatever you have already selected (e.g., run1). Also make sure the output control box is checked and assigned a unit number of 22. Click on OK.
Now click on Model>MODFLOW>basic package. Enter a title for your simulation. Note that it is in this screen that you choose the time unit for the simulation, for which the default option is days. Click OK. Now click on the BCF package. Click on the box "compute leakance (VCONT)" and remove the check mark. (For this simulation we have input the correct values for leakance and do not want them to be re-calculated by the model.) Select “unconfined aquifer” for layer 1 and “confined aquifer” for layers 2 and 3.
Now click on the output control package. Remove the check mark on the disable printing box in the last line of the screen. Choose the head print format to be 9G13.6. Choose the wrap format. The drawdown print format can be 11G10.3.
Click on solver package. Review the information on this screen. There is no need to change anything.
Click on initial heads. Check the values for default heads and make sure they are set to some value that is reasonable (like 100).
Click on print formats. Note that you can select how much information is printed in the output file. No need to change this screen.
Create the input data files and run the model. Now click Model. Note that “Use MODFLOW” is checked at the bottom of the screen. Now click on Model>MODFLOW>create data sets. GW Vistas will create the data files and save them in the nopump directory. Display the error/warning file.
Now click on Model>MODFLOW>Run MODFLOW. GWVistas will run the model. A simulation screen will appear. Wait until the screen disappears and GW Vistas gives you a message that the simulation is complete. Output files will be written to the nopump directory.
When the "import model results" screen appears, check to make sure that the root file name for the head (hds) file and the cell-by-cell flow (cbb) file is correct and then place a check mark next to the cell-by-cell flow file (there should already be a check mark for the head file.) If the root file name isn't correct, click on the browse box and select the correct file name (e.g., run1.hds and run1.cbb).
You should now see a contour plot of the heads. Click on the layer number next to the cube to view the plots for layers 2 and 3.
Save your file again. Click on file>save as. Save your gwv file in the a:\nopump directory. You can now exit GwVistas by clicking the x at the very upper right hand corner of the screen. You may choose to take a break and come back to the lab later to complete the exercise, or you can continue with run 2 at this time.
Run 2. Steady-State Flow with Pumping (pump).
Import the files from Run 1. Insert your disk into the a: drive and access GW Vistas. Click on File>open. A screen will appear. Click on browse and choose the nopump directory that contains the gwv file you created previously (e.g., a:\nopump\run1.gwv). Click OK. After you exit this screen, your grid should appear.
Click on Model>paths to models. Enter the working directory to be a:\pump. The input and output files generated by GW Vistas as well as a new gwv file will be written to that directory.
Enter the pumping well. You now need to enter the well. Click the BCs option. Check "well" and then choose insert>single cell to insert the well at the correct location in the grid by positioning the arrow over the appropriate cell. When the well screen appears specify the correct flow rate for the well. Remember that discharge rates are entered with a minus sign (i.e., -50000.)
Save the file. Click on “File>save as” and enter a file name for the new gwv file that contains the well. You can call it a:\pump\run1.gwv or you may wish to change the root file name to run2 so that the file will be called a:\pump\run2.gwv. Click on OK.
Run the model. Go through all the steps as before, checking the MODFLOW packages. Look at the root file name in the Model>MODFLOW>packages window and the working directory in the Model>paths to models screen to make sure the appropriate names are entered. Also be sure that the well option in the Model>MODFLOW>packages window is checked and is assigned a unit number of 12. One important new item to change is in the initial heads window. Check the box marked "Set heads from file", then click on the browse box. Select the heads from run 1 (e.g., a:\nopump\run1.hds) as the starting heads for the run 2 simulation. This will allow GW Vistas to compute drawdown as the difference between the heads from run 1 and run 2. Remove the check from the disable printing of output in the output control screen.
Go through the steps as before to create the data files and run the model. When the import results screen appears, check to make sure that you will import results from the pump directory. If necessary, click on the browse button and choose the head file (hds), drawdown file (ddn) and the cell-by-cell flow file (cbb) from the pump directory. Check the boxes to import all three of these files.
Viewing the results. Click on plot>what to display. Check head contours and "display velocity vectors". The contour map with velocity vectors should appear.
Now change the contour interval. Click on Plot>contour>parameters (plan). Choose a 0.5 ft contour interval. You can also experiment with labeling intervals and font size. The add option allows you to add a title.
Now click on View>full>screen. Then click on File>page set-up and remove all checks from the clip windows box. Then click on file>print preview. If you like the way the plot looks, click on Print and collect your plot from the printer at the rear of the lab. You have just generated a plot for layer 1. Now move to layer 2 by clicking on the cube at the left hand side of the screen that shows the layers. Generate contour plots for layers 2 and 3 and complete question #1 on page 6 of this handout.
Save your gwv file. Save the gwv file by clicking on file>save as. Make sure you save the gvw file in the pump directory. You may now take a break and leave GW Vistas or continue.
Analyzing Results. If you have taken a break and left GWVistas, you will have to load up the gwv file from run 2 with the pumping well by clicking file>open and choosing the appropriate gwv file. You also need to load in the output from the model. Do this by clicking Plot>Import results. The Import results screen appears; click the hds, ddn and cbb boxes. The head contour plot should appear. Move the arrow over the grid and note that the head value appears in the lower right hand corner of the screen. Refer to Question #1 – see below.