27th Congress of the International Association for Hydraulic Reseach (ASCE),
San Francisco, California, August 10-15, 1997.
Water for a Changing Global Community, Groundwater: An Endangered Resource, Proceedings of Theme C
1 Colorado School of Mines, Golden, Colorado, USA
2 U.S. Army Waterways Experiment Station, Vicksburg, Mississippi, USA
New geostatistical simulation tools that allow better definition of hydrofacies were coupled with previous techniques to reduce uncertainty associated with estimates of hydraulic parameters. Confidence intervals on the estimated hydraulic conductivities are 20 to 36 percent smaller when the new simulation tools are used.
Uncertainty concerning the subsurface arrangement of hydrofacies precludes accurate prediction of groundwater flow and contaminant transport. Combined use of disparate types of data reduces the uncertainty associated with solutions to groundwater flow and transport problems relative to the uncertainty level attained using only one type of data or only hard data. Hard data are data with negligible uncertainty, such as direct observations of lithologic units in core, while soft data are data with non-negligible uncertainty such as indirect measurement gathered in a geophysical survey, or expert opinion regarding geologic trends. In previous work (McKenna and Poeter, 1994), uncertainty associated with hydraulic parameter estimation of units within the Fountain and Lyons Formations was reduced by: 1) generating stochastic realizations of hydrofacies geometries conditioned to both hard and soft data, incorporating both geologic observations and geophysical survey results, through multiple-indicator, conditional, geostatistical simulation; 2) eliminating, hydrofacies realizations that were determined to be implausible based on results of inverse flow modeling; and 3) estimating hydrologic parameters using inverse flow modeling. Use of soft data produced a substantially different hydrofacies distribution, which was believed to be more representative of the site, than those produced with hard data alone, because forward groundwater flow modeling, using estimates of hydraulic conductivity (K) from field testing, yielded smaller head residuals for the "soft" realizations.
In this paper, uncertainty is further reduced for the same data set by utilizing new geostatistical simulation tools (Wingle and Poeter, 1996, in preparation a,b) that allow better definition of units through alleviation of stationarity constraints and semivariogram model limitations, and through use of a class, versus threshold, geostatistical simulation technique.
The site is located on the Colorado School of Mines campus in Golden CO. Site data include: core and chip data from eighteen boreholes; borehole geophysical logs; eight, two-dimensional, cross-hole, seismic tomogram sections; geologic knowledge, hydraulic test data, and observations of hydraulic head. The site includes the Fountain and Lyons Formations which consist of channel and overbank deposits that have undergone variable diagenesis, resulting in more hydrofacies than would be encountered in unconsolidated deposits.
The disparate types of available data are integrated to yield a coherent hydrofacies classification through use of discriminant analysis and soft data geostatistical simulation techniques, improving definition of the complex hydrofacies and increasing knowledge of their spatial correlation. Indicator values represent hydrofacies with initial hydraulic conductivity values (K) estimated from field and laboratory tests (Table 1) or estimated using inverse flow modeling.
|1||conglomerate of Lyons Formation||1x10-3a|
|2||well cemented, poorly sorted, sandstone w/conglom lense||1x10-3b|
|3||sandstone of Lyons Formation||1x10-2a|
|4||well cemented, poorly sorted, sandstone||5x10-3c|
|5||moderately consolidated, porous material w/some clay content||5x10-6a|
|6||combined well sorted, silty sandstone; poorly consol, poorly|
sorted sandstone; w/ fine gravel conglom and siltstone lenses
|7||weathered poorly consol, moderately porous, clay poor sandsy||1x10-3a|
|8||well fractured area of any previous hydrofacies||8e|
|a||no data, value is inferred from lithologic character|
|b||based on two packer tests|
|c||based on one packer tests|
|d||based on seven packer tests|
|e||based on one packer tests|
Although the site includes two geologic formations, it was initially simulated as a single zone (McKenna and Poeter, 1994) because a zonal simulation tool was not available. In this study, zonal kriging (Wingle and Poeter, 1996) is used in combination with independent directional semivariograms (Wingle and Poeter, in prep a) and class indicator simulation (Wingle and Poeter, in prep b) to model the CSM survey field. The boundary between the zones was clearly delineated on the seismic tomograms as areas east (Lyons) and west (Fountain) of the formation contact that dips 40o ENE. Zonal kriging is needed because the data populations and their spatial statistics are different in the two formations. For example, semivariograms for indicator 4 in the Fountain and Lyons Formation have ranges of 72 and 18 feet in the North-South direction, respectively, while in the East-West direction the ranges are; 27 and 48 feet, respectively. Not only are the ranges for the same indicator substantially different, but the principle anisotropy directions are 90oapart.
Geostatistical simulations were conducted on a regular grid of 81 columns representing 405 feet in the X direction, 61 rows representing 305 feet in the Y direction, and 72 layers representing 144 feet in the Z direction, with grid blocks 5x5x2 feet. Eighty multiple indicator conditional stochastic simulations of the site are generated, forty using the previous simulation methods and another forty using the new simulation tools. The only difference between the two sets of realizations are due to the different zone definitions and their associated semivariogram models.
More information about the software is available on the following web site: http://www.mines.edu/fs_home/wwingle/uncert/. For specific, updated details on access to, and use of the software, communicate with: email@example.com.
INVERSE FLOW MODELING
The groundwater gradient is essentially west to east across the site, thus fixed head boundaries are imposed on the west and east boundaries with zero flux boundaries along the north and south to simulate easterly flow. Recharge in the vicinity of the site (Kiusalaas, 1992) is essentially zero. The full scale geostatistical realizations of 355,752 elements are too large to solve the inverse flow model for multiple realizations with the available facilities in a reasonable time frame, thus the geostatistical realizations were upscaled to cell dimensions of 20x20x2 feet (x, y, z), incorporating 16 of the fine-scale cells into one coarse cell and reducing the number of cells by more than 90%.
Parameter estimation (i.e. calibration) is accomplished with a non-linear regression algorithm implemented in MODFLOWP (Hill, 1992). The K's of hydrofacies 6 and 8 were known with greater confidence relative to other hydrofacies, therefore they were not estimated. Sensitivity to K of hydrofacies 5 was low, consequently estimated values did not deviate from the prior estimation, thus for all intents and purposes the K of hydrofacies 5 was not estimated. "Acceptable" realizations were selected based on the following criteria: 1) the regression must converge; 2) hydrofacies 5 (moderately consolidated, clayey material) must have the lowest K of all the hydrofacies; 3) hydrofacies 8 (fracture zone) must have the highest estimated K of any hydrofacies; and 4) the "best fits" are selected by, after applying criteria 1 through 3, retaining the half of the solutions remaining that have the smallest sum-of-squared-head residuals. The first rule may be most controversial due to the control that a model user has over convergence parameters. In this case liberal parameters (e.g. a large number of iterations, loose tolerance) are applied to quickly find zonations that deviate substantially from the actual field conditions. Earlier studies indicate that realizations eliminated by this criteria do deviate from actual zonations (Poeter and McKenna, in review).
All realizations yielded values of K higher than hydrofacies 5 and lower than hydrofacies 8, for all of the other hydrofacies. Seven percent of the single zone and 22 percent of the two zone realizations failed to converge on stable parameter values. After parameter estimation, weighted head residuals were lower for the two-zone simulations (mean 0.15 ft and standard deviation 0.11 ft), relative to the single-zone simulations (mean 0.22 ft and standard deviation 0.14 ft). Although, estimates of K from two-zone realizations are similar to estimates obtained using single-zone realizations, the confidence intervals are 20 to 36 percent smaller (Figure 1).
Insensitivity to K of hydrofacies 5, coupled with the prior information provided on the parameter value caused the regression to estimate a value equal to the prior information for every realization, thus we treat K of hydrofacies 5 as if it were a "fixed" value. The mean K value of hydrofacies 4 for both ensembles of realizations was approximately an order of magnitude lower than the value estimated from field data, consequently we conclude that the one available measurement in hydrofacies 4 was not representative of the effective K value for that unit. The remaining estimated K values (hydrofacies 1,2,3, and 7) are closer to the field estimates (estimated K are within 5-30% of an order of magnitude from the value determined by packer tests) when the two-zone geostatistical realizations are used in the flow model, than when the single-zone realizations are used (within 40-100% of an order of magnitude from the field estimates). Of the hydrofacies for which K values were estimated, the estimated K's for hydrofacies 2 were closest to those measured in the field, regardless of which geostatistical ensemble was utilized, suggesting the two packer tests conducted in hydrofacies two reflected the effective hydraulic conductivity of that unit.
To make more accurate estimates of the parameters it is desirable to have observations of both heads and flows, such that absolute values (rather than relative values) of K are estimated. In this study, only 28 head observations are available. Flow rates are not available at this site, so the solution is constrained by the assumed knowledge of K for hydrofacies 5, 6, and 8. The estimates of K for the remaining hydrofacies are calculated relative to the values of hydrofacies 5, 6, and 8.
Use of geostatistical simulation tools that allow zonal, directional, and class indicator kriging were coupled with established soft-data simulation and inverse-model culling techniques to reduce uncertainty associated with estimates of hydraulic parameters below previously established levels using soft-data simulation and inverse-model culling techniques alone. Head residuals were reduced 32% when two zones were used to generate geostatistical realizations of the hydrofacies units. Confidence intervals on the estimated hydraulic conductivities are 20 to 36 percent smaller when the new simulation tools are used.
This work results from research conducted under the Groundwater Modeling Program sponsored by the Army Environmental Center. Permission was granted by the Chief of Engineers, U.S. Army Corps of Engineers, to publish this information.
Hill, M.C., 1992, A computer program (MODFLOWP) for estimating parameters of a transient, three-dimensional, groundwater flow model using nonlinear regression, United States Geological Survey, Open File Report 91-484, 358 pp.
Kiusalaas, N.J., 1992, Estimation of Groundwater Recharge Using Neutron Probe Moisture Readings Near Golden, CO, Colorado School of Mines ME Thesis ER- 4191, 186 pp.
McKenna, S.A. and E.P. Poeter, 1994, Simulating Geological Uncertainty with Imprecise Data for Groundwater Flow and Advective Transport Modeling in: eds. Yarus J.M. and R.L. Chambers: Stochastic Modeling and Geostatistics: Case Histories and Practical Examples, AAPG Special Publication, Computer Application in Geology, no. 3.
McKenna S.A. and E.P. Poeter, 1995, Field example of data fusion in site characterization, Water Resources Research, vol. 31, no. 12, pp 3229-3240.
Poeter, E.P. and S.A. McKenna, invited chapter in review, Combination of Geologic Information and Inverse Parameter Estimation for Improving Groundwater Modeling, Chapter of SEPM Special Publication: Concepts in Hydrology and Environmental Geology.
Wingle, W.L. and E.P. Poeter, 1996, Evaluating subsurface uncertainty using conditional indicator simulation and zonal kriging, Proceedings of Uncertainty '96 (ASCE), Madison Wisconsin.
Wingle, W.L. and E.P. Poeter, in preparation a, Independent directions semivariogram models versus anisotropic factors in kriging.
Wingle, W.L. and E.P. Poeter, in preparation b, Class versus threshold approach to indicator kriging.