Published in Vadose Zone Journal 3:59-74 (2004)
© 2004 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA
SPECIAL SECTION: UNDERSTANDING SUBSURFACE FLOW AND TRANSPORT PROCESSES AT THE IDAHO NATIONAL ENGINEERING & ENVIRONMENTAL LABORATORY (INEEL) SITE
Regulatory Modeling for the Idaho National Engineering and Environmental Laboratory's Subsurface Disposal Area and Conceptual Model Uncertainty Treatment
Swen Magnuson*
Bechtel BWXT Idaho, LLC, P.O. Box 1625, MS 2107, Idaho Falls, ID 83415
* Corresponding author (smm{at}inel.gov).
Received 4 February 2003.
 |
ABSTRACT
|
|---|
The Subsurface Disposal Area (SDA) at the Idaho National Engineering and Environmental Laboratory (INEEL) contains low-level, hazardous, and transuranic wastes and is the subject of ongoing environmental assessments. The buried waste is in shallow land burial overlying a 180-m-thick vadose zone comprised of fractured basalts and sedimentary interbeds. Applied simulation studies to predict future contaminant concentrations from this buried waste have been conducted as part of ongoing regulatory-driven environmental assessments for the SDA. An overview of the most recent numerical simulation study is presented along with a treatment of conceptual model uncertainty. The main premise of this paper is that conceptual uncertainty likely dominates parametric uncertainty, and this dominance should be considered when designing characterization studies. Given that the current modeling is uncalibrated due to inconsistent monitoring data in the vadose zone, it is only appropriate to consider conceptual model uncertainty. Conceptual model uncertainty was addressed by modeling different scenarios. Scenarios presented in this paper are the assignment of infiltration boundary conditions, a simplistic approximation to represent facilitated transport, and alternative representations of the moisture characteristic curve for fractured basalts. Results of this subsurface transport modeling have been used in part to guide additional field characterization activities. Eventually, predicted concentrations from this applied model will be used to estimate human health risks for a potential future receptor as part of the environmental assessment process.
Abbreviations: CERCLA, Comprehensive Environmental Response, Compensation, and Liability Act • INEEL, Idaho National Engineering and Environmental Laboratory • SDA, Subsurface Disposal Area
|
INTRODUCTION
|
|---|
THE SUBSURFACE DISPOSAL AREA was established in 1952 and contains low-level, hazardous, and transuranic wastes in shallow land burial in the southwestern portion of the INEEL. As part of an environmental assessment for an ongoing Comprehensive Environmental Response, Compensation, and Liability Act (CERCLA) evaluation, contaminant migration through the hydrologically complex subsurface beneath the SDA was simulated. Holdren et al. (2002) presented a model for simulating flow and transport through the vadose zone beneath the SDA. An overview of this recent numerical simulation study is presented along with a treatment of conceptual model uncertainty.
The environmental assessment for the SDA is driven by regulatory requirements. The modeling performed for this assessment strives to have a technical basis but recognizes limitations imposed by uncertainties. Uncertainty in modeling results is inherent, especially in applied simulation studies. The uncertainty in simulation results derives from two primary sources: parametric uncertainty and conceptual uncertainty. Parametric uncertainty results from the spatial variability in hydrologic and transport properties and an inability to characterize this variability perfectly, or even if it were perfectly characterized, an inability to describe this variability within a numerical model. Conceptual uncertainty results from inadequate knowledge of whether correct processes or external influences in the real system are appropriately simulated in the model. Examples of conceptual uncertainty would be appropriateness of a boundary condition such as the timing and amount of water applied at land surface or the appropriateness of the constitutive (relative permeability–saturation–pressure) relationships used to describe variably saturated water movement. In regulatory modeling applications, uncertainty in predicted future concentrations and their associated human health risks affects confidence in the ability to make remedial decisions. In lieu of accounting explicitly for uncertainty, conservative assumptions and parameterization of hydrologic and transport properties are often used that maximize the movement of contaminants toward potential receptors.
Field characterization activities are often performed with the goal of reducing uncertainty in modeling applications. Conceptual uncertainty likely has a greater impact on uncertainty in flow and transport simulations than parametric uncertainty (Neuman, 2002). This is not necessarily a new topic. Carrera (1993) also addressed conceptual uncertainty and its influence on transport simulation and the difficulty of quantifying its effect. However, the topic is worthy of revisiting as applied modeling studies are conducted frequently. The underlying premise in this paper is that characterization activities should be oriented first at reducing conceptual uncertainty before targeting parametric uncertainty. This simulation study investigates the impact of including selected processes and external influences on flow and transport phenomena at the SDA as part of investigating conceptual model uncertainty. Specifically addressed were the use of spatial varying infiltration rates, facilitated transport of actinides, and constitutive relations used to describe water movement in fractured basalt. Although the simulations included flow and transport in both the vadose zone and within the underlying aquifer, the emphasis in this paper is on the vadose zone.
|
OVERVIEW OF SITE GEOLOGY AND CONCEPTUAL MODEL
|
|---|
The SDA is 180 m above the regional Snake River Plain Aquifer. The vadose zone beneath the SDA is composed of fractured basalts and unconsolidated sediments. Thick layers of fractured basalt are interspersed with thin, mostly continuous, sedimentary interbeds that were deposited through aeolian and fluvial processes during periods of volcanic quiescence. Figure 1
shows the lithologic positions and nomenclature by which the interbeds are identified. The conceptual model for the subsurface consists of a heterogeneous porous media that has both isotropic and anisotropic structures. Surficial sediments and sedimentary interbeds are modeled as isotropic with varying thickness and elevations. Three sedimentary interbeds are included in the conceptual model and are termed the AB, BC, and CD interbeds following the nomenclature of Anderson and Lewis (1989). These interbeds occur at approximate depths of 10, 34, and 70 m below land surface. These three interbeds are the best characterized, with the CD interbed being the most continuous in the vicinity of the SDA. The deeper interbeds, although present, are less characterized because of fewer wells to those depths. These deeper interbeds are conservatively neglected.
The surficial sediments and AB interbed are modeled as homogeneous and isotropic. The BC and CD interbeds are modeled as heterogeneous and isotropic. The BC and CD interbeds, although deeper than the AB interbed, are better characterized because of efforts to evaluate possible preferential pathways across the interbeds because of spatial variability of hydrologic and transport properties.
Fractured basalts are modeled as homogenous equivalent porous media with a low effective porosity and an anisotropic high permeability compared with the sediments. The basalt matrix is treated as if it were impermeable, with flow only occurring in the fracture network. This treatment of the fractured basalt resulted in fast movement of water controlled primarily by gravity, with no propensity for lateral movement until an interbed was encountered.
Sources of water that affect contaminant behavior in the vadose zone are infiltration of meteoric water and occasional lateral recharge from the nearby Big Lost River system (see Fig. 2)
during periods of high flow. The Big Lost River is approximately 3 km northwest of the SDA and flows sporadically depending on regional snowpack accumulation. Excess runoff water from the Big Lost River can optionally be diverted to isolated playas that serve as spreading areas. The diversion system was constructed in 1958 to provide flood protection for facilities downstream on the INEEL. These spreading areas range in distance from 1 to 1.5 km southwest and west of the SDA. Nimmo et al. (2002) conducted a tracer study within the spreading areas during a Big Lost River runoff event and observed that spreading area water migrated laterally within the vadose zone to a location beneath the SDA. This source of lateral recharge is included in the simulations and is approximated as a constant lateral flux assigned at depth above the CD interbed. Infiltration of meteoric water is modeled as being spatially variable and dependent on topography and degree of soil disturbance. This spatially variable pattern of infiltration was assumed to take effect across the entire SDA at the beginning of waste emplacement in 1952. Currently only a soil cover exists over the SDA to restrict infiltration of meteoric water, so the assigned spatially variable infiltration is considered to continue indefinitely into the future. In addition to these two primary sources of water, there have been three occasions where runoff from the local basin in which the SDA is positioned entered onto and flooded portions of the SDA. These occurred in 1962, 1969, and 1982, and their estimated magnitudes of additional infiltrating water were included in the simulations.
View larger version (26K):
[in this window]
[in a new window]
|
Fig. 2. Location of the Big Lost River and spreading areas in relationship to the Subsurface Disposal Area (SDA).
|
|
Movement of water and contaminants locally within the aquifer is controlled by the regional flow in the aquifer, with occasional fluctuations in flow direction caused by the Big Lost River discharges to the spreading areas. However, flow within the aquifer was simulated as steady state given the long duration of the simulations (i.e., hundreds to thousands of years) and the assumption that the Big Lost River diversion would not be maintained indefinitely.
Processes included in the conceptual model for dissolved-phase transport were advection, dispersion, diffusion, adsorption, and radioactive decay. Given the relatively low surface area of the fractures and the fast movement of water through them, it was assumed that sorption was negligible in the fractured basalt. Sorption was included only in the surficial and interbed sediments. Linear reversible isotherms were used to describe sorption onto the unconsolidated sediments. The partition coefficients for all contaminants were based on best-estimate values developed by Dicke (1997). Facilitated transport mechanisms (e.g., colloidal transport) may affect contaminant migration in the SDA subsurface. Column studies with SDA interbed sediments have been conducted that show very small fractions of the Pu and Am may move in a facilitated manner at the SDA (Fjeld et al., 2000). Facilitated transport was included conceptually by treating small fractions of waste as mobile.
|
SIMULATOR DESCRIPTION
|
|---|
The simulations were performed using the TETRAD simulator (Vinsome and Shook, 1993; Shook, 1995). The TETRAD code uses a block-centered finite-difference numerical approach to solve the governing equations for conservation of mass and momentum. TETRAD has multiphase, multicomponent simulation capabilities and can mimic the behavior of any number of components in aqueous, gaseous, oleic, and solid phases. The modeling in this study was limited to dissolved phase aqueous transport. TETRAD has capabilities for local grid refinement, which were used to match the level of discretization in the model to available information. Since there are more shallow wells than deep wells at the SDA, this resulted in local grid refinement being implemented in the upper portions of the vadose zone model.
The general conservation equation solved by the TETRAD simulator for accumulation, flux, decay, or degradation, and sources for any component i can be written as
 | [1] |
