Navigation:  Quick Start Tutorials >

PEST with Pilot Points

Print this Topic Previous page Return to chapter overview Next page

 

This exercise demonstrates some of the advanced and exciting opportunities for model calibration and uncertainty analysis using PEST in Visual MODFLOW Flex. This exercise is based on the problem described in "Using Pilot Points to Calibrate a MODFLOW/MT3D Model", by John Doherty (Watermark Numerical Computing), and has been adjusted to work with the PEST workflow inside Visual MODFLOW Flex.

 

Objectives

By the time you have finished this exercise you will have:

 

used pilot points as a means of characterizing the spatial distribution of an aquifer hydraulic property;

used PEST’s advanced regularization functionality in conjunction with geostatistically-based regularization constraints;

used the Visual MODFLOW Flex GUI to:

oBuild input files for PEST

oRun PEST

oAnalyze the results

oSave adjusted parameters as new model inputs

 

Before you start

You are encouraged to familiarize yourself with the concepts and applications of PEST prior to using in Visual MODFLOW Flex. The time spent on this will make your experience with parameter estimation much more productive, and will likely help you to overcome any difficulties you may experience the first time you run PEST.

 

In addition, if you are not familiar with the Visual MODFLOW Flex graphical environment, please take a few minutes to review the Program Overview section.

 

Note:

You must have a Pro or Premium license in order to use the PEST module in Visual MODFLOW Flex.

Required Files

Several files are required for this exercise, which should be included with the Visual MODFLOW Flex installation.  

 

These files are available in your public "My Documents" folder, typically:

C:\Users\Public\Documents\VMODFlex\Tutorials\PEST

 

If you cannot find these files, please download the PEST Tutorial project from our website, and unzip to a desired folder on your computer.

Introduction

Launch Visual MODFLOW Flex

Click [File]>[Open Project]

Navigate to your Public 'My Documents' folder, and locate '\Documents\VMODFlex\Tutorials\PEST', and open the 'pest-tutorial.amd' file  

[Open]

In the Numerical Workflow window, locate the workflow tree (this is shown on the left in the window below)

If you're not already on this step, click on the 'Select Run Type' in the numerical workflow tree.

Click on the [PEST Run] button

 

 

A new PEST Workflow window will load as shown below

       Note: if this does not appear right away, select this tab from the list of workflows across the top of the window, highlighted in                image below)

 

Define Observations

The first step is to define the observations you want to include for calculating the objective function, and assign weights to various observations.

In this example, there are 21 observation wells where heads were measured at several intervals.

All head locations and associated time varying heads are selected by default, and will be used with a default weight = 1.  

No changes will be needed for this exercise

Define Property Parameters

Click [] (Next Step) to proceed to the 'Define Property Parameters' step.

The 'Define Property Parameters' window will then appear:

 

 

At this step, select which parameters you want to include in the PEST run.

In the table at the top, select which property parameters you want to include; in the table at the bottom, select which property zones you want to include. For this exercise, all Kx property zones will be included (these are selected by default)

For each parameter, you can specify to "Tie" it to another parameter. You can also specify the Transformation option (by default, all Conductivity parameters are set to Log transformation)

In the table at the bottom, for each property zone, you must specify minimum and maximum values; defaults are provided.

In addition, the value from each zone is also displayed to assist in defining reasonable minimum/maximums.

Enter 1 for the Minimum for each zone

Enter 300 for the Maximum for each zone.

 

 

Click [] (Next Step) to proceed to the 'Define Pilot Points' step.

Define Pilot Points

The next step is to 'Define Pilot Points' as shown below

 

 

Pilot points are XY points with an initial value for each parameter. Pilot points can be imported from .TXT file, XLS, .SHP file, or assigned manually by digitizing in the 2D environment. An example of the pilot points on top of conductivity zones is shown below

 

 

The general steps are as follows:

 

Select 'pilot-points-zone1' data object from the Data Explorer (tree)

Click [] (Insert button) at the top of the 'Define Pilot Points' window to add these points.

Select which parameter zones the points represent under the Parameter Zones grid in the upper right section of the window.

For this set of pilot points, Kx-Zone1 is fine (this should be selected by default)

Your display should now appear as shown below:

 

 

