Search tips
Search criteria 


Logo of nihpaAbout Author manuscriptsSubmit a manuscriptHHS Public Access; Author Manuscript; Accepted for publication in peer reviewed journal;
J Air Waste Manag Assoc. Author manuscript; available in PMC 2010 March 10.
Published in final edited form as:
J Air Waste Manag Assoc. 2009 April; 59(4): 447–460.
PMCID: PMC2835946

Development and Application of a Three-Dimensional Finite Element Vapor Intrusion Model

Kelly G. Pennell, assistant professor (research), Ozgur Bozkurt, doctoral candidate, and Eric M. Suuberg, professor


Details of a three-dimensional finite element model of soil vapor intrusion, including the overall modeling process and the stepwise approach, are provided. The model is a quantitative modeling tool that can help guide vapor intrusion characterization efforts. It solves the soil gas continuity equation coupled with the chemical transport equation, allowing for both advective and diffusive transport. Three-dimensional pressure, velocity, and chemical concentration fields are produced from the model. Results from simulations involving common site features, such as impervious surfaces, porous foundation sub-base material, and adjacent structures are summarized herein. The results suggest that site-specific features are important to consider when characterizing vapor intrusion risks. More importantly, the results suggest that soil gas or subslab gas samples taken without proper regard for particular site features may not be suitable for evaluating vapor intrusion risks; rather, careful attention needs to be given to the many factors that affect chemical transport into and around buildings.


Contaminated site closure criteria are often determined by risk assessments that evaluate the human health risks associated with contaminated groundwater and soil. The risks associated with vapor intrusion into inhabited structures have received increasing attention in recent years. Vapor intrusion involves indoor air contamination as a result of chemical volatilization from subsurface sources into structures located near contaminated sites. These vapor intrusion risks have been the subject of controversy, partly because there is a lack of scientifically based guidance on how to properly characterize the vapor intrusion pathway. For instance, Donohue1 notes that there is still a need to establish the value of specific vapor intrusion characterization techniques, and that vapor intrusion does not always represent a problem even in the presence of high contaminant vapor concentrations in a soil.

Gaining a better understanding of how vapors migrate from source to structure requires systematic scientific investigation and progress would be of great value in removing a critical source of uncertainty from the policy debate. Developing tools for evaluation of specific site features that may play an important role in characterizing vapor intrusion pathways is the focus of this research. The paper presented here describes the development of a quantitative model that can help guide both site investigation and mitigation efforts.

Vapor intrusion has been the focus of research, debate, and public policy for decades. Until the early 1990s, most of the discussion centered on intrusion of radon and methane. More recently, the discussion of vapor migration has shifted to concern for volatile organic compounds (VOCs), semi-volatile organic compounds, and in some cases metals (e.g., mercury) that are present as contaminants at contaminated sites.

During the radon era, an important landmark study in the area of vapor intrusion was performed by Nazaroff et al.2 The results of this study demonstrated that radon transport into a building from the soil was normally pressure driven. This finding confirmed the probable pathway for radon intrusion into buildings was entry through cracks or joints that connected the indoor air with the soil gas beneath the structure. On the basis of this conceptual understanding of the vapor intrusion process, several mathematical models were developed to predict radon transport through the subsurface into structures.26 More recently, Wang and Ward7,8 also have applied a three-dimensional (3-D) computation fluid dynamic (CFD) radon transport model.

Most of the initial radon models were simple one-(or two-) dimensional models that only considered advective (i.e., actual pressure-driven) transport of radon vapors into buildings through a single, well-defined pathway (e.g., a foundation crack). These models were later expanded to include diffusion (as a secondary transport mechanism), multiple cracks, and some simple geological variations in close proximity to the building foundation. The most advanced 3-D CFD radon models were developed by Loureiro et al.5 (modified by Revzan6) and the model later developed by Wang and Ward.7,8 These 3-D models illustrate that the quantitative prediction of transport of radon into a structure requires accounting for a multitude of effects and variables, and some considerable sophistication in the modeling geometry.

In contrast with the fairly advanced 3-D transport models developed by Loureiro et al.5 and Revzan6 in the early 1990s to describe radon transport, simpler one-dimensional (1-D) models to describe the transport of VOCs in the subsurface were first proposed in 1990 by Jury et al.9 and in 1991 by Johnson and Ettinger.10 Despite their spatial simplicity, both of these models included processes not previously included in most radon models, such as biodegradation9 (i.e., reduction of VOCs or semi-volatile organic compound concentrations in the soil by bacterial action, something that has no analogy in radon modeling) and contaminant migration from a source that is located a given distance from the building10 (which is different from many radon models that assume radon emanates directly from the soil surrounding the foundation).

The Johnson and Ettinger (J&E) model10 is at present the best-known and referenced model in the area of vapor intrusion. It was proposed as a screening tool in the U.S. Environmental Protection Agency (EPA) Draft Vapor Intrusion Guidance.11 Prior and subsequent to issuance of that guidance, the J&E model has come under considerable scrutiny. Several useful articles have been written to clarify proper application of the J&E model as a screening tool for identifying the potential for vapor intrusion problems,12,13 rather than its use as a quantitative predictive model, and it has been judged to be quite conservative in many cases.14

Several other 1-D vapor intrusion models have also been developed and used for various purposes.1520 An EPA document21 and other recent publications22,23 summarize some of the earlier models. It is generally acknowledged that 1-D models cannot capture the full range of effects that play a role in actual vapor intrusion scenarios, but they are useful because they can be more simply solved than full 3-D models.

It is the significant theoretical background available from earlier radon intrusion models,28 what has been learned from the J&E model (and other 1-D VOC vapor intrusion models mentioned previously), and recent two-dimensional (2-D) approximations24 that have led to development of more robust 3-D models for vapor intrusion of other contaminants, such as that by Abreu and Johnson25,26 and the model developed as part of the research reported herein. It is through the use of these 3-D CFD models that actual soil gas concentration profiles around a structure are available. It should be noted that Robinson and Turczynowicz23 also developed a 3-D model, but it is restrictive in geometry (axisymmetric). It is solved by Laplace transform methods, which make it computationally more straightforward than the other 3-D models, but this comes at the expense of being able to capture the full range of potentially important features.

