VZJ sign up for etocs
HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
 QUICK SEARCH:   [advanced]


     


This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in Web of Science
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via Web of Science (5)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
Related Collections
Right arrow Volatile Organic Compounds
Right arrow Multicomponent Transport Models
Right arrow Organic Compounds
Published in Vadose Zone Journal 3:1262-1275 (2004)
© 2004 Soil Science Society of America
677 S. Segoe Rd., Madison, WI 53711 USA

ORIGINAL RESEARCH

Modeling Natural Attenuation of Multicomponent Fuel Mixtures in the Vadose Zone

Use of Field Data and Evaluation of Biodegradation Effects

Petros Gaganisa, Peter Kjeldsenb and Vasilis N. Burganosa,*

a Institute of Chemical Engineering and High Temperature Chemical Processes–Foundation for Research and Technology, Hellas (ICE/HT-FORTH), P.O. Box 1414, GR 26504 Patras, Greece
b Environment & Resources DTU, Technical University of Denmark, DK-2800 Kgs. Lyngby, Denmark

* Corresponding author (vbur{at}iceht.forth.gr)

Received 21 October 2003.



    ABSTRACT
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Data from a controlled fuel source emplacement field experiment were used to explore the potential for predictive modeling using an approximate model and parameters from laboratory-scale studies, and to quantify the effects of biodegradation on the natural attenuation of volatile organic compound (VOC) mixtures in the vadose zone. The application of a species grouping technique is investigated in an attempt to reduce the computational cost for predicting the fate of multicomponent mixtures with a minimum loss in accuracy. We show how several difficulties in vadose zone transport modeling may be reasonably well overcome by simplifications that are supported by field data and sensitivity calculations. The estimated case-specific field-scale biodegradation values for the majority of mixture compounds were found to fall within the range of biodegradation rate constants determined from experiments performed in columns that were filled with the field soil. A sensitivity analysis showed that using the upper or lower bounds of literature values for the biodegradation constants had a negligible effect on gas concentrations at the early stages of the contamination process, when contaminant concentrations were still high. This is a very encouraging result regarding the usefulness of model calculations in risk assessment procedures.

Abbreviations: NAPL, nonaqueous phase liquid • VOC, volatile organic compound


    INTRODUCTION
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
WHEN USING mathematical models to simulate the fate of multicomponent VOC mixtures in the vadose zone it is necessary to adopt certain simplifications since it is practically impossible to either characterize or mathematically describe the true complexity of the physical system. Specifically, the need for simplifications arises from limited knowledge of the processes themselves or their synergy, limitations of modeling tools, a dearth of data, and computational cost constraints (Gaganis and Smith, 2001). In spite of the model limitations and approximations, there is a strong need for predictive numerical tools. Even in simplified versions, models can be especially useful in the early stages of a contamination event when some first estimates for a preliminary risk assessment are needed.

The fate of VOC mixtures in the vadose zone is typically controlled by a combination of transport phenomena (i.e., advection, dispersion, diffusion) and phase change and reaction processes (i.e., volatilization, sorption, and biodegradation). Vadose zone modeling becomes particularly complicated when the number of organic compounds in the pollution source increases, and interactions among phases, interfaces, and chemical components are too significant to neglect. Other difficulties encountered in VOC mixture modeling include uncertainty in the evaluation of the species properties and their usually unknown temperature dependence. For example, the evaluation of sorption distribution coefficients is, in most cases, only approximate, and typically based on octanol–water partition coefficients (Kow) and measurements of the organic C of the soil. As revealed from field (Hers et al., 2000) and laboratory (Franzmann et al., 1999; Höhener et al., 2003) studies, biodegradation is one of the most important attenuation mechanisms of VOCs in the vadose zone. However, the biodegradation process is site specific. It is common for the biodegradation rate constants reported in the literature to vary by orders of magnitude for the same compound, thus making it difficult to choose a rate constant value for a specific modeling application. In the field, it is usually difficult to distinguish and quantify biodegradation among the different natural attenuation processes such as volatilization and sorption without the aid of numerical modeling. Laboratory experiments, on the other hand, are usually conducted under conditions different from those at the site under consideration. Consequently, direct adoption of such biodegradation values in field-scale models may be questionable.

A controlled fuel source emplacement field experiment was conducted recently at Airbase Værløse, Denmark. Its purpose was to investigate the potential of natural attenuation of a hydrocarbon mixture in the vadose zone, and to provide a better evaluation of the risk for groundwater contamination. Using selected data from the above field experiment, the objectives of the present study were to (i) identify the important factors that affect the natural attenuation and transport of VOC mixtures in the vadose zone, (ii) evaluate the effects of the biodegradation process and quantify the biodegradation rates at the site, (iii) validate and assess the applicability of a compound grouping methodology (Gaganis et al., 2002) at the field scale for reducing computational requirements in the case of a large number of constituents, and (iv) explore the potential for predictive modeling using an approximate model and parameters from the literature and laboratory studies. An important aspect of the Værløse field experiment is the known amount, composition, geometry, and time of source emplacement. Consequently, the results we present here should be particularly important given that these parameters are among the most important yet most uncertain factors in such studies. These factors are known to strongly influence the reliability of the results of most investigations of actual spill sites.