Now repeat these steps to define the pilot points for Kx Zone2

 

Select 'pilot-points-zone2' data object from the Data Explorer (tree)

Click [] (Insert button) at the top of the 'Define Pilot Points' window to add these points.

Select Kx-Zone2 for these pilot points (if it is not already selected)

 

Repeat these steps for Kx Zone3 and Zone4

 

Select 'pilot-points-zone3' data object from the Data Explorer (tree)

Click [] (Insert button) at the top of the 'Define Pilot Points' window to add these points.

Select KxZone3 for these pilot points

 

Select 'pilot-points-zone4' data object from the Data Explorer (tree)

Click [] (Insert button) at the top of the 'Define Pilot Points' window to add these points.

Select Kx-Zone4 for these pilot points

 

When you are finished, your display should appear similar to the figure below.

 

 

In the table at the bottom of the window, you can adjust some parameters for specific pilot points. The main use of this table is to specify which pilot points (if any) are Fixed (Hard) and their initial values. Fixed pilot points are those locations where you are confident in the measured parameter value (eg. pumping/slug test locations), and you want these values to remain Fixed during the PEST run.

 

Above this table, there is a combo box where you select which parameter zone should be shown.  VMOD Flex allows you to combine multiple pilot point sets (eg. hard/soft) for parameter zones.

 

In this example, leave the defaults as is.

 

Click [] (Next Step) to proceed to the Define Kriging Variograms.

Define Kriging Parameters

The use of pilot points in characterizing the spatial distribution of a hydraulic property must be accompanied by a mechanism whereby hydraulic property values assigned to pilot points are spatially interpolated to the cells of the finite difference grid. Spatial interpolation is accomplished using the Kriging algorithm. Kriging is a method of spatial interpolation based on geostatistics. The cornerstone of geostatistics is the variogram. A variogram describes the extent to which hydraulic property values (or any other type of data) pertaining to any two points are likely to be different from each other as a function of the distance between those points.

 

One of the benefits of using Kriging as a basis for spatial interpolation is that the factors by which hydraulic properties at pilot points are multiplied before summation to obtain the hydraulic property value at a particular grid cell are independent of the actual hydraulic property values at the pilot points. Hence a set of “Kriging factors” pertaining to each of the cells of the finite difference grid can be calculated in advance of the actual interpolation process.Since the interpolation process is undertaken again and again as the model is run repeatedly by PEST, the fact that it is not necessary to repeat calculation of the Kriging factors on each occasion that the model is run can result in large savings in the time required to complete the overall parameter estimation process.

 

In this example, there are 4 hydraulic property zones as shown below:

 

 

Zone 1: river alluvium

Zone 2: creek alluvium

Zone 3: western basalt

Zone 4: eastern basalt

 

Each zone is represented by a geostatistical structure. Each of these structures cites one variogram (though it could cite up to five). “Structure1” will be used to characterize zone 1 of our model domain (i.e. the river alluvium), “Structure2” will be used to characterize zone 2 (ie. the creek alluvium) whereas “Structure3” will be used to characterize zones 3 and 4. Note that the variogram assigned to these latter zones is quite unimportant; because there is only one pilot point assigned to each of them, all cells within these zones will be assigned the one interpolated value (same as the respective pilot point) irrespective of the variogram.

 

In VMOD Flex, a default variogram is generated for each parameter zone, with type 2 (Exponential).  

 

 

However, the structures should be modified.

 

Click 'Zone1' under 'Kx' in the tree. The display will appear as below:

 

 

 

In the editor on the right side of the display, for 'Transform', set this to log.

Repeat these steps, setting log 'Transform' for Kx in Zones 2, 3, and 4.

 

Thus any variogram cited in each of these structures must pertain to the spatial distribution of the logarithm of the pertinent hydraulic property. This is in accordance with the fact that most studies cited in the groundwater literature which treat transmissivity and/or hydraulic conductivity as a regionalised variable indicate that its distribution is better described by a log variogram than by a variogram based on native property values.

 

For this example, the default parameters for the variograms (all will use Exponential) is sufficient. However, for Zone2, we will define a value of 2.0 for the Anisotropy, with the direction of anisotropy coinciding with the direction of the creek. Alignment in the direction of the creek is based on the premise that channel structures within this old creek valley will make it more likely for hydraulic property similarity to prevail in this direction than in a direction at right angles to it.

 