where 
Wi/t is the accumulation term that consists of net changes in concentration of component i in any phase including the adsorbed phase, 
Ni is the flux of component i, Ri is the change in concentration arising from decay of component i, and qi represents sources or sinks of component i.
The accumulation term can be written as
 | [2] |
where 
is the porosity; Sj are the phase saturations ("w" aqueous, "g" gaseous, "o" oleic); 
j are the phase molar densities; wi, yi, and xi are the mole fractions of i in the aqueous, gaseous, and oleic phases, respectively; s is the solid phase density; Vi is the mole fraction of i adsorbed on the solid phase, and Mi is the molecular weight of i. A generalized adsorption relationship is available in TETRAD that allows for adsorption onto the solid phase from any of the other three phases.
The flux term in Eq. [1] is comprised of an advection and dispersion term for each phase given by
 | [3] |
In Eq. [3] the uj are the phase advective fluxes, given by the multiphase version of Darcy's Law:
 | [4] |
where k is the intrinsic permeability, krj is the relative phase permeability, µj is the phase viscosity, Pj is the phase pressure, and g is gravitational constant in vector form. Dij is the phase-dependant dispersion tensor comprised of molecular diffusion modified by porosity, phase saturation, and tortuosity and mechanical dispersion consisting of phase dispersivities modified by directional components of advective phase fluxes (Bear, 1972).
The reaction term in Eq. [1] accounts for decay or degradation of component i and is written as
 | [5] |