To date, there has been little discussion about how relevant site features, such as impermeable surfaces (parking lots, sidewalks, foundations of detached garages, etc.), or backfill material beneath a foundation, could affect vapor intrusion. The research being performed by the authors is specifically aimed at understanding how various site features can affect vapor intrusion, and some examples are presented in this paper. Problems in risk assessments that may result from soil gas concentrations measured during typical site characterization activities are also considered here (and this will also be discussed in future papers by the authors of this manuscript).

Although one VOC 3-D vapor intrusion model is already available to investigate vapor intrusion fate and transport processes,25 better progress in understanding the many complexities of vapor intrusion can be achieved by testing multiple models against one another, until all can be validated with properly conducted field tests, which are difficult and expensive to conduct. In addition to serving as a comparison for previous modeling efforts by others, the model described herein is based upon readily available commercial computational platforms, which offer access to a broader range of users.

A commercially available CFD package, Comsol Multiphysics, is used for the subject model. The model differs from the other 3-D models because it uses a finite element code as compared with a finite difference scheme.25,26 Finite element approaches often provide more flexibility than other numerical solving methods because they incorporate non-structured gridding, which allows for modeling complex geometries. Thus the flexibility of the modeling platform and general availability of the model to any user (because of its commercial nature) are strengths of the model described herein. Despite these strengths, it should be noted that advocates of a particular method (finite difference vs. finite element) can easily support the strengths of their preferred solution method; however, it is the authors’ opinion that the finite element model incorporated in this research is a user-friendly approach that allows various site features to be modeled with the relative ease associated with the use of commercial codes.


General Model Structure

Figure 1 illustrates the conceptual framework for the model domain. Within the domain, advective soil gas flow is induced by pressure gradients in the soil. In the implementation presented here, these pressure gradients result from a “chimney effect” slightly depressurizing the building of interest (which sits in the middle of the domain) relative to atmosphere.27 This depressurization is transmitted into the soil through the assumed breaches (cracks) in the building foundation, which also serve as the pathways for vapor intrusion into the building. The pressure differential leads to the flow of soil gas into the structure via cracks in the foundation. In addition to advection, there is an upward diffusion of the contaminant from the subsurface source toward the surface of the ground and the structure. This diffusional transport occurs because of the concentration gradient that exists between the in-ground source and the atmosphere, which presumably contains negligible amounts of the contaminant vapor. A portion of the upward-diffusing contaminant is captured by and carried with the advective flow entering the structure. The other portion of the diffusing contaminant enters the atmosphere from the open surface of the ground surrounding the structure. It should be noted that for the scenarios presented herein, a negative pressure differential was assumed to control advection into the building; however, positive pressure differentials can also occur and result in indoor air entering the subsurface.28 The model presented herein is appropriate for both positive and negative disturbance pressures.

Figure 1
Schematic of model domain.

Table 1 lists the governing equations used as the basis of the model. Figure 2 outlines the solution strategy used for solving these equations.

Figure 2
Overview of modeling process.
Table 1
Summary of model equations.

In Figure 2, steps 1–3 define the relevant geometry and geology for the model. Step 1 requires gathering relevant soil and building data and defining the scope of the case being modeled. Steps 2 and 3 involve rendering the case solvable using the finite element method by dividing the domain of interest into small volume elements, the relationships between which are defined by the working equations shown in Table 1.

As part of steps 2 and 3, variable mesh sizes are used to provide more elements in areas of interest (i.e., where steep concentration and velocity gradients are anticipated). Typical element growth rates for finite element models range from 1.1 to 2, with 1.1 considered fine and 2 considered coarse. Because vapor intrusion involves severe differences in scale with regard to the site dimensions compared with crack dimensions (as discussed below), an element growth rate near 2 was used for the research discussed herein. The overall element dimensions were constrained to a minimum size of 0.012 m and a maximum size of 3.5 m. These growth parameters are useful in generating meshes that can be solved using personal computers; however, users should carefully evaluate the solutions to ensure accurate model results. For an introduction to the finite element method, the reader is directed to Reddy.31

A typical meshed model domain for the research discussed herein is shown in Figure 3. To provide adequate mesh density in the area of steep concentration and velocity gradients, an embedded surface was created to surround the immediate vicinity of the foundation (step 2). The volume between the embedded surface and the foundation walls was meshed finely and is shown in Figure 3 as the darkened area near the center of the domain, around what would be the foundation of the structure (not explicitly shown in that figure). To ensure that the model was converging to an accurate solution, the mesh density in this volume was increased during model iterations until a stable solution was obtained, as outlined in Figure 2.

Figure 3
Cross-section A–A’ (from Figure 1) showing meshed domain.

The model solves for advective flow (solution of which is shown as step 4), using the soil gas continuity equation (eq 1) subject to a user-defined disturbance pressure at the foundation-soil interface. This disturbance pressure is applied at the soil gas entrance into a crack, and it is the influence of the internal pressure in the structure that plays an important role in determining the flow of gas in the soil around the structure. Strictly speaking, the pressure immediately outside of a foundation crack is not identical to the indoor pressure, but it is closely related through eq 2 of Table 1 (see step 6).

The result of the calculation in step 4 is a soil gas flow field and pressure field around the structure. An example of the calculated pressure field for scenario 1 (a building in the middle of an open field) is included as Figure 4. After this calculation, the user should verify that the specified boundary conditions are indeed satisfied. Normally, this means there should not be significant lateral soil gas flow into or out of the domain, because the domains are normally chosen to be large enough that there is no influence of the structure disturbance pressure on soil gas flow far from the structure. This is the basis of a common “no-flow” boundary condition assumption at the edges of the domain.

Figure 4
Color pressure profile for single building (scenario 1).