The variogram parameters can be adjusted by selecting a Variogram from the tree as shown below.

 

Click 'Variogramkx2' under the Variograms node in the tree.

Enter 2 for the 'Anisotropy'.

Enter 45 for the 'Bearing'. This value allows you to make the variogram anisotropic in a certain direction; 'Bearing' is an angle of rotation.

The display should appear as below:

 

 

Click [] (Next Step) to proceed to the Select Run Type.

Select Run Type

At this step, choose the type of PEST Run. If you want to run PEST, then some additional options will be needed, such as define regularization and adjust the PEST control file.  If you want to run Sensitivity Analysis, this can also be launched.

 

 

Click the [Parameter Estimation] button from the main window, as shown above.

The next step will be to choose the Regularization options.

Select Regularization

At this step, choose the type of Regularization to run.

 

 

Click [No Regularization] button from the main window

The next step will be to adjust the PEST control file.

Edit PEST Run Settings

The last step before running PEST is to view and adjust the PEST Control file.

 

 

If you are familiar with the PEST file format/structure, you can adjust the PEST Control file in this window, or copy into a text editor, make changes, and paste the adjusted contents back in this window.  A full explanation of the PEST control file is available in the PEST manual
(http://www.pesthomepage.org/getfiles.php?file=pestman.pdf)

 

For this exercise, the default values are fine.

 

Before starting the PEST Run, it is a good idea to check the PEST Input files.  PEST provides a utility to do this, called PESTCheck.

 

Click the [] button on the workflow toolbar

You should receive a confirmation that no errors were found.

Click [OK] to proceed.

Click [] (Next Step) to proceed to the Run PEST window.

 

Run PEST

Click [] (Run PEST) to start the PEST Run.

 

The pest.exe will load in a DOS command window, and show the progress as seen below.

 

 

Depending on the speed of your computer, the PEST run should take between 3-5 minutes. As PEST runs, you should see the objective function (phi) decreasing over each optimization iteration; pay attention to these values in the DOS window. PEST will run a total of 25 optimization iterations and a total of 1366 model runs. PEST should reach a final objective function (phi) value of approximately 2.06 E-02.

 

When PEST finishes, you should see a confirmation message in the main window, below the 'PEST Run Log' tab, as shown below.

 

 

After the PEST run completes, you can analyze the results.

 

Click [] (Next Step) to proceed to the 'Analyze Results' step.

Analyze Results

VMOD Flex presents the results of the PEST run, with one tab per output file.

 

 

Record file (.REC): contains parameter values, objective function, sensitivities, etc..

Sensitivities for Observations (.SEO): contains the observed and simulated values with the sensitivities  

Sensitivities for Parameters (.SEN): provides the information on parameter sensitivities

Residuals (.RES): contains the adjusted calculated vs. observed values and residuals

 

The results from these files can be exported into Excel for charting.

 

If the results look reasonable, you can save the adjusted Kx parameter zonation as inputs for a new model; this is explained in the next section.

 

Click [] (Next Step) to proceed.

Save PEST Parameters as New Inputs

After reviewing the PEST output, if the adjusted parameter values seem reasonable, you can save these parameters as inputs for a new model run.

 

 

 

Click on the [Create New Model Run with PEST Results] button.

 

VMOD Flex will save the adjusted model parameters in a new model run within the same project. This new model run will appear in the Model Explorer below the most recent model run. A new workflow window will also appear with this model run.

 

 

Click on 'Define Properties' from the workflow tree (if not already selected)

From the Toolbox, select 'Kx' as shown below; you should then see a color flood of the Kx values.

.

 

You can mouse over the 2D display to see the range of Kx values. Click on the [] button on the toolbar to show color shading with contour lines.

 

You must Translate and Run this new model run in order to see the updated MODFLOW results using the adjusted Kx parameters from PEST.

 

*****This concludes the 'Model Calibration Using PEST with Pilot Points' tutorial.*****

 


Page url: http://www.waterloohydrogeologic.com/help/vmod-flex/index.html?vm_pest_tutorial.htm