where mi is the total aqueous mass of i. The first term on the right-hand side of Eq. [5] accounts for i decaying with a rate constant Ai
into component , whereas the second term on the right-hand side is the formation of i from destruction of component 
with a rate constant Ai.
The final term on the right-hand side in Eq. [1] is the source–sink term, qi. This term accounts for the addition or extraction of component i through wells or boundary conditions.
Two analytic expressions for relationships between capillary pressure or water potential, moisture content, and relative permeability were utilized in these simulations. For sediments, van Genuchten (1980) equations were used to describe the relationship between water potential and moisture content. The Mualem (1976) equation was used to describe the relationship between saturation and relative permeability. The three-phase van Genuchten constitutive equations, as adapted by Parker et al. (1987) with slight modifications to the normalized saturation terms, implemented in TETRAD were
 | [6] |
and
 | [7] |
where Pc_ow and Pc_go are the interphase capillary pressure between the oleic-aqueous and gaseous-oleic phases, respectively; 
ij is the surface tension between phases i and j; w is the density of water; 
is the van Genuchten fitting parameter related to the inverse air-entry pressure; ß is the van Genuchten fitting parameter related to the pore-size distribution, and 
= 1 – 1/ß. The normalized saturation terms in Eq. [6] and [7] are given by 
w = (Sw – Swr)/(1 – Swr) and g = Sg/(1 – Swr), where Swr is the irreducible aqueous-phase saturation. In this simulation study there was never any oleic phase present so go was set to zero, and ow and gw were set to one. The effect was to make Eq. [6] represent the capillary pressure between the aqueous and gaseous phases. This capillary pressure is the same as the water potential in the absence of osmotic pressure. Water potential expressed as centimeters of water will be used for the remainder of this paper.
The relative permeability terms for the aqueous and gaseous phases implemented in TETRAD are
 | [8] |
and
 | [9] |