In addition, the solution should satisfy continuity (i.e., the airflow into the domain should equal the airflow out of the domain). For the cases modeled herein, the air enters the domain from the atmosphere and exits the domain through the crack, as shown in Figure 1. To evaluate whether or not continuity is achieved, weak constraints are included when solving the model. The use of weak constraints solves for one or more Lagrange multipliers, which serves to balance the finite element projection and provide accurate determinations of fluxes at boundaries.32 The additional Lagrange multiplier allows the user to verify continuity and accurately determine the soil gas flow through the crack (or characteristic entrance region [CER], as discussed below). Without weak constraints, accurate fluxes cannot be readily extracted from the modeling software because the entire boundary flux is not calculated (even when the model has converged at an accurate solution). A more detailed discussion of weak constraints is included in the text by Reddy.31

If the boundary conditions or continuity for air are not satisfied, then the chosen finite elements may be too large, the element size should be decreased, and eq 1 should be re-solved. A properly converged solution should not be sensitive to further decreases in element size. This iterative process is typically most dependent on meshing (steps 2 and 3) and secondarily on domain size (step 1). If hardware limitations do not allow adequate meshing or a large enough domain, the calculations can be reconfigured to take advantage of axes of symmetry in the domain. For the scenarios presented herein, it was not necessary to use axes of symmetry.

Step 5 refers to the modeling of the region around the crack and will be separately discussed below.

Step 6 accounts for the pressure drop that occurs as a result of the soil gas passing through the foundation crack. Equation 2 is a standard expression for flow between parallel plates and is used to calculate the pressure drop caused by gas flowing through the crack into the house. This pressure drop is the difference between the applied disturbance pressure at the foundation-soil interface and the indoor pressure. For the simulations presented herein, the pressure drop through the crack was negligible because its magnitude was less than 5% of the applied disturbance pressure for all cases. As indicated in Figure 2, if the pressure drop calculated using eq 2 was more substantial, the model can be (and should be) rerun using a corrected user-defined disturbance pressure (and iterating to a correct pressure on both sides of the crack).

Once the imposed pressure field around the house has been successfully evaluated, the results of the solved soil gas continuity equation (eq 1) are then coupled with the chemical transport equation (eq 3), which allows for both advective and diffusive transport of contaminant species in the soil (steps 7 and 8 of the flowchart in Figure 2). The input of actual contaminant properties is deferred until this point because the solution up to this step is not contaminant-specific. It is important to note that soil gas consists mainly of air, which is what allows for the decoupling of the solution of the continuity equation (eq 1) from the solution of the chemical transport equation (eq 3). This is valid whenever the contaminant species represent only a small fraction of the soil gas concentration, which is typically the case for contaminants of interest. In such a case, the advective flow is dominated by the airflow, and the small amount of contaminant is carried along by this flow. Naturally, this assumption should be checked, but it will always be true for the cases presented here.

After step 8, the user should examine the results for any instabilities, which are present as negative concentrations of contaminant. Instability is usually a consequence of steep gradients not being properly handled by too large of finite element meshes. If instabilities are found, steps 2–8 are repeated until the entire solution is stable. Finally, the successfully converged coupled eqs 1 and 3 are used to produce 3-D pressure, velocity, and chemical concentration fields surrounding the structure. With these quantities available, flow of soil gas through a crack and into the structure can be calculated. As mentioned previously, the model should be solved using Lagrange multipliers (weak constraints), so that accurate fluxes at boundaries can be determined. On the basis of the calculated gas flow and contaminant concentration at the crack, the indoor air concentration is analytically determined (eq 4) as step 9.

Although the above procedure sounds laboriously iterative, an experienced user can converge solutions with limited time investment. The most time consuming aspect is domain generation and assigning appropriate mesh dimensions and densities. Once the domain and meshing have been defined, model run time typically takes less than 1 hr on an ordinary personal computer with 4 GB of RAM. For scenarios such as those presented herein, domain and mesh generation typically require an experienced user to invest less than 1 wk of time.

Consequences of Many Different Length Scales on Modeling

The calculations above are performed for a structure with typical characteristic dimensions of tens of meters, located in a domain (i.e., on a parcel of property) characterized by dimensions of tens or hundreds of meters. It was noted above that the entry of contaminants often takes place through foundation cracks in the structure. These cracks are typically millimeters in size (i.e., many orders of magnitude smaller than the structure or domain of interest). A finite element method requires dividing the whole domain of interest into finite elements of predefined sizes, for which the calculations are individually performed. Although computational packages permit varying the size of these elements throughout the domain, selection of element sizes must be carefully made to successfully converge a solution.

It is impractical to use finite elements of submillimeter dimensions (needed to describe the details of a very small crack) when dealing with a domain that is tens or hundreds of meters in size. To do so would require many very small elements to describe the geometry of the full domain, and this in turn would be computationally intense and require long computation times. On the other hand, using too coarse of a scale (finite elements that are too large) can cause convergence difficulties (i.e., inaccurate solutions result because the elements do not capture details important to the solution).

To overcome the disparate length scale issue and still maintain accurate solutions, an approach was adopted that compromised detail around the smallest length scale, defined by the foundation cracks, such that the model was tractable using an ordinary personal computer. Because cracks are inherently variable in nature, it makes little sense to implement a high level of detail when describing them. To this end, in solving eqs 1 and 3 near the foundation, the concept of a CER was adopted. The CER is simply a plane immediately adjacent to a crack at the soil-foundation interface that leads into the house, as illustrated in Figure 1. The CER, which is wider than the crack, is only a computational convenience to reduce the number of grid points needed in the vicinity of a crack, thereby allowing more points to be used in the rest of the domain while retaining a necessary level of detail near the foundation. For the scenarios presented herein, 4 GB of RAM was adequate for developing reliably meshed models for a CER 10 cm in width around the perimeter of the foundation.

