Improved Approaches to Calculating Dilution and Attenuation Factor in Contaminant Hydrogeology ()
1. Conceptual Model for Developing Site-Specific Dilution and Attenuation Factor
Soil to groundwater contamination is a major global concern and occurs when a contaminant moves through the vadose zone into the underlying aquifer. Sources of pollution include but are not limited to solid waste management units (SWMU), waste disposal, leaking storage tanks, and fertilizer applications. U.S. Environmental Protection Agency (USEPA, 1996) provided the general framework and is still being used to determine appropriate soil screening levels, which dictate the concentration of a constituent of potential concern (COPC) in the unsaturated soil below which the given contaminant does not present a health concern from subsequent contaminant leaching into groundwater. The various fate and transport mechanisms such as dilution, adsorption, and degradation over the course of contaminants moving through unsaturated and saturated zones to receptors are represented by a critical parameter in contaminant hydrogeology, i.e., dilution and attenuation factor (DAF) (USEPA, 2023).
DAF is defined as the ratio of contaminant concentration in soil leachate to the concentration in groundwater at the point of withdrawal. Many environmental guidance and regulatory programs use DAF to estimate the impact of unsaturated zone mass discharge on the underlying groundwater (American Society for Testing and Materials, 2022; Texas Commission on Environmental Quality, 2022; Newell et al., 2022). A higher DAF indicates a greater degree of dilution and attenuation of contaminants along the migration flow path. Figure 1 shows a conceptual site model (CSM) of the potential pathways of a COPC from a SWMU to a receptor well as the point of exposure. The migration generally consists of three distinct stages:
Leach downward through the unsaturated zone (vadose zone) to reach the groundwater table;
Mix with laterally flowing groundwater within top of the saturated zone (aquifer); and
Transport laterally in the saturated zone (aquifer) to the downgradient receptor well.
Determination of the concentration at which the COPC might be found in the receptor well as a result of a release from the SWMU requires a determination of how much the released COPC is attenuated and diluted at each step, which is addressed through the calculation of appropriate DAFs in the unsaturated zone, the mixing zone, and the saturated zone. Mathematically, the COPC concentration at the receptor well can be expressed by:
(1)
where:
Cwell = concentration at water supply well;
Cmix = concentration leaving mixing zone;
Csource = aqueous concentration in the source area such as an SWMU;
Cwt = concentration entering water table at the bottom of the unsaturated zone;
DAFsaturated = DAF in the saturated zone (groundwater) in which the COPC exits the mixing zone and migrate to the receptor well;
DAFmix = DAF in the mixing zone where leachate from vertical infiltration mixes with the laterally flowing groundwater;
DAFunsaturated = DAF in the unsaturated zone through which the COPC in the source area leaches to the groundwater table.
(2)
Equation (2) describes the dilution and attenuation processes in three distinct zones. These processes should be evaluated for each specific site where one or two processes may dominate. Precipitation is the main source of water recharging the unsaturated zone through the process of infiltration. The infiltration involves both liquid and gas phases because the pores are partially filled with water and partially filled with soil gas. The multiphase nature in the unsaturated zone gives rise to capillary effects, causing each fluid phase to have differing local fluid pressures. Capillary forces affect the soil water content and the infiltration rate. The amount of water that infiltrates through the subsurface, in turn, has a direct impact on the amount of chemical mass that is transported in the aqueous phase toward groundwater, thus the concentration at groundwater table, Cwt.
Not only DAF can help predict the concentration of COPC at an extraction well, but it can also help determine the target soil leachate concentration or the target source concentration, which is calculated by:
Cw = (Maximum Contaminant Level [MCL] or another health-based regulatory limit) × DAF.
The target soil leachate concentration is then integrated into the following two questions to calculate the soil screening level (SSL) (USEPA, 1996):
For inorganic COPCs:
(3)
For organic COPCs:
(4)
where:
Cw = target soil leachate concentration (mg/L), which is calculated by Cw = MCL × DAF;
MCL = maximum contaminant level or another health-based limit;
Kd = Soil-water partition coefficient (L/kg);
Koc = Soil organic carbon-water partition coefficient (L/kg);
foc = Organic carbon content of soil (kg/kg);
θw = Water-filled soil porosity;
θa = Air-filled soil porosity;
H' = Henry’s law constant (dimensionless);
ρ = Dry soil bulk density (kg/L).
While other parameters can be determined with geotechnical analysis of soil samples, DAF is the parameter that needs to be calculated from other data. USEPA (1996) recommends a default value of 20 if no site-specific data are available for a 0.5-acre source area but also allows development of site-specific DAF. It is perceivable that contaminant transport through the saturated zone toward the receptor well will dilute and attenuate its concentration. We recommend that DAFsaturated be conservatively assumed to be 1. This paper focuses on approaches to calculating DAFunsaturated and DAFmix. We propose two improved methods to calculate the DAF.
2. Overview of Commonly Used Methods in Calculation of DAF in Unsaturated Zone
2.1. Analytical Solution
When limited data is available for DAF calculation, analytical solutions can be an alternative approach. The one-dimensional chemical transport taking advection, dispersion, retardation, and biodegradation into account can be described by the following partial differential equation (Javandel et al., 1984):
(5)
where:
D = hydrodynamic dispersion coefficient (cm2/year);
v = water infiltration rate (cm/year);
C = contaminant concentration (mg/L);
x = distance along flow path (cm);
t = time (year);
λ = chemical decay constant (year−1);
R = chemical retardation factor (unitless).
The retardation factor R is a ratio between water flow rate and contaminant transport rate and is expressed as:
(6)
where:
Kd = soil-water partition coefficient (mL/g);
ρ = soil bulk density (g/cm3);
θ = soil porosity (cm3/cm3).
Under uniform flow conditions the analytical solution to the differential Equation (4) is:
(7)
The above solution is valid under the following initial and boundary conditions:
Initial Conditions (t = 0): C = C0 when 0 < x < A0 and C = 0 when x > A0, where A0 is the initial contamination thickness (Figure 1).
Boundary Condition (t > 0):
when x approaches infinite.
Figure 1. Conceptual model for developing site-specific DAF.
Although both chemical adsorption and chemical decay reduce chemical migration and increase dilution, there are significant uncertainties in estimating these parameters. For conservative purpose, no chemical adsorption and decay are assumed, i.e., R = 1; λ = 0. Equation (7) is then simplified into the following equation (Enfield et al., 1982):
(8)
Equation (8) can be solved to determine the maximum COPC concentration leaching into the groundwater when the peak of the COPC pulse (or plume) arrives at the groundwater table. In any instance, the initial concentration C0 is unknown in a disposal unit. Instead, soil concentration is available. The following soil leachate equation can be used to estimate C0:
(9)
where:
Cs = contaminant concentration in soil (mg/kg);
C0 = contaminant concentration in leachate in disposal unit (mg/L);
Kd = soil/water partition coefficient (L/kg);
θw = water-filled soil porosity;
θa = air-filled porosity;
H' = Henry’s law constant;
ρb = dry soil bulk density (kg/L).
The time it takes for the front of the plume to reach the groundwater table is:
where:
A is = thickness of the unsaturated zone (Figure 1).
The dispersion coefficient is the product of water flow rate v and dispersivility α:
Substituting t and D into Equation (9) yields the following reduced form for the peak concentration at the water table:
(10)
Cpeak represents the maximum concentration at the groundwater table. Subsequently, the smallest DAFunsaturated is calculated by:
(11)
As shown in Equation (11), the minimum DAF in the unsaturated zone is a function of three variables:
Initial depth of contamination at the SWMU;
Depth to the water table;
Dispersivity coefficient.
Gelhar et al. (1992) suggested that the observed dispersivity under field conditions was on the order of 10% of the flow length. Therefore, Equation (11) can be further simplified as:
(12)
Alternatively, Equation (10) can be used to determine if the soil to groundwater pathway is complete. For example, if Cpeak is less than the federal drinking water standards, then the soil to groundwater pathway is incomplete.
2.2. Dilution Method
USEPA (1996) presented four water balance models for dilution in the mixing zone of the aquifer for DAF calculations. Although written in different terms, all four models can be expressed by the following form:
(13)
where:
K = aquifer hydraulic conductivity (m/year);
i = hydraulic gradient (m/m);
Sd = mixing zone depth (m) (Figure 1);
I = infiltration rate (m/year);
L = length of source parallel to groundwater flow (m) (Figure 1).
The following equation is recommended to calculate the mixing zone depth:
(14)
where:
B = aquifer thickness (m) (Figure 1).
Such a development of the SSLs considers only the dilution of contaminant concentration through mixing with groundwater in the aquifer that can be assumed to be unconfined, homogeneous, and isotropic. For fractured aquifers, assumption justifications are to be provided for use of this approach. Because the hydraulic conductivity of fractured rock is typically much smaller than that of the overlying alluvium, most of the mixing may have occurred on top of the bedrock in the alluvium. This approach also has the following additional conservative assumptions:
The required input parameters that are required in calculation of DAFmix in Equations (13) and (14) are site-specific. The source area length can be based on historical records or measured directly. Aquifer hydraulic conductivity and hydraulic gradient, and aquifer thickness can be calculated from aquifer testing and groundwater monitoring data at the SWMU.
The infiltration rate is not readily measurable, especially for the infiltration rate at the groundwater table. Although the infiltration rate can be estimated with in-situ monitoring techniques including lysimeter, tensiometer, and ring infiltrometer (U.S. Department of Defense Environmental Security Technology Certification Program [ESTCP, 2022], 2022), they are more often estimated from water or mass balance. Several methods of estimating the infiltration rate are discussed in the following sections.
Water balance method
Water balance models couple climatic and hydrological data with a simplified model for estimating the infiltration rate. HELP (Hydrologic Evaluation of Landfill Performance) model is a commonly used water balance model. HELP is a layered, water budget (moisture routing) model for hydrologic evaluation of landfill performance, but can be applied more generally to evaluate infiltration rate with the following information:
Precipitation;
Evapotranspiration based on leaf area index, growing season length, evaporative depth, wind speed, and humidity;
Soil properties such as porosity, water content, saturated conductivity;
Factors affecting run-off including vegetation cover and slope.
Empirical method
Infiltration rate can be estimated from rainfall measurement data. In arid and semiarid regions, Woods (1999) suggests that the infiltration ranges from 2% to 4% of average precipitation and often is focused in playas, arroyos and topographic depressions. Data from 101 study sites compiled by American Petroleum Institute (1996) are plotted in Figure 2 (Newell et al., 2022).
Figure 2. Annual groundwater recharge as percentage of annual precipitation.
The best linear regression of these data (red line) indicates that relative recharge as a percent of annual precipitation can be found by multiplying annual precipitation in mm/year by 0.00017, i.e.:
where:
I = mean annual net infiltration (cm/year);
P = mean annual precipitation (cm/year).
The study sites are biased toward sandy soil in arid and semiarid regions. Connor et al. (1997) also provide empirical estimates of the infiltration rate for silt and clay, respectively:
For silt soil:
For clay soil:
Environmental tracer method
Conservative environmental tracers can be used to estimate the infiltration rate. Some example tracers are tritium, bromide, and chloride (Dassi, 2010). This method has been recognized by some researchers as the most successful method for estimating recharge in arid regions (Allison et al., 1994; Phillips, 1994). The assumptions made for this estimate are:
The natural tracer in the groundwater originates from precipitation;
The natural tracer is conservative in the system;
The tracer mass flux has not changed over time.
The following equation is used to calculate the infiltration rate:
(15)
At each study site the meteoric tracer concentration in precipitation does not change much. However, the tracer concentration in groundwater or unsaturated zone may change from one SWMU to another, depending on their spatial relation to preferential flow paths and the selected water samples. At Fort Wingate of New Mexico, the infiltration rate is calculated by yearly precipitation multiplied by the ratio of chloride concentration in rainfall over chloride concentration in groundwater (Henry et al., 2016). The calculated infiltration rate was 0.0000178 m/year.
Equation (15) oversimplifies the infiltration processes in the unsaturated zone where multiple layers of soil are present in the vertical profile. To accommodate different soil types at different depths, the accumulative tracer method is recommended to determine the recharge rate. In this method, the tracer concentration and water contents at a certain depth in a specific profile are accumulated to determine the multi-year average recharge capacity of the vertical infiltration. The water flux recharged by precipitation (and water discharges) can be expressed in the following equation for a unit area:
(16)
where:
I = Infiltration rate (mm/year);
TMp = yearly tracer mass in precipitation (mg/year);
WMs,z = accumulated amount of water in the unsaturated zone from surface to depth z (mm);
TMs,z = accumulated tracer mass (mg) in the unsaturated zone from surface to depth z.
The accumulated amount of water and tracer mass can be calculated using the flowing equations:
(17)
(18)
where:
n = number of soil layers;
hi = thickness of layer i (mm);
wi = mass water content of soil in layer i (%);
bi = dry density of soil in layer i (g/cm3);
ρi = density of water in layer i (g/cm3), which is set as 1 g/cm3;
Ci = tracer concentration of soil water in layer i (mg/l);
θi = volume water content of soil in layer i (%).
The accumulative environmental tracer method addresses the influence of soil characteristics on infiltration. Based on the soil characteristics, soil samples can be collected at discrete depths of the unsaturated zone and be analyzed for soil density, moisture content and tracer concentration. These data can then be used to calculate the infiltration rate. This method applies to sites where the unsaturated zone is relatively thick, and the soil profile consists of different layers.
van Genuchten method
The hydraulic conductivity of soil varies with saturation, with the maximum value occurring when the soil is saturated. The unsaturated hydraulic conductivity is a function of water content. Because the hydraulic gradient in the vadose zone is approximately 1, the maximum infiltration rate equals the unsaturated hydraulic conductivity, which can be calculated by:
(19)
where:
= unsaturated hydraulic conductivity;
Ksaturated = hydraulic conductivity in vertical direction at saturation;
= relative hydraulic conductivity at soil water content level of θ.
A closed-form analytical solution to the relative hydraulic conductivity is provided by van Genuchten (1980) for various soil water contents:
(20)
where:
θ = soil water content;
θs = soil water content at saturation;
θr = residual soil water content;
N = van Genuchten parameter.
The values of θ, θs, θr, Ksaturated, and N can be provided by laboratory analysis of soil samples collected at each site. The water content data can also be measured with in-situ instrumentation such as time domain reflectometry or neutron probes. The van Genuchten parameter N is derived from a water retention curve, as shown in Figure 3, of a soil sample.
If Ksaturated is known, the shortest travel time for a contaminant to reach the groundwater table can be calculated by:
(21)
Figure 3. Example predicted water retention curve and data points for a soil sample.
where:
T = travel time from source to groundwater table (year).
The travel time for contaminant percolating through the unsaturated zone could provide additional insight on the significance of impact to groundwater. If there is a thick unsaturated zone and calculations indicate relatively long travel time, there may be limited potential for significant contaminant flux to groundwater. This is particular true when there are mechanisms such as biodegradation and sorption that will result in attenuation of contaminants.
2.3. Numerical Models
Many numerical models, using either finite difference or finite-element method, have been developed and used to simulate the fate and transport of COPCs in the unsaturated zone. These models are applicable to evaluating the soil to groundwater pathways, estimating DAF and then SSLs based on simplifications and assumptions. Numerical models require extensive data input, and much of the data may not be available for some project sites. Therefore, applicability of a numerical model to a SWMU depends on the site-specific scenario and comparison of the data available (or potentially available) against the input requirements for the model. Table A1 of Appendix summarizes five commonly used numerical models and their data requirements. Although each of the five models has been reported for simulation of water and contaminant transport in the unsaturated zone, they emphasize different transport mechanisms in calculating the DAF.
The five unsaturated models evaluated herein can calculate the leachate concentrations entering groundwater (Cwt in Figure 1). The leachate concentrations are needed to develop the site soil SSLs and to estimate groundwater concentrations at the receptor well. These concentrations at the receptor point are then compared with the acceptable groundwater concentration (e.g. MCL). If they do not exceed the acceptable groundwater concentrations, it might be believed that there is no complete pathway at the site. The conclusion of such a comparison is based on the following assumptions:
A CSM is reasonably developed.
Site-specific data required by the models are properly collected and meet data quality standards.
Unsaturated zone model for COPC migration is properly used.
All five numerical models (SESOIL, HYDRUS, CHAIN 2D, MULTIMED_DP, and FECTUZ) are capable to estimate the soil to groundwater SSLs. When there is a concern for the uncertainty of the SSL estimates due to variability of the input parameters, sensitivity analyses can be conducted using the built-in Monte Carlo simulation in the MULTIMED_DP and FECTUZ, whereas the sensitivity analysis using a set of variable input parameters is recommended for SESOIL, HYDRUS and CHAIN 2D models.
While a few tools are available for solute transport modeling through the unsaturated zone, limitations and uncertainties in these models must be recognized, as summarized in Table A1 of Appendix. Furthermore, these models deal with porous media and are not readily applicable to fractured bedrocks. Because of uncertainties in hydraulic conductivity and non-linear relationships between hydraulic conductivity and soil suction, estimate of seepage rate is always approximate. Significant uncertainties are associated with water flow modeling in the arid climate where precipitation and evapotranspiration can vary dramatically daily.
The analytical and numerical models are understandably more applicable to the alluvium. For the unsaturated fractured rock, it can be treated as a separate layer, and the DAF is often conservatively assumed to be 1. The data requirements for numerical models are extensive. Some of the data may not be obtainable in practice to calculate the site-specific DAF.
3. Probability Method in Calculation of DAF
As an alternative to the dilution method, a probability method is proposed to estimate the order of magnitude of the DAFmix. It should be noted that the probability method is limited to the data used by USEPA in developing the SSL Guidance (USEPA, 1996) based on the following assumptions:
An infinite COPC source, with no chemical adsorption to soil and no chemical degradation.
The nearest drinking water wells are within 100 feet downgradient edge of the SWMU.
Wells are assumed to be screened within 15 to 300 feet beneath the groundwater table.
Because the USEPA adopted an infinite source and no chemical adsorption to soil, the derived groundwater DAFs implicitly exclude any dispersion/dilution within the unsaturated zone. The assumption that the receptor is within 100 feet downgradient edge of the waste unit excludes any lengthy transport processes in the saturated zone. Thus, the DAF values represent primarily the effects of mixing and dilution in the aquifer underlying the SWMU and its vicinity of less than100 feet. Because potential receptor wells are typically located miles downgradient of any SWMU, it is reasonable to assume that the DAFs derived by USEPA (1996) mostly describe the processes that would be included in DAFmix.
USEPA derived groundwater DAFs for a wide range of climatological and hydrogeological conditions encompassing the country. In order to address the widely varying conditions, the USEPA used a Monte Carlo framework coupled to a chemical fate and transport model. The framework was implemented by selecting a source area and then randomly selecting inputs for the fate and transport model repeatedly to produce a distribution of DAF values for a given contaminant release area. This procedure was repeated for a range of source areas from 0.02 to 69 acres leading to a family of DAF distributions for the different source sizes. USEPA reported the 85th, 90th, and 95th percentile lowest values from these distributions in table 5 of its SSL Guidance (USEPA, 1996). Although USEPA considered a range of source areas in its SSL guidance, it did not develop relationships between the parameters of DAF distributions (i.e., mean and standard deviation) and size of the source area. These relationships are needed to implement the probabilistic framework utilized in our calculation of DAFmix.
The three DAF percentile values reported by USEPA (85th, 90th, and 95th) were used to estimate the mean and standard deviation of the probability distribution of DAF as a function of source size. Because the lower bound of a lognormal distribution is zero, whereas the minimum value of DAF is 1, the transformed variable (DAF−1) was fit to the lognormal distribution (Gradient Corporation, 2013):
(22)
where:
Y ≤ y = cumulative probability of any value y from the distribution of the random variable Y;
y = transferred variable (DAF−1);
µ = mean of ln(DAF−1);
σ = standard deviation of ln(DAF−1).
Figure 4 shows the close correlations between the USEPA-derived DAF distribution and the calculated from the lognormal distribution for the 85th, 90th, and 95th percentiles.
(a)
(b)
(c)
Figure 4. Correlation between probability-based DAF and USEPA calculated DAF. (a) Correlation between probability-based DAF and USEPA calculated DAF at 85th percentile, (b) Correlation between probability-based DAF and USEPA calculated DAF at 90th percentile, (c) Correlation between probability-based DAF and USEPA calculated DAF at 95th percentile.
By sequentially fitting the mean and coefficient of variation (CV is the standard deviation divided by the mean, or σ/µ) to each percentile and area, we derived a best fit polynomial for µ and CV as a function of source area. The resulting polynomial equations for each are given below:
(23)
(24)
where:
x = log10(area) for source area in acres;
µ = mean of ln(DAF−1);
CV = coefficient of variation of ln(DAF−1).
Under the lognormal distribution, the DAFmix can be calculated for any source size and percentile levels by the following equation:
(25)
where:
Zα = Z score of a standard normal distribution (mean of 0 and standard deviation of 1) corresponding to α percentile.
Figure 5 shows the relationship between probability-based DAF and source sizes. The Z score can be found in typical statistics textbooks. Table 1 presents the Z scores for the 85th, 90th, and 95th percentiles.
Table 1. Z-scores in response to percentile for a standard normal distribution.
Percentile |
85th |
90th |
95th |
Z score |
−1.0364 |
−1.282 |
−1.645 |
Figure 5. Probability-based DAF versus source size.
4. Application of Probability Method to Calculating DAF
Figure 6 shows the study site where the two proposed approaches were applied. Groundwater occurs in fractured andesite bedrock at a depth of 37 m below ground surface under unconfined conditions. The vadose zone includes coalescent alluvial fan deposits of approximately 24 to 30 m thick, and the upper 6 to 12 m of fractured bedrock above the water table.
Figure 6. Layout of case study site.
The SWMU was a wastewater impoundment. Constituents dissolved in waste-water disposed to the SWMU were transported downward into the vadose zone as water infiltrated through any compromised areas of the lined impoundments. Following infiltration, wastewater and its dissolved constituents would percolate primarily downward through the vadose zone, with limited lateral spreading due to heterogeneities in soil texture and the absence of any laterally continuous soil horizons that could present a barrier to infiltration. During downward migration through the vadose zone, transport of some dissolved constituents would be retarded by sorption to mineral surfaces and, for hydrophobic organic constituents, sorption to natural organic carbon present in the native soil. The wastewater would ultimately reach the groundwater zone and recharge the aquifer. Hence, constituents that were historically discharged to the wastewater impoundments may have contributed to groundwater contamination.
Although there are groundwater monitoring wells in the vicinity of the SWMU, as shown in Figure 6, none of the soil borings encountered groundwater. The SSLs should be developed to evaluate observed soil concentrations and the potential for these concentrations to adversely impact groundwater underlying the soil. The site-specific DAF is the most important parameter to develop the SSLs.
Table 2. Summary of geotechnical properties.
Table 3. Summary of input parameters and calculated DAF values.
Dilution Method |
Probability Method |
Parameter |
Definition |
Value |
Source for the Value |
Parameter |
Definition |
Value |
Source for the Value |
K |
Aquifer hydraulic conductivity (m/year) |
70.7 |
Slug tests of monitoring wells located adjacent to the study site. |
Source length and source width |
Source area (acre) |
0.98 |
Figure 6 , field measurements |
i |
Hydraulic gradient (m/m) |
0.059 |
Figure 6 , potentiometric surface map from groundwater level measurements in monitoring wells. |
x |
Log10(source area) |
−0.00942 |
Equation (23) |
D |
Mixing zone depth (m) |
14 |
Calculated from Equation (14). The aquifer thickness used for these calculations is 73 m. |
μ |
Mean value |
16.73 |
Equation (23) |
I |
Infiltration rate (m/year) |
0.0067 |
Calculated from van Genuchten method |
CV |
Coefficient of variation |
0.58 |
Equation (24) |
L |
Source length parallel to groundwater flow (m) |
132 |
Figure 6 , field measurements |
σ |
Standard deviation |
9.68 |
μ × CV |
W |
Source width perpendicular to groundwater flow (m) |
30 |
Figure 6 , field measurements |
|
|
|
|
DAF |
Calculated DAF |
67 |
Equation (13) |
DAF |
Calculated DAF |
812 at 85% probability percentile 76 at 90% probability percentile 3 at 95% probability percentile |
Equation (25) |
Soil samples were collected from 10 soil borings for chemical analysis. Based on soil sampling results, the site-specific COPCs included trichloroethene, tetrachloroethene, and their daughter products. Additional soil samples were also analyzed for physical and hydraulic properties, which are referred to in this paper as geotechnical properties. Results of analysis of geotechnical soil samples are presented in Table 2.
Based on the geotechnical soil analysis results, especially the values of hydraulic properties, the van Genuchten method was used to calculate the unsaturated hydraulic conductivity or the infiltration rate. It is assumed that the soil beneath the site is draining under gravity conditions. The advancement of a sharp wetting front that might cause steep pressure head gradients was not observed. Soil water contents consistently less than 10% provided evidence that large pressure head gradients were unlikely to exist at the time of sample collection. In addition, the samples were collected in the months immediately following the monsoon season, when soil moisture conditions are expected to be near their highest. As a result, the maximum unsaturated hydraulic conductivity of 6.7E−03 m/year (Table 2) was used as a conservative estimate of the infiltration rate beneath the site.
Table 3 compares the input parameters, values, justification for values, and the calculated DAFs from both the dilution and probability methods. The calculated DAF from the probability method varies with the probability percentile. At 90% probability percentile, the calculated DAF is 76, while the calculated DAF are 812 and 3 at 85% and 95% probability percentile, respectively. The calculated DAF is 67 from the dilution method, which is comparable to the probability-based DAF of 76 at 90% probability percentile. Based on the calculated site-specific DAFs, the SSLs were developed for the COPCs and used for evaluation of potential migration pathways from soil to groundwater.
5. Conclusion
DAF is a critical parameter in contaminant hydrogeology to evaluate complete or incomplete pathways from soil to groundwater. Both analytical and numerical methods have been used to calculate site-specific DAFs based on assumptions that may not be justifiable. This paper proposed two empirical methods, dilution method and probability method, to calculate DAFs from variables that are readily available from field measurements and geotechnical analysis of soil samples. In the dilution method, the critical parameter is the infiltration rate, which is often not available. The van Genuchten method is recommended, and the infiltration is represented by the unsaturated hydraulic conductivity, which can be calculated from hydraulic properties measured in geotechnical laboratories on soil samples. In the probability method, the critical parameter is the acreage of the SWMU size. Both methods focus on the DAF in the mixing zone at the groundwater table but include parameters reflective of the transport processes in the vadose zone. These two methods were applied to an actual SWMU, and the site-specific DAF was calculated, 76 from the dilution method and 67 from the probability method at 90% probability percentile. Because the calculated DAF are 812 and 3 at 85% and 95% probability percentile, respectively, the calculated DAF from the dilution method is comparable to the probability-based DAF at 90% probability percentile. Based on the calculated DAFs, the SSLs were developed for the COPCs and used for evaluation of potential migration pathways from soil to groundwater.
Author Contribution
Xiangquan Li: conceptualization; methodology; supervision; writing review and editing. Chunchao Zhang: conceptualization; methodology; visualization; writing review and editing. Wanfang Zhou: investigation; methodology; data acquisition; draft preparation.
Data Availability
The data presented in this study are available from the communicating author.
Appendix
Table A1. Summary of select numerical models for the unsaturated zone.
Model |
Main features |
Data Requirements |
Soil properties |
Site characteristics |
COPC properties |
Other data |
SESOIL |
One-dimensional model to simulate flow and transport from the land surface to the water table using finite difference method. First-order chemical degradation, biodegradation, cation exchange, hydrolysis, equilibrium partitioning to soil, and volatilization. Including hydrologic, sediment, and COPC fate cycles. Consisting of up to four soil layers with each layer dividing into 10 uniform sub-layers. Using either monthly or annual data. |
Number of layers and sublayers Thickness of layers Pore disconnectedness index Effective porosity Organic carbon content Silt, sand, and clay fractions Soil loss ratio pH of each layer Hydrolysis constants (acid, base, or neutral) |
Mean air temperature Mean cloud cover fraction Mean relative humidity Total precipitation Mean storm duration Number of storm events Bulk density Intrinsic permeability |
Cation exchange capacity Freundlich exponent Solubility in water Air diffusion coefficient Henry’s law constant Organic carbon adsorption ratio Soil adsorption coefficient Molecular weight Biodegradation rates (liquid, solid) |
Spill index COPC load Mass removed or transformed Index of volatile diffusion Index of transport in surface runoff Erodibility factor Practice factor Manning coefficient |
HYDRUS |
1-D model to simulate solute and heat transport in variably saturated media using finite-element method. Incorporating diffusion, hydrodynamic dispersion, linear equilibrium reactions between the liquid and gaseous phases, nonlinear non-equilibrium partitioning (sorption) between the solid and liquid phases, and first-order decay/degradation. Accounting for water uptake by plant roots with a sink term. Allowing temporal variations in flow, and heat and solute transport along boundaries. Simulating finite sources and hysteresis in the water movement. Soil properties being described by the van Genuchten parameters. |
Number of soil materials Depth of soil layers Saturated water content Residual water content Saturated hydraulic conductivity Soil bulk density van Genuchten retention parameter, α van Genuchten retention parameter, β Rescaling factors for hydraulic properties |
Uniform or stepwise rainfall intensity Volumetric fraction of solid phase Volumetric fraction of organic matter Thermal conductivity Heat capacities of solid phase, organic matter, and liquid phase Number of solutes Contaminant concentrations in soil |
Molecular diffusion coefficient Dispersivity Freundlich isotherm coefficients Freundlich isotherm exponents First order rate constants Decay coefficient |
Potential transpiration rate Osmotic coefficient Pressure head at 50% transpiration Root density as a function of depth Power function in stress- response function |
MULTIMED_DP |
1-D steady-state flow with a semi-analytical solution to simulate contaminant migration from a SWMU with an option for unsaturated zone transport. Seasonal variability in precipitation and evapotranspiration in inputs. Modeling the following processes: advection, dispersion, linear or nonlinear sorption, volatilization, hydrolysis, biodegradation, and first-order chemical decay. Addressing finite or infinite sources. |
Number of physical flow layers Thickness of each layer Number of porous materials van Genuchten retention parameter, α van Genuchten retention parameter, β Soil bulk density Residual water content Temperature of layer |
Area of waste disposal unit Length scale of facility Width scale of facility Saturated hydraulic conductivity Recharge rate Reference temperature for air diffusion Air entry pressure Depth of unsaturated zone |
Duration of pulse Source decay constant Initial contaminant concentration at landfill Longitudinal dispersivity pH of layer |
|
FECTUZ |
1-D fate and transport model to simulate migration of contaminants from a SWMU through the unsaturated zone to an unconfined aquifer. Allowing for finite or infinite sources which may vary with time. Modeling linear and nonlinear adsorption and first order decay. |
Soil bulk density Saturated water content Saturated hydraulic conductivity Residual water content van Genuchten retention parameter, α van Genuchten retention parameter, β Fraction of organic carbon |
Thickness of unsaturated zone Uniform thickness for discretized soil layers Uniform infiltration rate except for surface impoundments Constant source or Decaying source or finite pulsed source |
Organic carbon partition coefficient Freundlich isotherm coefficients Dispersivity Decay coefficient (dissolved) Decay coefficient (adsorbed) |
|
CHAIN 2D |
2-D model to simulate variably saturated flow, contaminant transport, and heat transport using finite element method. Accounting for water uptake by plant roots with a sink term. Incorporating soil anisotropy. Including prescribed head, gradient, flux boundaries, or free drainage for boundary conditions. Modeling following processes: advection, dispersion, conduction, convection, non-linear non-equilibrium reactions between solid and liquid phases, linear equilibrium, reactions between liquid and gaseous phases, and two first-order decay reactions: one independent of other solutes, and one in sequential chain decay reactions. |
2-D cell discretization Saturated water content Saturated hydraulic conductivity Soil bulk density van Genuchten retention parameter, α van Genuchten retention parameter, β Residual water content |
Transpiration rate Evaporation/ infiltration rates Initial contaminant concentrations in soil Contaminant species initial and boundary conditions Initial head conditions Location and rates of pumping/ injection wells Seepage faces, tile drains |
Ionic or molecular diffusion coefficient in water and gas phases Longitudinal and transverse dispersivities First order decay coefficient in liquid, solid or gas phase Zero order rate constant in liquid, solid or gas phase Adsorption (Freundlich) isotherm coefficients Source Decay |
Root density as a function of depth Power function in stress- response function Pressure head where transpiration is reduced by 50% |