where Sgr is the irreducible gaseous-phase saturation.
Determining constitutive relationships for fractures continues to be an active research area. A Corey-type analytical formulation was used to describe the constitutive relationships for the fractures. The formulae implemented in TETRAD are Pc_gw = Agw(1 – Sw)Bgw and kri = Ai[(Si – Sir)/(1 – Sir)]Bi, where Agw, Bgw, and Bi are fitting parameters.
|
BASE MODEL DEVELOPMENT
|
|---|
A three-dimensional domain was established for the vadose zone model. A plan view showing the horizontal extent of the vadose zone simulation domain is shown in Fig. 3
. This domain was discretized as shown in Fig. 4
using conformable gridding to match the variable surface topography and the variable interfaces between kriged lithologic contacts (Leecaster, 2002). Two locally refined zones were used to further discretize the model to match varying levels of available lithologic data on contacts between basalt flows and sediments. These refined zones nest within the base model and represent progressively shallower regions of the domain, matching the density of available well control information. The outline of the SDA is shown in Fig. 4 at the top of the base and refined domains for orientation. The conformable gridding approach resulted in variably sized vertical dimensions for all grid blocks. This allowed refinement of vertical grid block dimensions at lithologic interfaces in the model.
View larger version (47K):
[in this window]
[in a new window]
|
Fig. 4. Conformable grid for (A) base domain and (B and C) locally refined domains. Surficial sediments are shown in yellow, interbeds in green, and fractured basalt in blue.
|
|
Hydrologic properties for the sedimentary lithologic units were derived primarily from site-specific characterization and inverse modeling studies (McElroy and Hubbell, 1990; Baca et al., 1992; Bhorgese, 1988; Martian, 1995; Magnuson, 1995). Table 1 shows the assigned homogenous properties for the surficial sediment and the AB interbed. Porosity and permeabilities for the BC and CD sedimentary interbeds were spatially variable and were kriged based on interbed samples that were hydrologically characterized (Leecaster, 2002). The spatially variable porosity and permeability for the BC interbed are shown in Fig. 5
. The model also includes a low-porosity, low-permeability feature at the top of both the BC and CD interbeds. This is partly based on observations of a clay layer that is sometimes encountered at these depths during drilling and was also necessary to increase simulated water contents and water potentials into the observed range. The parameters for the van Genuchten constitutive equations for all the sediments in the model were = 1.066 m–1, ß = 1.523, and Swr = 0.142. These values were obtained from averaging the results of characterization activities on core samples from the surficial sediments (Baca et al., 1992).
View this table:
[in this window]
[in a new window]
|
Table 1. Hydrologic properties and source for surficial sediments, AB interbed, and the imposed low porosity and low permeability feature at the top of the BC and CD interbeds.
|
|
View larger version (46K):
[in this window]
[in a new window]
|
Fig. 5. Kriged (A) porosity and (B) permeability for most refined grid for the BC interbed with core sample values posted at their well locations.
|
|
The hydrologic properties for the fractured basalt portion of the vadose zone were taken from Magnuson (1995). In that study, an inverse model was applied to monitoring results from a large-scale infiltration test conducted approximately 1.5 km south of the SDA (Dunnivant et al., 1998). The empirical Corey-type parameters derived from Magnuson (1995) are given in Table 2.
Infiltration at land surface within the SDA was treated as spatially variable and was a function of topography and degree of soil disturbance. Three constant rates representing areas with high (24.1 cm yr–1), medium (3.68 cm yr–1), and low (0.64 cm yr–1) annual infiltration rates were assigned as shown in Fig. 6
. The equivalent average infiltration calculated from these spatially variable rates was 8.5 cm yr–1. These spatially varying estimates were based on inverse modeling to site-specific calibrated thermal neutron monitoring in a network of neutron access tubes (Martian, 1995). These spatially varying rates take into account the effect of topography, soil disturbance, and vegetation, with topography being the dominant factor. A uniform low infiltration rate of 1 cm yr–1 annually was assigned outside the SDA (Cecil et al., 1992).
View larger version (44K):
[in this window]
[in a new window]
|
Fig. 6. (A) Spatially varying infiltration rates assigned inside the SDA and (B) assignment locations for simulated waste stream source releases.
|
|
The influence of Big Lost River water discharged to the spreading areas was included in the model by adding additional water at depth above the CD interbed. This was necessary since the horizontal extent of the vadose zone simulation domain did not include the spreading areas. The additional spreading area water was introduced at a constant rate beginning in 1965, the first year Big Lost River water was diverted to the spreading areas. The amount of additional water was adjusted so that the western half of the CD interbed under the SDA had increased moisture contents compared with simulations without the additional spreading area water. As an example of simulation results, Fig. 7
shows the maximum simulated water saturation at each location in both the BC and CD interbeds at a point in time before the introduction of the additional water above the CD interbed. The variation in saturations across the interbeds shows the extent to which infiltrating water gets focused at depth in the simulation.
View larger version (39K):
[in this window]
[in a new window]
|
Fig. 7. Simulated maximum water saturation in the (A) BC interbed and (B) CD interbed immediately before the introduction of the additional spreading area water above the CD interbed.
|
|
Dissolved-phase transport of a suite of contaminants of potential concern was simulated with this vadose zone model. The flow and transport simulations in this paper utilized a separate source release model (Holdren et al., 2002) that estimated time-dependent release of contaminants from buried locations. This source release model has a separate suite of uncertainty issues that are not addressed here. The source release model does utilize the spatial distribution of water indicated in Fig. 6. The infiltration rates assigned to the individual grid blocks representing each of the waste streams are averaged to obtain a single value for each waste stream. The resulting infiltration for the waste streams ranges from 0.6 cm yr–1 for Pad A to 11.7 cm yr–1 for Trenches 1 through 10. Mass fluxes at the bottom of the vadose zone model were input to an aquifer model. Simulated aquifer concentrations were then used to estimate hypothetical risk to a person using water from the aquifer in the vicinity of the SDA. Although there are uncertainties related to the aquifer model implementation, they are not discussed in this paper because the focus is on the vadose zone aspects of the model.
Although there are some contaminant concentrations from analysis of water samples from the vadose zone monitoring network, including recently installed monitoring points within and beneath the buried waste, consistent detections and trends are not evident (Olson et al., 2003). This lack of consistent trends precluded calibrating the simulations to observed contaminant advance. Therefore, the current flow and transport model is based on best judgment using information on lithology and hydrologic properties. Without calibration, the best use of the simulator is to test conceptual uncertainty, which is the topic of the remainder of this paper.
|
CONCEPTUAL UNCERTAINTY SIMULATION RESULTS
|
|---|
Given the limited characterization and field monitoring information on which this model is derived there are numerous conceptual uncertainties that could be considered. Addressed in this paper are spatial variability in infiltration, facilitated transport of contaminants, and hydraulic properties of fractures. The conceptual uncertainties investigated in this paper represent a wide range of impacts on the simulation results. Holdren et al. (2002) discussed additional conceptual uncertainties, such as the impact of the additional water at depth in the model from the spreading areas, the effect of possible gaps in the sedimentary interbeds, and the impact of variations in assigned partition coefficients.
Infiltration Rates
Infiltration of meteoric water with subsequent percolation downwards through the vadose zone is the dominant process affecting dissolved-phase transport of contaminants in the vadose zone at the SDA. To investigate the effect on predicted aquifer concentrations, two methods of assigning infiltration at the surface of the SDA were explored. The first method employed the spatially variable infiltration rates shown in Fig. 6 that account for topography, soil disturbance, and vegetation. The applied infiltration rates of 0.64, 3.68, and 24.1 cm yr–1 reflected local topography and degree of soil disturbance. The second method employed a simpler approach and applied a uniform rate of 8.5 cm yr–1 that was the weighted average of the spatially varying rates applied inside the SDA. This simpler method neglected the effect of local topography and soil disturbance. Each of these different infiltration assignments were linked with the source release model so that the areas of the domain representing waste saw differing amounts of assigned infiltration between the base simulation and this conceptual uncertainty simulation.
Figure 8
shows the maximum simulated 238U aquifer concentrations as a function of time when the spatially variable and uniform infiltration rates were employed. The peak simulated concentration with the uniform infiltration rate is approximately half an order of magnitude higher than the peak concentration obtained using the spatially varying infiltration rate. The difference in concentration results is caused by more water passing through those areas with waste streams that contain U in the uniform infiltration simulation. The majority of the 238U bearing waste is in Pad A, Pit 4, Pit 6, and Pit 10. In the spatially varying infiltration case, these areas were assigned infiltration rates of 0.6, 2.9, 3.7, and 2.0 cm yr–1, respectively. In the uniform case, each area was assigned a rate of 8.5 cm yr–1. This greater amount of water results in a quicker release of 238U from the source release model which is reflected in the peak aquifer concentration occurring slightly sooner in the uniform infiltration simulation than the spatially varying infiltration simulation. At later times, the aquifer concentration is lower with the uniform infiltration assignment because the contaminant mass has been flushed through the system. For evaluating risk, it is generally the peak concentrations that are most important in the decision making process.
Facilitated Transport
Facilitated transport via sorption onto colloids or organic complexation has long been suspected as a possible mechanism to explain sporadic, very low concentrations of Pu detected in the aquifer near the SDA. These aquifer concentrations are generally just at or slightly above the routine analysis method detection limit of about (0.02 pCi L–1). Some of these low level detections have been confirmed with the same analytic method through reanalysis of a duplicate sample, indicating that Pu may indeed be present (Holdren et al., 2002). However, analysis with more sensitive analytical methods has not confirmed Pu in the aquifer (Roback et al., 2000). Column studies with SDA interbed sediments to evaluate if facilitated transport is a viable mechanism have been completed. The studies have shown that a very small fraction of Pu may become mobile (Fjeld et al., 2000).
To test the effect of conceptual uncertainty associated with facilitated transport, a series of simulations were performed to evaluate if small fractions of Pu released from the buried waste on an annual basis would reasonably explain the aquifer measurements. In the source release model, the annual mobile fractions of plutonium were assumed to be 0 (base case), 1 x 10–6, 1 x 10–4, and 1 x 10–2. For the non-zero mobile fractions, this was equivalent to saying all the mass was released in 1000000, 10000, and 100 yr, respectively. Once this mass was released from the source release model, it was also treated as if it were completely mobile in the vadose zone, which was different from the base case where the partition coefficient for Pu was assigned a site-specific value of 5100 mL g–1 (Dicke, 1997) for the sedimentary features. Figure 9
shows the maximum simulated aquifer 239Pu concentrations as a function of time for the base and range of mobile Pu fraction cases. Simulations were also performed for other actinides, which showed similar results. In each mobile fraction case, the 239Pu concentration rapidly rises, even before the current date that is indicated on Fig. 9. The simulation with the 1 x 10–2 mobile fraction shows a rapid decrease as the source is rapidly depleted and transported through the vadose zone. The decline stops approximately 2000 yr into the simulation, because this was actually a multiple contaminant run that included disposed 243Am decaying into 239Pu. The 243Am decay and ingrowth stabilizes the 239Pu concentration and stops its continued decrease. The 1 x 10–4 and 1 x 10–6 mobile fraction simulations both stabilize without reaching the earlier peak because they have not had all their mass released during the time period shown on Fig. 9. By comparison, the base case with no mobile fraction shows much lower simulated concentrations in the aquifer because the high sorption coefficient retains most of the mass in the source release model and what mass is released is transported very slowly.
Comparisons of the simulated aquifer concentrations to the low-level concentrations detected in the aquifer show that none of the simulated mobile fraction results are plausible. The simulations predict that 239Pu should be widespread and easily detected in the aquifer at the current date, which is not the case. This representation of facilitated transport is simplistic and did not include processes that could filter colloids being transported across the interbeds. Inclusion of attenuating processes would yield slower breakthrough and lower simulated concentrations in the aquifer. However, model calibration would still be unachievable because monitoring for Pu contaminants typically yields nondetect results. Nonetheless, the public interest in Pu transport remains high, even though modeling and monitoring results show transport is not likely.
Variably Saturated Fracture Flow Moisture Characteristic Curve
The majority of the subsurface at the SDA is fractured basalt, which is conceptualized as a high-permeability, low-porosity continuum in the simulations. The constitutive relations describing moisture content, water potential, and hydraulic conductivity for the base case environmental assessment simulations use a Corey-type constitutive relation curve. The parameters for this Corey-type curve were determined empirically through inverse modeling of a wetting front advance from the Large Scale Infiltration Test conducted in 1994 at the INEEL (Magnuson, 1995).
Testing conceptual uncertainty related to the appropriateness of this Corey-type curve for describing water movement in fractured basalt was accomplished by simulating a wetting front moving down inside the SDA and comparing against tensiometer measurements from a network in the SDA. Two other constitutive curves, one using alternate parameters for the Corey curve and one using a van Genuchten curve, were also used to simulate water movement in the fractured basalt to determine which of the three constitutive curves was most representative and which was most conservative with respect to simulating transport. The three fractured basalt constitutive relationship curves are shown in Fig. 10
. The alternate Corey curve was parameterized to reduce the water potential indicated at low saturations. It is not likely to be physically plausible that water potentials of about –300 cm of water could be sustained in fractures. To reduce the water potential at low saturation, the Agw parameter was changed from 325 cm of water to 122 cm of water. The van Genuchten curve parameters for the fractured basalt constitutive relationship curve were = 10.0 m–1, ß = 2.5, and Swr = 0.001. These parameters were selected to result in low water potentials for most water saturations. The permeability curves show large differences between the Corey-type and the van Genuchten constitutive relationships, with the van Genuchten curve being appreciably lower over almost the entire water saturation range.
A series of nested tensiometers in Well 76-5, which is within the SDA and is located along the northern boundary, recorded the movement of three successive wetting fronts (see Fig. 11)
to a depth of at least 11.6 m resulting from snow melting events at the surface. These water potentials were measured in Well 76-5 in both fractured basalt and a thin interbed at a depth of 9.4 m. Figure 12
shows a construction diagram for Well 76-5 with the nested tensiometers. Of the three uppermost tensiometers, two are located within fractured basalt and one is located within the interbed at a depth of 9.4 m. Also shown in Fig. 12 is a one-dimensional simulation grid that was used to represent this site. The remainder of the fractured basalt hydrologic properties and the sediment hydrologic properties assigned to this domain were the same as those in the Base Model Development section above.
Figure 13
shows simulated water potentials with the three 1999 wetting events imposed at land surface for each of the fractured basalt constitutive relationships. Initial conditions for each simulation were obtained by imposing a low background infiltration rate of 1 cm yr–1 and simulating until water potentials were constant throughout the domain. The initial water potential at the 9.4-m depth were greater that those at the 6.7- and 11.6-m depths since the 9.4-m depth was within an interbed and was assigned different hydrologic properties.
View larger version (31K):
[in this window]
[in a new window]
|
Fig. 13. Simulated water potentials in response to surface wetting events with three fractured basalt constitutive relationships.
|
|
The three wetting events were imposed numerically by maintaining saturated conditions at the upper boundary for a duration of 1 d each time. In each simulation result, the series of wetting fronts can be seen passing down through the profile during the late winter and early spring of 1999. In the base case Corey-type curve simulations, the wetting fronts essentially are independent with no cumulative buildup of moisture. The same is true for the alternative Corey simulation. In the latter, the objective of having reduced capillary pressures compared with the base case was achieved. The simulated results using a van Genuchten constitutive relationship best mimic the monitoring data shown previously in Fig. 11, with the combined effect of multiple wetting fronts showing progressively wetter conditions with time and the lowest simulated capillary pressures. Of the three curves chosen for this simulation, the van Genuchten constitutive curve shows the most promise for accurate representation of water movement in the fractured basalt and would be the best candidate for further calibration efforts to improve the agreement between simulated and observed water potentials. Since the objective of this effort was to determine which constitutive relationship was most representative and which was most conservative, further calibration was not attempted.
The simulated water fluxes out the bottom of the simulation domain for all three curves are shown in Fig. 14
. The Corey-type curves show a faster and more complete advance of the wetting front through the simulation domain. From a conceptual uncertainty perspective, this demonstrates that a conservative approach was used to simulate water movement in the fractured basalt portions of the vadose zone for the base case environmental assessment simulations. As contaminant monitoring data become available against which to calibrate simulations, the constitutive curves used for the fractured basalt may need to be modified to obtain simulation results that accurately represent flow and transport in the vadose zone.
View larger version (16K):
[in this window]
[in a new window]
|
Fig. 14. Simulated flux out bottom of simulation domain with three fractured basalt constitutive relationships.
|
|
|
CONCLUSIONS
|
|---|
A model was developed to simulate flow and transport of contaminants as part of a CERCLA site assessment for the SDA on the INEEL. A limited set of conceptual uncertainties for the model were assessed by implementing conceptual variations in different infiltration rates at the surface, including facilitated transport, and using alternative methods to describe hydraulic properties of variably saturated fractured basalts comprising the majority of the subsurface. The water and contaminant fluxes out the bottom of the vadose zone model were input into an aquifer model to estimate aquifer concentrations. The sensitivity of the results to these conceptual variations has been used to guide additional characterization activities at the SDA. Of the sensitivities discussed in this paper, the primary sensitivity was to the possibility of facilitated transport occurring, although none of the scenarios in which facilitated transport was included were physically plausible. Lesser sensitivities were shown to the amount and method of assigning infiltration at land surface. The hypothesis that conceptual uncertainty dominates parametric uncertainty was utilized as the basis for focusing on conceptual uncertainty. Given the current state of the model development for representing flow and transport beneath the SDA, this hypothesis is still believed to be valid.
Improvements in modeling for the SDA will result primarily from calibration to monitoring data within the SDA waste pits and in the deeper vadose zone. Until representativeness of the model is established through calibration, conceptual uncertainty will continue to dominate the simulation results. Once the model is adequately calibrated to observed dissolved-phase transport through the vadose zone, quantification of the parametric uncertainty and its impacts on predicted concentration can become the focus of characterization efforts.
|
ACKNOWLEDGMENTS
|
|---|
This work was supported by the U.S. Department of Energy, Assistant Secretary for Environmental Management, under DOE Idaho Operations Office Contract DE-AC07-99ID13727. The author thanks the two anonymous reviewers for the insightful comments on this report.
|
REFERENCES
|
|---|
- Anderson, S.R., and B.D. Lewis. 1989. Stratigraphy of the unsaturated zone at the Radioactive Waste Management Complex, Idaho National Engineering Laboratory, Idaho. Water-Resources Investigations Rep. 89-4065, DOE/ID-22080. USGS, Idaho Falls, ID.
- Baca, R.G., S.O. Magnuson, H.D. Nguyen, and P. Martian. 1992. A modeling study of water flow in the vadose zone beneath the Radioactive Waste Management Complex. EGG-GEO-10068. EG&G Idaho, Inc., Idaho Falls, ID.
- Bear, J. 1972. Dynamics of fluids in porous media. American Elsevier Publishing Co., New York.
- Bhorgese, J.V. 1988. Hydraulic characteristics of soil cover, Subsurface Disposal Area, Idaho National Engineering Laboratory. M.S. thesis. University of Idaho, Moscow.
- Carrera, J. 1993. An overview of uncertainties in modeling groundwater solute transport. In J.I. Kim and G. de Marsily (ed.) Chemistry and migration of actinides and fission products. J. Contam. Hydrol. 13:23–48.
- Cecil, L.D., J.R. Pittman, T.M. Beasley, R.L. Michel, P.W. Kubik, P. Sharma, U. Fehn, and H. Gove. 1992. Water infiltration rates in the unsaturated zone at the Idaho National Engineering Laboratory estimated from chlorine-36 and tritium profiles, and neutron logging. In Y.K. Kharaka and A.S. Meest (ed.) Proceedings of the 7th International Symposium on Water–Rock Interaction—WRI-7, Park City, UT. Balkema, Rotterdam, The Netherlands.
- Dicke, C.A. 1997. Distribution coefficients and contaminant solubilities for the Waste Area Group 7 baseline risk assessment. INEL/EXT-97–00201. Lockheed Martin Idaho Technologies, Idaho Falls, ID.
- Dunnivant, F.M., M.E. Newman, C.W. Bishop, D. Burgess, J.R. Giles, B.D. Higgs, J.M. Hubbell, E. Neher, G.T. Norrell, M.C. Pfiefer, I. Porro, R.C. Starr, and A.H. Wylie. 1998. Water and radioactive tracer flow in a heterogeneous field-scale system. Ground Water 6:949–958.
- Fjeld, R.A., J.T. Coates, and A.W. Elzerman. 2000. Final report, Column tests to study the transport of plutonium and other radionuclides in sedimentary interbed at INEEL. INEEL/EXT-01–00763. Clemson University, Department of Environmental Engineering and Science, Clemson, SC, for the Idaho National Engineering and Environmental Laboratory, Bechtel BWXT Idaho, LLC, Idaho Falls, ID.
- Holdren, K.J., B.H. Becker, N.L. Hampton, L.D. Koeppen, S.O. Magnuson, T.J. Meyer, G.L. Olson, and A.J. Sondrup. 2002. Ancillary Basis for risk analysis of the Subsurface Disposal Area. INEEL/EXT-02–01125. Bechtel BWXT Idaho, LLC, Idaho Falls, ID.
- Leecaster, M. K. 2002. Geostatistic modeling of subsurface characteristics in the Radioactive Waste Management Complex region. Operable Unit 7–13/14. INEEL/EXT-02–00029, Revision 0. Bechtel BWXT Idaho, LLC, Idaho Falls, ID.
- Magnuson, S.O. 1995. Inverse modeling for field-scale hydrologic and transport parameters of fractured basalt. INEL-95/0637. Lockheed Martin Idaho Technologies Company, Idaho Falls, ID.
- Magnuson, S.O., and D.L. McElroy. 1993. Estimation of infiltration from in situ moisture contents and representative moisture characteristic curves for the 30', 110', and 240' interbeds. EG&G Engineering Design File RWM-93-001.1. EG&G Idaho, Inc., Idaho Falls, ID.
- Martian, P. 1995. UNSAT-H infiltration model calibration at the Subsurface Disposal Area, Idaho National Engineering Laboratory. INEL-95/0596. Lockheed Idaho Technologies Company, Idaho Falls, ID.
- McElroy, D.L., and J.M. Hubbell. 1990. Hydrologic and physical properties of sediments at the Radioactive Waste Management Complex. EGG-BG-9147. EG&G Idaho, Inc., Idaho Falls, ID.
- Mualem, Y. 1976. A new model for predicting the hydraulic conductivity of unsaturated porous media. Water Resour. Res. 12:513–522.
- Neuman, S.P. 2002. Accounting for conceptual model uncertainty via maximum likelihood bayesian model averaging. In Proc. 4th Int. Conf. on Calibration and Reliability in Groundwater Modeling: A few steps closer to reality, Prague, Czech Republic. 17–20 June 2002. Geologica, Berkeley, CA.
- Nimmo, J.R., K.S. Perkins, P.A. Rose, J.P. Rousseau, B.R. Orr, B.V. Twining, and S.R. Anderson. 2002. Kilometer-scale rapid transport of naphthalene sulfonate tracer in the unsaturated zone at the Idaho National Engineering and Environmental Laboratory. Available at www.vadosezonejournal.org. Vadose Zone J. 1:89–101.[Abstract/Free Full Text]
- Olson, G.L., D.L. Koeppen, A.M. Parsons, P.D. Ritter, and A.J. Sondrup. 2003. FY 2002 Environmental monitoring report for the Radioactive Waste Management Complex. INEEL/EXT-03-00055, Rev. 0. Bechtel BWXT Idaho, LLC, Idaho Falls, ID.
- Parker, J.C., R.J. Lenhard, and T. Kuppusamy. 1987. A parametric model for constitutive properties governing multiphase flow in porous media. Water Resour. Res. 23:618–624.
- Roback, R.C., D.W. Efurd, M.T. Murrell, R.E. Steiner, and C.J. Duffy. 2000. Assessment of uranium and plutonium in the saturated and unsaturated zones beneath the Surface Disposal Area, INEEL. INEEL/EXT-01-00771, Rev. 0. Los Alamos National Laboratory, Los Alamos, NM.
- Rodriguez, R.R., A.L. Schafer, J. McCarthy, P. Martian, D.E. Burns, D.E. Raunig, N.A. Burch, and R.L. VanHorn. 1997. Comprehensive RI/FS for the Idaho Chemical Processing Plant OU 3-13 at the INEEL—Part A. RI/BRA Rep. DOE/ID-10534. Lockheed Martin Idaho Technologies Company, Idaho Falls, ID.
- Shook, G.M. 1995. Development of an environmental simulator from existing petroleum technology. INEL-94/0283. Lockheed Martin Idaho Technologies Company, Idaho Falls, ID.
- van Genuchten, M.Th. 1980. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci. Soc. Am. J. 44:892–898.[Abstract/Free Full Text]
- Vinsome, P.K.W., and G.M. Shook. 1993. Multi-purpose simulation. J. Pet. Sci. Eng. 9:29–38.
This article has been cited by other articles:
|
|
|
|
 
C. L. Duke, R. C. Roback, P. W. Reimus, R. S. Bowman, T. L. McLing, K. E. Baker, and L. C. Hull
Elucidation of Flow and Transport Processes in a Variably Saturated System of Interlayered Sediment and Fractured Rock Using Tracer Tests
Vadose Zone J.,
November 20, 2007;
6(4):
855 - 867.
[Abstract]
[Full Text]
[PDF]
|
|
|
TOP
ABSTRACT
INTRODUCTION