The CER could theoretically be chosen to be any width, but it is critical that the CER be similar enough in size as compared with the actual crack dimensions so that the model is representative of the crack entry scenario being modeled. The choice of CER dimensions impacts the solution process in step 5 of Figure 2. The choice of CER dimensions (and associated meshing) is determined by an iterative procedure that guarantees convergence to a true solution. After the soil gas continuity equation is solved (eq 1), the user should increase the mesh density near the CER (step 5). If the obtained velocity and concentration results change with decreasing finite element size, then the original mesh choice was not adequate and steps 2–4 should be repeated (with attention to areas with steep concentration or pressure gradients) until the solution is stable. Once this solution has converged, the user should verify that continuity is satisfied by comparing the volumetric flow entering the domain through the surface and the volumetric flow leaving the domain through the foundation crack. If not, steps 2–5 are repeated until continuity is achieved.

Boundary Conditions

The domain should be extended long enough in width and length so that the effects of any pressure field disturbances near the structure are negligible at the chosen boundaries. There is typically a no-flow boundary condition imposed at the boundary of the domain (i.e., soil gas should not be laterally flowing into or out of the domain). In this way, the contaminant concentration profile in the soil far from the structure is determined purely by upward diffusion of the contaminant, as it would be in an open field. In the cases considered here, a 100- by 100-m domain satisfied the condition that there was no significant flow induced at the boundary for the building depressurizations of typical interest.

No-flow boundary conditions (and no contaminant flux conditions) were specified at all solid foundation boundaries (except at the CER). The lower horizontal domain boundary (the water table) is the source of contaminant vapor and is subject to a no-soil gas-flow condition (meaning that no air flows through this boundary, but of course contaminant vapor is produced at the boundary). The groundwater surface serves as an infinite contaminant source. The contaminant concentration was defined as the nonpotable groundwater standard for trichloroethylene (TCE) in Rhode Island (0.54 mg/L). The pressure at top of the domain (soil surface) was defined as reference gauge pressure (p = 0), and the contaminant concentration there was defined as zero (although it could be set to any concentration to permit modeling in areas with elevated atmospheric concentrations). The pressure at the CER is set to be p = −5 Pa. This pressure can be modified to account for the pressure drop that occurs as the soil gas travels through the crack (step 6, Figure 2).

Contaminant mass flow rate across the CER (MCER) was set to be equal to the contaminant mass flow rate through the crack (Mck) as shown in eq 5. Jck was determined using eq 4 in Table 2.


The indoor air concentration is initially assumed to be negligible. This assumption is verified following solution convergence. For all of the scenarios presented herein, and most vapor intrusion cases, this assumption will be accurate. The exception will be when the indoor air concentration is greater than the concentration at the CER, such as may be the case when elevated indoor air concentrations are present relative to subsurface concentrations.

Table 2
Summary of model parameters.

It should be noted that the model presented herein is suitable for evaluating diminishing sources, periodic pressure fluctuations, and other transient and spatially variable boundary conditions. The effect of these variable boundary conditions will be the focus of future studies.

Scenarios Modeled

The model was exercised for the five different site-specific scenarios described below and illustrated in Figure 5. The model was exercised for each of these scenarios assuming a disturbance pressure of −5 Pa at the CER and a perimeter crack (5 mm wide) around the entire floor of the basement (similar pressure and crack dimensions have also been assumed by others25). In addition, a purely diffusive case (disturbance pressure of 0 Pa) was also modeled for each scenario.

Figure 5
Color soil gas contour plots (k = 10−11 m2) for modeled scenarios 1–5 (a–e, respectively).

Each of the scenarios was modeled assuming homogenous geologies using five different intrinsic soil permeability values (k = 10−10, 10−11, 10−12, 10−13, and 10−14 m2). By specifying different permeabilities, the effect of advective transport on vapor intrusion rates could be examined. Soil permeability is a function of many factors, including soil moisture, compaction, and soil type, among others. For the purposes of the research discussed herein, the soil properties assumed to be associated with each of the permeabilities, along with all other relevant chemical and physical model parameters, are included in Table 2. The porosity and moisture content values are within typical ranges as reported by Driscoll33 and consistent with parameters values used in previous 3-D modeling efforts.25 As discussed below (see Sensitivity Analysis), the effect of factors such as soil moisture appear to affect vapor intrusion rates mostly through changes in soil permeability, rather than the corresponding changes in chemical diffusivity values, which are dependent on air-filled and moisture-filled porosity values.

Scenario 1—Single Building (Figure 5a)

This scenario includes a single building (10 × 10 m) with a basement (2 m deep) located in the center of an open field.

Scenario 2—Parking Lot Surrounding Building (Figure 5b)

To simulate the effects of a parking lot around the building, an impervious surface was located around the entire building from scenario 1, extending 5 m from the building walls. Complete connection between the building walls and the impervious surface was assumed.

Scenario 3—Single Building with Detached Garage (Figure 5c)

To simulate the effects of a detached garage, an impervious square surface (5 × 5 m) was located 12.5 m center to center, or 5 m wall to wall from the scenario-1 building (on center). This configuration represents a garage on slab. It should be noted that no disturbance pressure was applied to the garage.

Scenario 4—Single Building with Gravel Sub-Base (Figure 5d)

For this scenario, 10 in. of gravel backfill material was defined immediately beneath the basement foundation for the building in scenario 1. It was assumed that backfill material had an intrinsic permeability of 10−7 m2

Scenario 5—Two Adjacent Buildings (Figure 5e)

For this scenario, two buildings with disturbance pressures of −5 Pa were located adjacent to each other. The buildings were identical to the building in scenario 1 and were located 14 m apart center to center or 4 m apart wall to wall.


Figure 5 illustrates the modeling results (for k = 10−11 m2) as soil gas concentration contour plots for the five different scenarios considered here. Shown as color-graded figures, these plots immediately provide qualitative information regarding the soil gas contaminant concentration. The plots show the expected trends—contaminant vapor concentration is highest at the surface of the groundwater source and decreases as it moves upward toward the surface of the soil. The contaminant concentrations throughout the domain depend on the boundary conditions at the groundwater source (which is the maximum concentration in the domain) and the atmosphere (which is the minimum concentration in the domain serving as a sink for the contaminant). At the furthest points from modeled structures, a smooth concentration gradient from source to sink exists, which corresponds to a purely diffusive process.