Following a brief description of the field experiment and the mathematical formulation, we present a parametric analysis for the purpose of identifying the key processes and parameters that are expected to play a significant role in numerical simulations of the attenuation process. Such a preprocessing stage is very significant for problems involving contamination of the vadose zone by multicomponent VOC mixtures inasmuch as it can substantially reduce the computational cost without sacrificing accuracy. This is particularly important when phenomenological parameters used in the model have both high uncertainty in their value and significant effect on the calculations. The results of this analysis serve as a guide to the field modeling, which employs the constituent averaging technique. Subsequently, a sensitivity analysis is used to identify the effects of potential uncertainties in key parameter estimates (e.g., biodegradation constants) and to assess the predictive capability of the numerical approach. Finally, the main conclusions are summarized in the last section, along with further remarks on the modeling of vadose zone contamination by VOC mixtures.


    MATERIALS AND METHODS
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
The hydrologic cross section of the Værløse site includes a sandy dark brown top soil of approximately 50-cm thickness overlying a 2- to 3.3-m layer of homogeneous glacial melt water sand, which, in turn, overlies a thin layer (0.5–1 m) of moraine sand–gravel. The moraine sand–gravel is underlain by a thick moraine clay layer (up to 50 m). The water table ranges from 2.5 to 3.5 m below ground surface. The saturated hydraulic conductivity as measured on core samples was 3.1 m d–1 for the depth interval of 30 to 40 cm below the soil surface, and 8.3 m d–1 for depths of 80 to 130 cm. The field experiment was started in the beginning of July 2001. A cylindrical source (diameter of 0.75 m) consisting of sand from the site mixed with the oil phase was placed between 0.8 and 1.3 m below ground surface. The source was covered to prevent direct leaching of the oil phase with infiltrating rainwater. The oil phase was an artificial hydrocarbon mixture of 14 volatile and semivolatile compounds (BTXs, n-, iso- and cycloalkanes) similar to jet fuel (Table 1). Samples from the source material were taken immediately after the installation and four times (Day 23, 113, 225, and 352) during the field experiment. The volume of oil contained in the source (10.2 l) at the beginning of the experiment was measured immediately after the source installation. Monitoring devices were placed mainly in one radial transect (direction of groundwater flow) with only control monitoring points in the other directions.


View this table:
[in this window]
[in a new window]
 
Table 1. Fuel mixture composition after source installation and physicochemical properties of individual compounds.

 
The monitoring network of the site and the source consisted of 107 soil gas probes (Fig. 1) , seven multilevel water samplers (four of them located on the main radial transect 0, 1, 5, and 10 m from the source), and six porous cups for sampling water from the unsaturated zone 1 and 2 m from the source. Within the main radial transect, the hydrocarbon composition of the pore gas was monitored to a distance of 20 m from the center of the source. The migration of the hydrocarbons and tracers in soil gas, pore water, and groundwater was monitored until July 2002, when the source was removed and the site was remediated. Experimental concentration contours for hexane and CFC-113 along the main radial transect at the Værløse site are presented in Fig. 1. The soil moisture profile was also continuously monitored in one location during the experiment using time domain reflectometry probes, while various other parameters were estimated by laboratory tests (e.g., saturated hydraulic conductivity, parameters {alpha} and n for the van Genuchten model, expected biodegradation rates). At all times, the 50-cm topsoil showed higher water content than the underlying sand. The organic C content in the upper 50-cm zone was 0.1 to 1.5% and was 0.03% in the glacial melt water sand. For a more detailed description and experimental results of the field experiment, the reader is referred to Christophersen et al. (2002), Broholm et al. (2003), and Kjeldsen et al. (2003).



View larger version (31K):
[in this window]
[in a new window]
 
Fig. 1. Location of soil gas probes and experimental concentration contours along the main radial transect at the Værløse site.

 

    MODELING
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Mathematical Formulation
The governing equations for three-phase flow in the subsurface, assuming an incompressible porous medium, incompressible liquid phases, and a compressible gas phase, may be written as

[1]

[2]

[3]

[4]
where the subscript p denotes a particular phase (organic, aqueous, or gas phase), {phi} is the porosity, {rho} is the fluid density, S is the fluid saturation, Kpij is the p-phase conductivity tensor, krp is the relative permeability of phase p, Kswij is the saturated conductivity tensor for water, nrp is the ratio of the p phase to the water viscosity, hp is the water height-equivalent pressure head of phase p, P is the p-phase pressure, g is the acceleration of gravity, {rho}w* is the density of pure water, Qp is the source (+) or sink (–) term, and t is time. When gas phase flow is considered, gas compressibility also enters the calculations through the relation (Katyal et al., 1991):

[5]
where hg is the water height-equivalent pressure head in the gas phase, {lambda} is the gas compressibility, assigned here a value of 1.17 x l0–6 g cm–4, {rho}og is the density of native soil air, and {rho}cw is the density of contaminants in the gas phase, which is determined from the solution of the transport model. Assuming equilibrium partitioning among phases, the transport of VOC mixture component {alpha} in the p phase inside a porous medium is given by

[6]

[7]