The influence of the building and other site features are also immediately apparent from Figure 5. There is an upwelling in the soil gas concentration profile beneath the structure. This characteristic of the concentration profile occurs regardless of advective flow around the structure. It is purely a consequence of “capping” caused by the building foundation and other impervious surfaces, which prevents direct diffusion of the contaminant to the sink (atmosphere). For contaminant molecules to diffuse upward directly below the structure toward the atmosphere (ignoring for the moment any possible entry into the house), they must “turn the corner” around the building foundation, requiring a longer diffusion pathway than if the structure were not there. Hence, the concentration cannot fall off as rapidly, causing the curvature of the concentration profiles near the foundation corner.

The capping effect of impervious surfaces explains why, in the absence of other factors that fundamentally change the nature of the problem, the soil gas concentrations beneath a structure will typically be much higher than soil gas sampled at the same depth some distance from the building (assuming vapor intrusion rates do not exceed chemical transport from the source). This is consistent with field observations in which subslab concentrations were detected at concentrations greater than soil gas concentrations measured exterior from the building.34

Comparing the single building scenario (Figure 5a) with the parking lot around the building scenario (Figure 5b) further highlights the importance of any impervious surface on the concentration profile. The more extensive the impermeable capping around a structure, the broader the high concentration “neck” beneath the building.

On the one hand, these concentration profiles suggest that the vapor intrusion problem is relatively insensitive to many factors such as permeability and nearby structures. All five scenarios show similar looking profiles. This is, however, why soil gas concentrations by themselves are not sufficient for characterizing vapor intrusion risks. It is not the concentration of the soil gas outside of the structure alone that determines indoor air concentration. Rather, it is the combination of soil gas concentration and flow rate of air carrying the contaminant into the structure that determines the indoor air concentration.

The quantitative results for the scenarios modeled herein are summarized in Table 3. From these results, it is clear that the indoor air concentration is a function of soil properties. As the permeability value decreased, the indoor air concentrations also decreased, regardless of the scenario. Each scenario resulted in different indoor air concentrations, but Table 3 reveals that for the scenarios modeled herein less than a factor of two separated the maximum and minimum indoor air concentrations for each permeability value.

Table 3
Summary of model results.

As shown in Table 3, indoor air concentrations vary by orders of magnitude; however, soil gas concentrations at the crack and the subslab were nearly identical for permeability values ranging from 10−11 to 10−14 m2. This is consistent with chemical transport being primarily diffusion dominated in the subsurface, which is consistent with EPA’s conceptual model for vapor intrusion.11,35 The exception to the diffusion-dominated paradigm is the highest permeability case (k = 10−10 m2), which displayed atmospheric dilution resulting in the lowest soil gas concentrations, although this permeability resulted in the highest indoor air concentration. For diffusion-dominated scenarios, soil gas concentration profiles are fairly insensitive to permeability. Concentrations are governed by air- and moisture-filled porosities, which are used to determine effective chemical diffusivities (see eq 3).

Attenuation Factors

The relationship between indoor air concentrations and subsurface soil gas concentrations is commonly expressed using attenuation factors (also known as alpha values). In the equations shown below, attenuation factors are typically less than 1, although in some cases attenuation factors can be greater than 1, which would indicate subsurface concentrations are lower than indoor air concentrations, and indoor or atmospheric sources may be responsible for the measured indoor air contamination. In general, lower attenuation factors indicate more attenuation.

αi=Indoor Air ConcentrationSoil Gas Concentrationcollected at location i

αgw=Indoor Air ConcentrationKH·Groundwater Concentration

Attenuation factors calculated using EPA’s target indoor air concentrations and site-specific soil gas data are often used for screening purposes. EPA’s draft vapor intrusion guidance11 recommends values of 0.001 and 0.1 for αgw and αsubslab, respectively, for generic screening values. These screening values are calculated using EPA’s target indoor air concentration in the numerator of eqs 6 and 7. When calculated in this way, attenuation factors less than the recommended values may warrant additional evaluation of vapor intrusion risks. When the numerator and denominators in eqs 6 and 7 are defined using site-specific data, the attenuation factor indicates how much attenuation is occurring between the subsurface and the indoor air. Again, the more attenuation, the lower the attenuation factor.

Figure 6 presents the calculated attenuation factors for the scenarios modeled herein. The αgw and αsubslab values are nearly identical. On the basis of the effect of impervious surfaces on the soil gas concentration profiles (discussed previously), this relationship is expected. The subslab concentrations are representative of the soil gas concentrations immediately beneath the foundation, but these concentrations are similar in magnitude to the source concentration at the groundwater table (see Table 3). Consequently, the groundwater and subslab attenuation factors are also similar.

Figure 6
Comparison of attenuation factors and mass flow rates through crack for TCE.

This similarity of groundwater and subslab attenuation factors is in contrast to the 2 orders of magnitude different attenuation factors EPA recommends for screening purposes11 (0.001 and 0.1 for groundwater and subslab, respectively). In addition, Dawson36 reported differences in αgw and αsubslab calculated using data from the EPA vapor intrusion database; however, the attenuation factors calculated from this dataset showed considerable scatter. Although the reason for this scatter in the empirical datasets is not well understood, it should be noted that only limited information is available for the data within the database; therefore analysis is fairly uncontrolled. Dawson36 suggests that spatial and temporal variability may be responsible for the variations observed in the individual attenuation factors, as well as nonrepresentative samples. Additional controlled studies are needed to fully evaluate attenuation values for real field sites.

As shown in Figure 6, qualitatively similar trends are present for the calculated attenuation factors and the mass flow (intrusion) rates of the contaminant into the building. Both metrics increase as indoor air concentrations increase. From eq 4, it is apparent that the indoor air concentration is primarily a function of the mass flow rate and the building exhaust rate (because Qck is typically small compared with AeVb). Therefore, for a given building exhaust rate, the indoor air concentration is solely dependent on the mass flow rate of the contaminant into the building. Although mass flow rate is dependent on soil gas concentrations, it is also greatly affected by soil gas flow rate (which is a function of soil type and building features). Although it is difficult to measure without the assistance of a vapor intrusion model (such as the tool presented herein), mass flow rate appears to be a good indicator of vapor intrusion potentials and is qualitatively similar to the commonly used attenuation factor.

Figure 7 provides additional evidence for why soil gas concentrations (specifically subslab concentrations) by themselves are not good indicators of vapor intrusion risks. In general, the modeled subslab concentrations are approximately 75% of the source concentration (at 8 m below ground surface [bgs]) regardless of scenario and/or permeability. Figure 7 depicts a relatively flat slope for subslab concentrations. The specific value for each subslab concentration with respect to scenario/permeability is shown in Table 3. Overall the subslab concentrations vary less than a factor of 2. However, these subslab concentrations are associated with nearly 2 orders of magnitude difference in indoor air concentrations. This observation is potentially important because regulations often rely on subslab concentrations during vapor intrusion investigations.34,35 For the cases modeled herein, subslab soil gas concentrations alone are not sufficient for characterizing health risks related to vapor intrusion.

Figure 7
Comparison of TCE subslab concentrations and indoor air concentrations.

As discussed previously, other 3-D models have been used to predict vapor intrusion potentials for buildings located in open fields (scenario 1).25 The results of the research presented herein are in general agreement as compared with other published results. For instance, the soil gas flow rate calculated using the model described herein is approximately 1.15 times greater than the soil gas flow rates presented by Abreu and Johnson.25 In addition, indoor air concentrations for scenario 1 for specific permeabilities were also in good agreement.

Spatial Variation in Soil Gas Profiles

Spatial variation in soil gas concentrations in the vicinity of the building of interest has been reported by many practitioners.3740 There are several factors that may be responsible for the spatial variability. Two possible factors investigated as part of this research are atmospheric dilution caused by building depressurization and the effect of impervious structures on chemical transport. These two factors are theoretically a function of soil gas flow rate (eq 1) and diffusivity (eq 3), respectively.

Effect of Soil Permeability

Figure 8 illustrates the soil gas concentrations that result for five different permeability values for scenario 1 (a single building in an open lot). As shown in Table 3, these different permeability values result in soil gas flow rates ranging from 7.9 × 10−4 m3/sec to 7.9 × 10−8 m3/sec. The effect of these flow rates on soil gas concentrations is minimal, except for the highest flow rate case, 7.9 × 10−4 m3/sec, which corresponds to a gravelly soil (k = 10−10 m2). For this exceptional case, the soil gas flow rate is so large that atmospheric gas is pulled into the subsurface and ultimately into the building. Although the soil gas concentration is diluted by this process, the indoor air concentration is greatest because of the combined effect of soil gas concentration and soil gas flow rates, which result in the highest mass contaminant flow rate of all of the permeabilities (Table 3).

Figure 8
Variability in TCE soil gas concentrations for scenario 1 (single building) at depths of (a) 2 m, (b) 4 m, and (c) 6 m.

Figure 8 also indicates that concentrations are generally greatest less than 5 m away from the building, but then the concentrations become constant (except the most permeable soil—k = 10−10 m2) as the distance is increased. This spatial trend in soil gas concentrations is related to the building structure itself, which serves as a cap and causes elevated soil gas concentration profiles near the structure. With lateral distance, as shown in Figure 5, these concentration profiles become constant.

It should be noted that for the cases modeled herein in which soil permeabilities were less than approximately k = 10−11 m2, the soil gas concentration profile is virtually identical whether a crack (or CER) is or is not present. This reinforces that the soil gas concentration profiles depicted in Figure 5 and Figure 8 are not affected substantially by advective transport through the crack. It is, however, important to emphasize that although a diffusive process largely determines the contaminant concentration profile in the soil, advection into the structure plays an important role in determining the indoor air contaminant concentration. Although diffusive processes dominate the global contaminant profile in the subsurface, advection and diffusion determine contaminant entry through the crack into the structure. The above discussion highlights important considerations for understanding the relationships between soil gas flow rates, soil gas concentrations, and indoor air concentrations.

Effect of Site Features

Figure 9 illustrates the soil gas concentration profiles for each of the scenarios modeled for permeability k = 10−11 m2. This single permeability was selected as an example case, because it has similar soil gas concentration profiles to the other permeabilities (except k = 10−10 m2) modeled herein.

Figure 9
Spatial variability in TCE soil gas concentration profiles for scenarios 1–5 (a–e, respectively; k = 1 × 10−11 m2).

Site-specific features included in the modeled scenarios affected the soil gas concentration profiles. For instance, elevated soil gas concentrations are present beneath impervious structures, such as the parking lot (Figure 9b) and adjacent buildings (Figure 9, c and e). The presence of a porous sub-base beneath a building appears to have a slight effect on the soil gas concentration profiles because of atmospheric dilution in the vicinity of the building. As shown in Figure 9d, the soil gas concentration near the domain boundary (35 m) is slightly higher than the concentrations 2.5, 5, and 7.5 m away from the building for depths of 2, 4, and 6 m. Although subtle, the effect of atmospheric dilution in the vicinity of the building is also present in Figure 5d.

Overall, variation in soil gas concentrations are most substantial for “shallow” depths (<4 m bgs) at locations close to the building. This is most apparent in Figure 9b, which has elevated soil gas concentrations beneath the parking lot surrounding the building. For instance, the soil gas concentration at the building wall 2 m bgs is nearly 3 times higher than the soil gas concentration detected 7.5 m from the building at the same depth. Figure 9b shows that the difference between the soil gas concentrations converge as the depth is increased.

For characterization purposes, a direct relationship between soil gas concentrations and indoor air concentrations is not obvious. Deep soil gas concentrations are nearly identical for all scenarios (Figure 9) and permeabilities (Figure 8), although indoor air concentrations vary substantially. The spatial variability of shallow soil gas concentrations does not allow direct comparison to indoor air concentrations. As mentioned previously, it is the combined effect of soil gas concentrations and soil gas flow rates that determines the indoor air concentration (eq 4); consequently, drawing conclusions about indoor air quality from only one of these parameters is imprudent and will likely result in inaccurate risk assessments.