[8]
where C{alpha}p (M L–3) is the concentration of the {alpha} component in p phase (o, oil; w, water; g, gas) expressed as mass per phase volume, D{alpha}pij (L2 T–1) is the dispersion tensor, Ddif{alpha}p (L2 T–1) is the diffusion coefficient of {alpha} in the p phase of the porous medium (estimated from bulk phase molecular diffusion coefficient using the tortuosity model of Millington and Quirk, 1959), Ddisppij (L2 T–1) is a mechanical dispersion coefficient, qpi (L T–1) is the p-phase Darcy velocity, µ{alpha}p (T–1) is an apparent first-order decay coefficient, Qap is a source (+) or sink (–) of mass of constituent {alpha} in the p phase due to interphase mass exchange (i.e., dissolution, volatilization) (M L–3 T–1, K{alpha}o is the equilibrium partition coefficient for the {alpha} component between water and organic liquid (Raoult's constant), and K{alpha}g is the equilibrium partition coefficient for the {alpha} component between water and gas (Henry's constant). Linear, instantaneous sorption of {alpha} from the aqueous phase onto the solid phase is accommodated as

[9]
where Kd (L3 m–1) is the linear Freundlich sorption coefficient, {rho}b (M L–3) is the porous medium bulk density, and C{alpha}s is expressed as mass of adsorbed component {alpha} per porous medium volume. The flow and transport equations were solved using the numerical solver of the MOFAT code (Katyal et al., 1991) to take also advantage of several powerful features, such as, the three-phase extension of the van Genuchten model (van Genuchten, 1980) to describe the relationships between phase permeabilities, saturations, and pressures at locations where all three phases appear simultaneously (Kaluarachchi and Parker, 1989, 1992).

Because temperature variations have a profound effect on volatilization, dissolution, and the biodegradation of VOC compounds, an external module was developed and coupled to MOFAT. The module updated those fuel compound properties that are most sensitive to temperature at selected times according to literature expressions (Table 2). Numerical simulations were implemented through a batch file:

  1. The simulation was terminated at time t1.
  2. Property values that are functions of temperature were updated using the new external module.
  3. The simulation continued using the model output at time t1 as initial condition.


View this table:
[in this window]
[in a new window]
 
Table 2. Temperature dependence of vapor pressure and Henry's Law constant of fuel mixture compounds.{dagger}

 
The frequency of repeating this procedure was dictated by the rate of variation of the temperature measured at the field with time.

Identification of Case-Specific Significant Processes
As mentioned above, attenuation of VOCs in the vadose zone involves a multitude of transport phenomena and reaction processes, each of which requires knowledge of parameter values that are typically site-specific and, in many cases, are characterized by moderate to high uncertainty. For this reason, it is important to identify first a limited number of processes that most significantly control the plume spreading. Doing this should not only simplify the complexity of the model structure, but also save considerable time and effort in estimating phenomenological coefficients that may have only a negligible effect on the dynamics of the process. Still, care must be exercised to ensure proper account for processes that may vary in significance at different stages of the plume spreading, as will be shown below.

In the case of volatile compounds, much evidence exists in the literature indicating that certain mechanisms are important to VOC attenuation in the vadose zone. For instance, vapor diffusion cannot be ignored, considering that gas phase diffusion coefficients are four to five orders of magnitude higher than diffusion coefficients in the liquid phase (Sleep and Sykes, 1989; Grathwohl, 1998). Field investigations have identified volatilization and biodegradation as the key processes that control VOC attenuation in the unsaturated zone (e.g., Lahvis and Baehr, 1996; Lahvis et al., 1999). Similar results were observed in laboratory experiments (e.g., El-Farhan et al., 1998). Biodegradation appears to be a key factor that may greatly reduce the spreading distance, the maximum gaseous concentrations, and the time to depletion for all biodegradable mixture components, especially for those with a high aqueous solubility and a low Henry's constant (e.g., benzene). A numerical analysis by Thomson et al. (1997) showed that water table fluctuations and seasonal variations in moisture content and temperature can have a significant impact on groundwater and gas concentrations. Specifically, the dissolution and volatilization rates can vary by almost 100% of their average values, whereas gas concentrations can vary by 50% due to the effects of temperature alone.

The main features (e.g., depth to the water table, location of contaminant source, boundary conditions) of the hydrogeological setting adopted in our numerical parametric analysis are identical to those of the field site. Two-dimensional Cartesian coordinates were adopted to incorporate the effects of regional groundwater flow. A number of simulations were performed for a contaminant mixture consisting of five compounds with equal mass fractions (i.e., benzene, methylcyclopentane, hexane, 1,2,4-trimethylbenzene, and decane), which were selected from those used in the field experiment. These compounds were selected to cover a broad range of property values of those in the actual nonaqueous phase liquid (NAPL) mixture. The parameter values used for the simulated compounds are given in Table 1. Biodegradation is described by means of first order-kinetics. The biodegradation rate coefficients used in our analysis were 0.20, 0.68, 0.46, 1.0, and 97.0 d–1 for benzene, methylcyclopentane, hexane, 1,2,4-trimethylbenzene, and decane, respectively. These values are within the range of biodegradation values reported in the literature for the above compounds. Equilibrium partitioning among phases was assumed in these simulations, which is justified by the fact that the evolution in source composition at the Værløse site appears to follow closely that calculated using Raoult's Law (Broholm et al., 2003).

In the present analysis we investigated the effects of the following additional factors, whose exclusion may allow significant simplifications in the hydrogeologic model: (i) groundwater flow, (ii) temporary variable barometric pressure, and (iii) advective transport in the gas phase. Numerical simulations indicated that the majority of the contaminants escape to the atmosphere following direct volatilization from the source (Fig. 2) . However, this process becomes less important at higher water content. As can be seen in Fig. 2, the contaminant mass dissolved in the aqueous phase (in both the vadose and saturated zones) represents only about 1% of the initial contaminant mass in the source, while the mass in the aqueous phase below the water table is at least one order of magnitude lower. Mass transfer to flowing groundwater, hence, does not have a significant effect on VOC transport in this case. This result is also consistent with the field observation (Kjeldsen et al., 2003) that few compounds were found at the water table, and then only for a short period of time.



View larger version (21K):
[in this window]
[in a new window]
 
Fig. 2. Temporal evolution of mass in the oil, aqueous, and gas phases, and contaminant mass escaping to the atmosphere for the base case of the parametric analysis.

 
It is generally accepted that gas-phase diffusion from the source is a key transport mechanism in the vadose zone. However, in our case, considering diffusion as the sole VOC transport process in the gas phase (i.e., ignoring gas advection) greatly underestimated travel distances, at times during 0 to 10 d, for the majority of the volatile mixture compounds (up to 50% at Day 6). The enhanced gas advection at early times appears to result from overpressure created in the source area as a consequence of evaporation of the volatile compounds. This explanation is consistent with the calculated overpressure of approximately 0.006 MPa (0.06 atm) in the source area at the Værløse site immediately after source installation. At later times, as the source is being depleted, gas phase advective flow becomes progressively less significant and the VOC transport process is dominated by diffusion. However, simulation of early time air-phase advection appears essential in our study for proper interpretation of our field data and for the inverse analysis aimed at estimating reliable biodegradation rate values. In addition to increasing the travel distance, gas advection due to overpressure in the source area also favors the escape of contaminants to the atmosphere at early times (see Fig. 2). This may considerably influence the transport process at later times and result in overestimation of the actual biodegradation rates.

Temporary variable barometric pressure (as measured at the field site, with mean value of 0.101 MPa [1 atm], and a standard deviation of 0.0008005 MPa [0.0079 atm]) at the ground surface in the model had only a minor effect on the gas concentrations. Figure 3 shows a comparison between breakthrough curves obtained for constant and fluctuating barometric pressure conditions, in terms of normalized gas concentrations (i.e., divided by the product of aqueous solubility and Henry's constant) at a horizontal distance 2.5 m from the source. The fluctuations in the plotted gas concentrations, are about 0 to 5% around the breakthrough curve obtained for a constant (0.101 MPa, 1 atm) barometric pressure.



View larger version (22K):
[in this window]
[in a new window]
 
Fig. 3. Temporal evolution of gas concentrations of two mixture components at a distance of 2.5 m from the source for both constant and fluctuating barometric pressure conditions.

 
The parametric analysis in this section indicates that reliable modeling of the Værløse field site necessitates a numerical model that accounts for several other processes in addition to those expected to be important for VOC transport, such as biodegradation, sorption, and partitioning among phases. The model must also accurately describe the soil types and water content profiles at the site, and account for advective transport in the gas phase as well as for the effects of temperature variation on VOC concentrations. Simulation of the temporal evolution of VOC mixture composition is also an important factor that significantly affects the extent of the gas-phase plumes. On the other hand, barometric pressure fluctuations and groundwater flow in the saturated zone were found to have an insignificant effect on VOC transport in this case. The latter result is particularly significant since it allows us to employ a two-dimensional model with cylindrical coordinates, at substantial reduction in computational cost without compromising the accuracy of our predictions.

Hydrogeological Model
The main features of the hydrogeologic simulation model include (i) quantified relationships between phase permeabilities, saturations, and pressures according to Kaluarachchi and Parker (1989)(1992); (ii) flow due to pressure and density differences and advective transport in the gas phase; (iii) equilibrium partitioning between the gas–aqueous, gas–NAPL, and aqueous–solid phases; and (iv) modeling of biodegradation, which takes place only in the aqueous phase, and is assumed to be a first-order process. As shown in the section above, the second feature may be especially important for VOC transport in the vadose zone. Furthermore, even if gas pressure gradients have a negligible effect on liquid phase flow, significant gas flow can occur either due to changes in fluid saturations with time, or due to pressure and density gradients in the gas phase resulting from phase partitioning. Regarding the biodegradation process, all relevant laboratory data indicated first-order kinetics for the particular organic components and soil found at the Værløse site (Höhener et al., 2003). Furthermore, field data showed that the pore air phase did not contain more than 2.5% CO2 and <17% O2 at all depths and all times, thus indicating that the supply of electron acceptors was not a limiting factor at the site. On the basis of this information and available laboratory data on biodegradation, a first-order biodegradation model was adopted in our study to simulate VOC transport in the vadose zone at the Værløse site. As mentioned, the assumption of equilibrium partitioning between phases appeared to produce good estimates of the actual field measurements of the source composition at the Værløse site (Broholm et al., 2003) as functions of time.

The field site exhibited relatively large temperature changes. The temperature at the source had been fairly constant (maximum fluctuation of about 1°C) at about 17°C for 60 d after the beginning of the experiment. Following this, however, the temperature steadily decreased to become slightly <5°C at Day 167. Among the parameters that most strongly affect VOC transport, the vapor pressure, Henry's Law constant, and aqueous solubility were found to be the most sensitive to temperature variations. Consequently, these three parameters were treated as temperature dependent in the model, whereas all other parameters (e.g., molecular diffusion coefficient, and gas and water viscosities) were assumed to remain constant throughout the experiment because of their lesser impact (Thomson et al., 1997). For the particular case considered here, the values of the temperature-dependent properties were updated using the new module that was coupled with the code, at times 75 and 130 d when the temperature change was >3°C. Table 2 shows the relationships that were used to express the temperature dependence of the vapor pressure and Henry's Law constant for the compounds under consideration. The aqueous solubility values at different temperatures were calculated as the ratio of the vapor pressure/Henry's Law constant. There are no quantitative data in the literature regarding the temperature dependence of biodegradation. However, the expected variations are probably small compared with the uncertainty in the values of the rate constant. Hence, biodegradation was treated in the model as being independent of temperature.

Cylindrical coordinates in a 20 by 5 m domain were used for the simulations (Fig. 4) . The model domain contained 24 grid-blocks in the horizontal direction and 16 in the vertical direction, and was divided into five homogeneous zones with different soil characteristics (Fig. 4). The available data could not justify a more detailed description of the heterogeneity at the site. The model domain was discretized in both directions using a variable nodal spacing of 15 to 20 cm near the source zone (where the largest concentration gradient occurred) and the water table, and gradually larger nodal spacing toward the outer boundaries (up to 2 m in the horizontal direction). This discretization scheme was found to provide a balance between accuracy in the calculations and savings in computational time. As indicated by field measurements in this particular study (Kjeldsen et al., 2003), and by results of the parametric analysis, ignoring the influence of groundwater flow did not have a significant impact on the transport process. Thus, groundwater flow was not incorporated in the model. The boundary conditions for fluid flow and contaminant transport were (i) constant atmospheric pressure, time dependent groundwater recharge, and constant concentration at the upper boundary (C = 0); (ii) temporally varying water table elevation according to the field measurements, which implies a time-varying specified hydraulic head below the water table, no-flow conditions above the water table, and zero-concentration gradient along the left and right boundaries; and (iii) no-flow and zero-concentration gradient at the lower boundary (Fig. 4). The soil properties assigned to the five zones with different soil types (Table 3) were obtained from independent laboratory measurements on soil samples, or estimated from the profiles of the grain-size distribution with depth at the site.



View larger version (22K):
[in this window]
[in a new window]
 
Fig. 4. Model domain and boundary conditions.

 

View this table:
[in this window]
[in a new window]
 
Table 3. Property values assigned to the five soil types used in the model.

 
Grouping Approach
Although the numerical model used in this work was able to simulate all site-specific important processes identified in the parametric study, it could only handle a maximum of five compounds. This is only a fraction of the many organic species making up kerosene (14) used in the field experiment. Similar restrictions hold for the majority of published approaches and codes for unsaturated zone modeling of VOCs, which account for either a single compound or only a limited number of compounds to cope with the complexities of VOC fate and transport (see Karapanagioti et al., 2003, for a relevant literature report). Because of this, we grouped the mixture compounds into effective pseudospecies, that is, groups of species sharing similar property values. Given the synergistic impact of the various processes on the temporal evolution of the concentration profiles, such a grouping must be implemented in a way that ensures minimum error in the calculation of individual species transport. A commonly used approach to treat multicomponent mixture transport is to consider the attenuation of a small number of "effective" species that are assigned property values equal to the mole fraction–weighted averages of the property values of all member compounds, usually of similar chemistry (e.g., Baehr and Corapcioglu, 1987; Klenk, 2000). The accuracy of this approach was examined by Karapanagioti et al. (2004) through detailed calculations for a multicomponent mixture. They concluded that such a simplification would not only result in different shapes of the concentration profiles but also overestimate vapor phase concentrations of the mixture components by one to two orders of magnitude.

An alternative grouping approach, termed the constituent averaging technique, was recently suggested by Gaganis et al. (2002) and is used here also. Gaganis et al. (2002) suggested an objective function for determining case-specific criteria for grouping individual compounds into "effective" species (composite constituents) and formulated a new algorithm for estimating their effective thermodynamic properties. A brief description of the approach is given below for completeness.

A composite constituent is defined to have the following concentration at all times:

[10]
where Cc{alpha} (M L–3) is the concentration of the composite constituent in phase {alpha} (o, organic phase; w, aqueous phase; g, gas phase) and Ci{alpha} is the concentration of member component i = 1, 2,..., l in phase {alpha}. As opposed to individual compounds, the thermodynamic properties of composite constituents can be time dependent due to temporal composition changes resulting from different mass transfer rates of the member components from the pure phase to the aqueous and gas phases. Gaganis et al. (2002) introduced an objective function for minimizing the error associated with the time dependence of the effective property values of the composite constituents through a selection of an appropriate grouping criterion, which is case specific since it is influenced by the hydrogeologic setting, the nature of the contaminants, the processes involved, and the mixture composition. The procedure for selecting the optimum grouping criterion includes (i) random selection of a small number n of different sets of m individual hydrocarbons (m = maximum number of constituents that can be cost-effectively handled by the numerical code) from the actual organic mixture, (ii) simulation of the temporal evolution of the mole fraction in the organic phase of individual compounds in each set, (iii) investigation of the temporal change of the composition of the composite constituents (organic phase) for alternative grouping schemes, and (iv) selection of the property that minimizes the following objective function:

[11]
where p is the property used to form the composite constituents; j = 1,..., n represents a set of m individual hydrocarbons; t is time; i = 1,..., k represents a composite constituent that consists of two individual components ai and bi; and Ri0 and Rit are the mole fraction ratios of ai and bi (Ri = xai/xbi) at times zero and t, respectively.

Because source composition changes are generally more rapid in the early stage of a contamination event, the simulations in this stage may need to be performed only for a fraction of total time. The property that minimizes the objective function J(p) is selected as the optimum grouping criterion since it is associated with the smallest composition variability with time. Furthermore, Gaganis et al. (2002) suggested an algorithm to relate the effective thermodynamic properties of the composite constituents to the properties of the individual members. Based on data from the lysimeter experiment presented in Pasteris et al. (2002), they showed that the effective aqueous solubility and Henry's Law constant of a composite constituent are best approximated by the mole fraction–weighted geometric and arithmetic averages of the solubility and Henry's Law constants of the individual components, respectively. They also demonstrated that the values of the sorption distribution coefficients and the biodegradation rate coefficients of the least sorbing and least biodegradable member hydrocarbons (smallest reported or measured sorption and biodegradation rate coefficients) should be assigned to the composite constituents. An evaluation of the accuracy of the grouping approach with respect to different grouping criteria, and substantiated by a comparison with simulations that include all mixture components individually, can be found in Karapanagioti et al. (2004).

Numerical Simulations and Parameter Estimation
The flow model was calibrated to the measured soil moisture profile by adjusting the recharge rate to the actual precipitation data. A comparison of simulated and measured soil moisture profiles at the beginning of the experiment (time = 0 d) and at Day 97 is shown in Fig. 5 .



View larger version (13K):
[in this window]
[in a new window]
 
Fig. 5. Experimental and simulated soil moisture profiles at the beginning of the experiment (time = 0 d) and at Day 97.

 
For the transport problem, one compound was selected from the fuel mixture as a target contaminant, and simulated as an individual constituent, whereas the other hydrocarbons in the mixture were simulated using a small number of composite constituents following the constituent averaging technique described above. In this case, the composite constituents were used to determine the effect of the mixture composition on the transport process of the selected target compound. Independent simulations for seven different sets of five individual hydrocarbons, randomly selected from the artificial kerosene components, were performed to select the optimum criterion for grouping the fuel mixture compounds into composite constituents. The initial mass fractions were scaled as if the entire mixture consisted of the selected five compounds only (Table 1). For these simulations, the values shown in Table 1 were assigned to the respective properties of the selected individual hydrocarbons. First-order biodegradation rate coefficients for the aqueous phase assigned to those compounds were calculated from the mean apparent first-order biodegradation rates applied to the gas phase, which were estimated by means of column experiments (Kjeldsen et al., 2003). Equation [11] was subsequently applied to different properties (i.e., vapor pressure, Henry's Law constant, aqueous solubility, octanol–water partition coefficient, biodegradation rate coefficient, and molecular gas diffusion coefficient). Our analysis revealed that the aqueous solubility was most instrumental in minimizing the objective function J(p). This result indicates that, for the particular field case examined here, selecting aqueous solubility as a grouping criterion leads to the smallest temporal variation in the composition of the composite constituents, and to the smallest temporal variation in their effective thermodynamic properties. The selection of aqueous solubility as the optimum grouping criterion is quite reasonable given the wide range of values among the fuel mixture components, and the fact that sorption and biodegradation take place only in the aqueous phase. The biodegradation process can create steeper concentration gradients around the contaminant source, resulting in enhanced dissolution of the organic phase. This effect is more pronounced for the compounds with larger mass in the aqueous phase.

Following the above grouping of the fuel compounds into composite constituents using the aqueous solubility as the grouping criterion, a number of simulations were performed using different fuel components as target contaminants to evaluate apparent biodegradation rates and assess the applicability of the constituent averaging technique to large-scale contamination problems. All parameters, except for the biodegradation rates, were assigned values that were either measured at the site or obtained from the literature (individual constituents) or calculated from member compound property values (Table 1) using the algorithm described above (composite constituents). The sorption coefficients were estimated from octanol–water partition coefficients (Table 1) and organic C fractions measured at the site using relationships reported in Schwarzenbach and Westall (1981). The longitudinal and transverse dispersivity values were assumed to equal 0.1 and 0.01 m, respectively (Gelhar et al., 1992). These parameters were assigned deterministic values to simplify the parameter estimation procedure and to assess the performance of the constituent averaging technique under conditions of reduced uncertainty and having a minimum number of adjustable parameters. The code requirements in random access memory (RAM) and storage were 640 kilobytes and approximately 3 megabytes, respectively. The computational time required for one simulation ranged from 30 min to 2 to 3 h on a Pentium IV, 1.9 gigahertz processor, depending on the involved convergence criteria.

Biodegradation rate constants of the simulated compounds were adjusted through calibration of the transport model, and subsequently compared with literature data where possible. The concentrations of hydrocarbons in the gas phase as measured at the Værløse site were used as calibration data for this inverse procedure. The minimization problem was solved using a bounded least-squares algorithm (McLaughlin and Townley, 1996):

[12]
where K is the biodegradation rate vector, m is the number of constituents (either single or composite) used in the simulations, Cgj is the experimental gas concentration value for the jth constituent, and fj(K) is the model prediction with regard to the jth constituent. The second term on the right-hand side is the prior probability density of the parameter vector calculated either from independent in situ measurements or from literature sources. The prior probability density is a probabilistic description of the range of values that is physically meaningful for a given parameter. The prior density pK(K) is a positive constant when K lies within the expected range of values, and 0 otherwise. This term enforces the constraints imposed by prior information on the parameters by assigning an infinite penalty on estimates that are outside the physically acceptable parameter range. The prior densities pK(K) of the first-order biodegradation rate constants of the simulated constituents were assigned log-uniform distributions with minima and maxima equal to one order of magnitude smaller and greater, respectively, than the average biodegradation values measured in laboratory experiments (Kjeldsen et al., 2003; Höhener et al., 2003).

In the first stage of model calibration, simulations were performed with biodegradation values randomly sampled from the prior parameter space. This stage identified the high probability area in the parameter space that constrains the global minimum of the minimization function, Eq. [12]. The second stage consisted of fine-tuning the model (using Eq. [12]) to locate the maximum likelihood biodegradation rate constants within the high probability area estimated in the previous step. The confidence intervals of the parameter estimates were calculated through a parameter space analysis. For a large number of measurements, assuming that (i) the uncertainties in the data are normally distributed and (ii) the model structure error is negligible, an approximate (1 – a) 100% confidence region for the parameter set can be estimated as (Weiss and Smith, 1998)

[13]
where p is the number of parameters, is the maximum likelihood parameter estimate, F(p, np) denotes an F distribution with p and (np) degrees of freedom, S2 = J()/(np), and n is the number of available measurements of the dependent variable.


    RESULTS AND DISCUSSION
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
Three sets of simulations were performed in this study using different fuel components as target contaminants. In all cases, model calibration resulted in gas concentration profiles that were in close agreement with those measured in the field. For the first set of simulations, toluene was selected as the target contaminant, while the remaining compounds were grouped (using aqueous solubility as a grouping criterion) into four composite constituents: (i) benzene, CFC-113, m-xylene and cyclopentane (composite constituent 1); (ii) 1,2,4-trimethylbenzene, methylcyclopentane, methylcyclohexane, 3-methylpentane and hexane (Composite Constituent 2); (iii) isooctane and octane (Composite Constituent 3); and (iv) decane and dodecane (composite constituent 4). Figure 6 shows the simulation results of the best match between the model and experiment for the five simulated constituents in terms of gas concentration at depth 1.05 m below the ground surface. Note that gas concentrations in this figure range over several orders of magnitude. For all simulated constituents, the first-order biodegradation model achieved reasonably good reproductions of the experimental concentration profiles. Results of simulations using the maximum likelihood biodegradation rate constants (best match of predictions with observations) are also compared with field measurements in terms of mole fraction evolution in the source zone, and total mass in the gas phase to assess the performance of the constituent averaging technique and to verify the validity of the model assumptions, such as equilibrium partitioning between phases.



View larger version (26K):
[in this window]
[in a new window]
 
Fig. 6. Reproduction of experimental data for the five simulated constituents in terms of gas concentrations at depth 1.05 m below the ground surface.

 
Figure 7 shows that the temporal evolution of the source composition (organic phase) as measured experimentally was satisfactorily approximated by the model, which was calibrated only to gas concentrations. These results suggest that, at least for our data, the assumption of equilibrium mass partitioning across the interfaces between the gas–aqueous, gas–NAPL, and aqueous–NAPL phases was reasonable. They also show that Eq. [10], which defines the composite constituents according to the constituent averaging technique, holds for both the oil and the gas phases.



View larger version (17K):
[in this window]
[in a new window]
 
Fig. 7. Simulated and experimentally determined mole fraction evolution of simulated constituents in the source (organic phase) with time.

 
Additional simulation results using hexane and CFC-113 as target contaminants were also consistent with the above findings. For the simulations that involved hexane as the target contaminant, the four composite constituents consisted of (i) toluene, benzene, CFC-113, m-xylene, and cyclopentane (composite constituent 1); (ii) 1,2,4-trimethylbenzene, methylcyclopentane, methylcyclohexane, and 3-methylpentane (composite constituent 2); (iii) isooctane and octane (composite constituent 3); and (iv) decane and dodecane (composite constituent 4). The four composite constituents in the case of simulating CFC-113 as the individual compound were: (i) toluene, benzene, m-xylene, and cyclopentane (Composite Constituent 1); (ii) 1,2,4-trimethylbenzene, methylcyclopentane, methylcyclohexane, 3-methylpentane, and hexane (Composite Constituent 2); (iii) isooctane and octane (Composite Constituent 3); and (iv) decane and dodecane (Composite Constituent 4). The experimental gas concentration contours for hexane and CFC-113 along the main radial transect at the Værløse site are presented in Fig. 1. For comparison purposes, the respective simulated concentration contours for the two simulations using hexane and CFC-113 as target contaminants and the maximum likelihood biodegradation rates are shown in Fig. 8 . In accord with previous simulations and studies, a comparison of Fig. 1 and 8 suggests that the assumption of first-order biodegradation kinetics was reasonable and satisfactorily reproduced the gas phase concentrations of the fuel mixture compounds in both cases. A further comparison of experimentally determined and simulated results is provided in Fig. 9 , which shows the total mass of the three target contaminants (toluene, hexane, and CFC-113) in the gas phase within the model domain as a function of time. In this figure, the mass calculated from measured concentrations is compared with the simulated results (mass of hydrocarbon in the gas phase) using the maximum likelihood biodegradation rates.



View larger version (30K):
[in this window]
[in a new window]
 
Fig. 8. Simulated concentration contours along the main radial transect at the Værløse site.

 


View larger version (16K):
[in this window]
[in a new window]
 
Fig. 9. A comparison of experimental data and model results in terms of the mass of three fuel mixture compounds in the gas phase as a function of time.

 
The estimated biodegradation values were compared with biodegradation rate constants obtained from independent measurements to assess the transferability of laboratory measurements of biodegradation rates to the field, and their usefulness in numerical studies. Figure 10 shows first-order biodegradation rates associated with the best match of the modeling results, the upper and lower limits of their respective 95% confidence intervals, and the range of values obtained from column measurements (mean values and 95% confidence intervals) using soil from the site (as reported by Kjeldsen et al., 2003, and Höhener et al., 2003). As can be seen in this figure, the biodegradation rates that provided the best match to the experimental data in this study were, in most cases, within the range of biodegradation rate values obtained from column experiments for soils from the Værløse site. In fact, biodegradation rates inferred from the field experiment were close to the lower bound of the laboratory column experiments, except for benzene. We note that in the present work, the parameters that usually introduce the highest uncertainty in field modeling, including the spill history (i.e., the total and relative amounts of mixture compounds present in the organic phase and their spatial distribution in the soil) were determined relatively accurately, so the uncertainty was substantially reduced. This fact enhances the reliability of the estimated biodegradation values in this study, which compared well with the experimentally measured values. The agreement could be further improved if one considers that the temperature at the field site was considerably lower than the temperature (26°C) at which the laboratory measurements were taken (Kjeldsen et al., 2003), given that the biodegradation rate is expected to increase with temperature.



View larger version (18K):
[in this window]
[in a new window]
 
Fig. 10. Estimated biodegradation rate values, 95% univariate confidence intervals, and range of biodegradation rate constants (error bars) obtained from independent laboratory measurements (column tests using soil from the site).

 
Sensitivity Analysis
Additional simulations were performed to investigate the sensitivity of the attenuation calculations to the biodegradation rate constant, which is among the most uncertain, yet most significant parameters for plume spreading. Our analysis here explores the predictive capabilities of the proposed model and its potential to assist decision making processes. Figure 11 shows a comparison between the experimental gas concentration contours of toluene and model predictions if the biodegradation rates are assigned the upper and lower bounds (Fig. 10) of the laboratory measurements. Figure 12 presents the predicted gas concentration contours of hexane for the maximum and minimum biodegradation rates as defined by the 95% percentile of the experimental values reported in the literature. The impact of biodegradation on VOC transport in the vadose zone can be also seen in Fig. 12. At early stages, increasing the biodegradation constant by a factor of 4 appears to affect only negligibly the gas phase concentration of hexane. This may be attributed to the effect of advective transport caused by pressure gradients created in the source area as a consequence of evaporation of the compounds after installation of the source. The lack of concentration sensitivity to the precise value of the biodegradation rate constant at the early stage of the contamination event, when concentration levels are still high, is very convenient for risk assessment and decision making. At later stages, the same increase of biodegradation constant leads to one to two orders of magnitude lower concentration levels. Nevertheless, this strong sensitivity refers to concentrations that already are small compared with those at the source and, hence, the impact on decision making may be relatively small. In general, the adoption of conservative values within the range of biodegradation constant measurements from column experiments provided better agreement between calculated concentration profiles and those measured at the site for the majority of the mixture compounds.



View larger version (33K):
[in this window]
[in a new window]
 
Fig. 11. Experimental gas concentration contours of toluene compared with model predictions for the minimum and maximum biodegradation rates as defined by the 95% percentile of the values reported by lab experiments.

 


View larger version (34K):
[in this window]
[in a new window]
 
Fig. 12. Model predicted gas concentration contours of hexane for the minimum and maximum biodegradation rates as defined by the 95% percentile of the values reported by lab experiments.

 

    CONCLUSIONS
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 
This study suggests that many problems in transport modeling within the vadose zone may be reasonably well addressed by simplified numerical models or approaches once the dominant factors and transport processes are identified. The appropriateness of the applied approximations and simplifications, however, should be judged on the basis of analyses of field data and case-specific sensitivity calculations. There may be an advantage in conducting preliminary parametric studies using a detailed model to identify important case-specific factors and to guide the formulation of a simpler, but still reliable model for use in risk assessment calculations.

Using data from a large-scale field experiment at Airbase Værløse, Denmark, our simulations show that an overpressure is created in the source area, which gives rise to significant advective transport in the gas phase. This advective transport greatly affected the spreading of the VOCs in the vadose zone, especially at the early stages of the field experiment. At later times, as the source was being depleted, the impact of ignoring gas phase advection in the simulation model weakened and the VOC transport process became more diffusion dominated. However, early-time gas advection may also be responsible for appreciably enhancing the escape of contaminants to the atmosphere. This can considerably influence the transport process at later times. Large temperature variations may also greatly affect VOC volatilization, dissolution, and biodegradation of mixture compounds. To deal with the large seasonal variations of temperature at the Værløse site, an external module was developed and coupled to the numerical code MOFAT, which then allowed for an update of the most temperature-sensitive properties of the fuel mixture components. These properties were adjusted using expressions for their temperature dependence obtained from the literature.

In all simulated cases, the constituent averaging technique provided close approximations to the temporal evolution of mixture composition, and of the effect of mixture composition changes on transport and natural attenuation of the selected mixture compounds. This approach appears to provide an attractive way to reduce computational costs when dealing with complex VOC mixtures in the vadose zone. The assumptions of equilibrium mass partitioning and first-order kinetics resulted in satisfactory estimates of the mass transfer rates between gas–aqueous, gas–NAPL, and aqueous–NAPL phases, and of the biodegradation rate of VOCs in the vadose zone.

The biodegradation constants that were estimated by the model were found to fall within the range of biodegradation rate constants determined by column experiments for the majority of mixture compounds. Such column experiments can provide valuable information to a more informative probabilistic description of biodegradation in numerical modeling. Unfortunately, biodegradation rates at the field scale are generally highly uncertain parameters, even when experimental data are available. Hence, biodegradation rate constants are attractive candidates for treatment as stochastic parameters. We are currently pursuing such an approach.


    ACKNOWLEDGMENTS
 
Funding for this project was provided by EC (project GRACOS-EVK1-CT-1999).


    REFERENCES
 TOP
 ABSTRACT
 INTRODUCTION
 MATERIALS AND METHODS
 MODELING
 RESULTS AND DISCUSSION
 CONCLUSIONS
 REFERENCES
 





This Article
Right arrow Abstract Freely available
Right arrow Figures Only
Right arrow Full Text (PDF) Free
Right arrow Alert me when this article is cited
Right arrow Alert me if a correction is posted
Services
Right arrow Similar articles in this journal
Right arrow Similar articles in Web of Science
Right arrow Alert me to new issues of the journal
Right arrow Download to citation manager
Right arrow reprints & permissions
Citing Articles
Right arrow Citing Articles via Web of Science (5)
Right arrow Citing Articles via Google Scholar
Google Scholar
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
Right arrow Search for Related Content
PubMed
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
GeoRef
Right arrow GeoRef Citation
Agricola
Right arrow Articles by Gaganis, P.
Right arrow Articles by Burganos, V. N.
Related Collections
Right arrow Volatile Organic Compounds
Right arrow Multicomponent Transport Models
Right arrow Organic Compounds


HOME HELP FEEDBACK SUBSCRIPTIONS ARCHIVE SEARCH TABLE OF CONTENTS
The SCI Journals Agronomy Journal Crop Science
Journal of Natural Resources
and Life Sciences Education
Soil Science Society of America Journal
Journal of Plant Registrations Journal of
Environmental Quality
The Plant Genome