Sensitivity Analysis

Because diffusivity and permeability are functions of soil type (see Table 2), a sensitivity analysis was performed to demonstrate the degree of importance for soil diffusivity versus soil permeability on the resulting indoor air concentrations. First, constant soil permeability (k = 10−11 m2) was modeled while changing diffusivity values that were expected based on variations in soil porosity due to increases in moisture content. Effective diffusivity values (see eq 3, Table 1) were calculated assuming the soil had moisture contents of 10, 20, 50, and 90%. Then, a constant diffusivity value (8.68 × 10−7 m2/sec) was modeled and permeabilities were varied from 10−10 to 10−14 m2. The results are shown in Figure 10 and illustrate that indoor air concentrations are most sensitive to changes in soil permeability (Figure 10b), not contaminant diffusivity (assuming unsaturated conditions) (Figure 10a). Therefore, changes in soil permeability due to increased moisture content, compaction, etc., are far more likely to affect vapor intrusion rates because of changes in permeability, which affect soil gas flow rates, as opposed to changes in diffusivity caused by corresponding changes in soil porosity. In addition, the figures show that the model results are fairly insensitive with respect to the selection of diffusivity values (on the basis of corresponding residual water contents). Therefore, considering how soil moisture variations affect permeability appears more important for vapor intrusion than accounting for changes in effective diffusivities as a result of changes in air- and moisture-filled porosities.

Figure 10
Sensitivity analysis for (a) constant intrinsic permeability and (b) constant effective diffusivity.


Vapor intrusion risks are a result of vapor-phase contaminant transport into a building. When characterizing vapor intrusion risks, special consideration should be given not only to measuring chemical concentrations, but also to understanding how those subsurface concentrations may or may not result in indoor air pollution on the basis of the potential for soil gas flow to occur. Higher soil gas concentrations do not necessarily imply higher rates of gas transport, and it is the gas transport rate that is the key consideration. The results presented herein show that for the scenarios modeled, soil gas concentrations universally increased when approaching the source and varied only slightly with soil permeability values over 4 orders of magnitude. In contrast to these similarities in soil gas contaminant concentration profiles, indoor air concentrations varied substantially, indicating soil gas concentrations by themselves do not serve as good indicators of vapor intrusion risks.

The presence of site-specific features such as sub-base material beneath a foundation or an impervious surface surrounding a building can also affect soil gas concentration profiles. For certain scenarios, increases in soil gas flow rates and decreases in soil gas concentrations at the foundation crack are observed (scenario 3—porous sub-base material), but for other scenarios, soil gas concentrations increase near the building (scenario 4—parking lot surrounding building). For cases such as these, in which the effect of site features is not easy to predict, modeling is a valuable tool to guide decision-making in risk assessment and site investigation. The 3-D model described herein allows for site-specific features to be considered and is a useful tool for better understanding how vapor intrusion risks vary from site to site and building to building. The step-by-step modeling approach outlined herein can be used to investigate various vapor intrusion scenarios to better understand and characterize vapor intrusion risks.


Among the emerging issues of national concern are the potential health risks associated with vapor intrusion into structures that might exist or be built on or near sites contaminated with volatile organic compounds or semi-volatile organic compounds. Chemical vapors released from soil and/or groundwater contamination can migrate into buildings and pose a long-term hazard to human health. The issue of vapor intrusion has recently received considerable attention in the environmental community. Unfortunately, tools to assist vapor intrusion site characterization efforts are thus far limited in number. The model described herein is proposed as such a tool, with the belief that it can help guide site investigations and improve vapor intrusion risk predictions.


The project described was supported by grant no. P42ES013660 from the National Institute of Environmental Health Sciences. The content is solely the responsibility of the authors and does not necessarily represent the official views of the National Institute of Environmental Health Sciences or the National Institutes of Health. The authors gratefully acknowledge the thoughtful input provided by Professor Merwin Sibulkin of Brown University, who passed away before the publication of this manuscript.


1. Donohue MK. Special Letter from the Managing Editor: Report in the DOE/PERF Vapor Pathway Workshop, January 27–28, 2004 at Brea, CA. Environ. Foren. 2004;5:109–122.
2. Nazaroff WW, Lewis SR, Doyle SM, Moed BA, Nero AV. Experiments on Pollutant Transport from Soil into Residential Basements by Pressure-Driven Airflow. Environ. Sci. Technol. 1987;21:459–466. [PubMed]
3. Nazaroff WW. Predicting the Rate of 222Rn Entry from Soil into the Basement of a Dwelling Due to Pressure-Driven Air Flow. Radiation Protect. Dosimetry. 1988;22:199–202.
4. Garbesi K, Sextro RG. Modeling and Field Evidence of Pressure-Driven Entry of Soil Gas into a House through Permeable Below-Grade Walls. Environ. Sci. Technol. 1989;23:1481–1487.
5. Loureiro CO, Abriola LM, Martin JE, Sextro RG. Three-Dimensional Simulation of Radon Transport into Houses with Basements under Constant Negative Pressure. Environ. Sci. Technol. 1990;24:1338–1348.
6. Revzan KL, Fisk WJ, Gadgil AJ. Modeling Radon Entry into Houses with Basements: Model Description and Verification. Indoor Air. 1991;2:173–189.
7. Wang F, Ward IC. The Development of a Radon Entry Model for a House with a Cellar. Build. Environ. 2000;35:615–631.
8. Wang F, Ward IC. Radon Entry, Migration and Reduction in Houses with Cellars. Building Environ. 2002;37:1153–1165.
9. Jury WA, Russo D, Streile G, Abd HE. Evaluation of Volatilization by Organic Chemicals Residing below the Soil Surface. Water Resource Res. 1990;26:13–20.
10. Johnson PC, Ettinger RA. Heuristic Model for Predicting the Intrusion Rate of Contaminant Vapors into Buildings. Environ. Sci. Technol. 1991;25:1445–1452.
11. OSWER Draft Guidance for Evaluating the Vapor Intrusion to Indoor Air Pathway from Groundwater and Soils. Washington, DC: U.S. Environmental Protection Agency; 2002. EPA530-D-02-004.
12. Johnson PC. Identification of Critical Parameters for the Johnson and Ettinger (1991) Vapor Intrusion Model. Washington, DC: American Petroleum Institute (API) No. 17; API; 2002.
13. Johnson PC. Identification of Application-Specific Critical Inputs for the 1991 Johnson and Ettinger Vapor Intrusion Algorithm. Ground Water Monitor. Remed. 2005;25:63–78.
14. Hers I, Zapf-Gilje R, Johnson PC, Li L. Evaluation of the Johnson and Ettinger Model for the Prediction of Indoor Air Quality. Ground Water Monitor. Remed. 2003;23:62–76.
15. Devaull GE. Indoor Vapor Intrusion with Oxygen-Limited Biodegradation for a Subsurface Gasoline Source. Environ. Sci. Technol. 2007;41:3241–3248. [PubMed]
16. Parker JC. Modeling Volatile Chemical Transport, Biodecay, and Emission to Indoor Air. Ground Water Monitor. Remed. 2003;23:107–120.
17. Turczynowicz L, Robinson N. A Model to Derive Soil Criteria for Benzene Migration from Soil to Dwelling Interior in Homes with Crawl Spaces. Human Ecolog. Risk Assess. 2001;7:387–415.
18. Mills WB, Liu S, Rigby MC, Brenner D. Time-Variable Simulation of Soil Vapor Intrusion into a Building with a Combined Crawl Space and Basement. Environ. Sci. Technol. 2007;41:4993–5001. [PubMed]
19. Krylov VV, Ferguson CC. Contamination of Indoor Air by Toxic Soil Vapors: the Effects of Subfloor Ventilation and Other Protective Measures. Build. Environ. 1998;33:331–347.
20. Turczynowicz L, Robinson N. Exposure Assessment Modeling for Volatiles—towards an Australian Indoor Vapor Intrusion Model. J. Toxicol. Environ. Health A. 2007;70:1619–1634. [PubMed]
21. Tillman FD, Jr, Weaver JW. Review of Recent Research on Vapor Intrusion. Athens, GA: U.S. Environmental Protection Agency; 2005. EPA/600/R-05/106.
22. Hers I, Zapf-Gilje R, Evans D, Li L. Comparison, Validation and Use of Models for Predicting Indoor Air Quality from Soil and Groundwater Contamination. Soil Sediment Contam. 2002;11:491–527.
23. Robinson NI, Turczynowicz L. One and Three-Dimensional Soil Transportation Models for Volatiles Migrating from Soils to House Interiors. Transport. Porous Media. 2005;59:301–323.
24. Lowell PS, Eklund B. VOC Emission Fluxes as a Function of Lateral Distance from the Source. Environ. Progress. 2004;23:52–58.
25. Abreu LD, Johnson PC. Effect of Vapor Source—Building Separation and Building Construction on Soil Vapor Intrusion as Studied with a Three-Dimensional Numerical Model. Environ. Sci. Technol. 2005;39:4550–4561. [PubMed]
26. Abreu LD, Johnson PC. Simulating the Effect of Aerobic Biodegradation on Soil Vapor Intrusion into Buildings: Influence of Degradation Rate, Source Concentrations. Environ. Sci. Technol. 2006;40:2304–2315. [PubMed]
27. IAQ Reference Guide: Tools for Schools. 3rd ed. Washington, DC: U.S. Environmental Protection Agency; 2005. EPA 402-K-95-001.
28. McHugh TE, Deblanc PC, Pokluda RJ. Indoor Air as a Source of VOC Contamination in Shallow Soils below Buildings. J. Soil Sediment Contam. 2006;15:103–122.
29. Massmann JW. Applying Groundwater Flow Models in Vapor Extraction System Design. J. Environ. Engineer. 1989;115:129–149.
30. Millington RJ. Gas Diffusion in Porous Media. Science. 1959;130:100–102. [PubMed]
31. Reddy JN. An Introduction to the Finite Element Method. 2nd ed. New York: McGraw Hill; 1993.
32. Multiphysics Modeling Guide, Version 3.4. Burlington, MA: Comsol; 2007. p. 300.
33. Driscoll FG. Groundwater and Wells. 2nd ed. St. Paul, MN: Johnson Filtration Systems; 1986. p. 67.
34. Final Guidance for Evaluating Soil Vapor Intrusion in the State of New York. Albany, NY: New York State Department of Health; 2006.
35. Brownfields Technology Primer: Vapor Intrusion Considerations for Redevelopment. Washington, DC: U.S. Environmental Protection Agency; 2008. EPA 542-R-08-001.
36. Dawson H. EPA Vapor Intrusion Database: Preliminary Analysis of Attenuation Factors. Proceedings of the 16th Annual AEHS West Coast Conference on Soils, Sediment and Water; Amherst, MA. Association for Environmental Health and Sciences; 2008.
37. Dawson H. Revised Guidance on Exterior Site Specific Screening. Proceedings of the 16th Annual AEHS West Coast Conference on Soils, Sediment and Water; Amherst, MA. Association for Environmental Health and Sciences; 2006.
38. DiGiulio D. Assessment of Vapor Intrusion in Homes near Raymark Superfund Site—Discussion of Findings Related to EPA’s Guidance. Proceedings of the 16th Annual AEHS West Coast Conference on Soils, Sediment and Water; Amherst, MA. Association for Environmental Health and Sciences; 2006.
39. Wertz W, Anders K. Near-Building and Sub-Slab Sampling at the Endicott (NY) Site, Implications for Site Screening Approaches. Proceedings of the 16th Annual AEHS West Coast Conference on Soils, Sediment and Water; Amherst, MA. Association for Environmental Health and Sciences; 2006.
40. Eklund BM, Simon MA. Concentration of Tetrachloroethylene in Indoor Air at a Former Dry Cleaner Facility as a Function of Subsurface Contamination: a Case Study. J. Air & Waste Manage. Assoc. 2007;57:753–760. [PubMed]