Objectives

The purpose of this thesis is to develop a two-component model to predict the mole fraction distributions of VOC and the migration of evaporation front of NAPL phase after the removal of the storage tank in a homogeneous unsaturated soil. To the best of our

knowledge, the existing models in simulating the natural evaporation of VOC in unsaturated soil have the problem domains all with fixed boundaries. In contrast, the present model considers a moving boundary to describe the downward move of evaporation front of NAPL phase. The front, which migrates downward with time, is treated as a lower boundary while the land surface is used as the upper boundary for the region above the front. For the region below the front, the front becomes the upper boundary and the lower boundary is located at some distance below the land surface. Based on the present model with these boundaries, VOC mole fraction distributions and the front location are solved by the finite difference method with a moving grid approach. This numerical model can predict the mole fraction distributions between two components and the movement of the front, analyze evaporation time of VOC, and assess the influences of initial mole fraction, soil porosity as well as chemical volatility on VOC migration in the unsaturated soil. Moreover, a

single-component model is obtained from simplifying the two-component model and solved analytically using Boltzmann’s transformation. Then, an analytical expression for the moving speed of the front can be developed from the solution and used to assess the time of

vanish of NAPL at a specific location below the land surface.

CHAPTER 2 METHMETICL MODEL

This chapter presents the mathematical models and related solution procedures for VOC transport in homogeneous and unsaturated soils.

2.1. Mathematical model: Two-component case

Figure 1(a) shows a leaking storage tank located on the top of land surface and filled with VOC. Consider that the VOC leaking from the tank has four different phases (namely, gas, liquid, adsorbed and residual NAPL phases) and distributes uniformly in the unsaturated zone. Each phase has a saturated or equilibrium concentration and constant volumetric content. The VOC has an initial NAPL saturation S0 in the soil. The saturation of each phase denotes as the volume percentage in the soil pore and the sum of saturation of each phase equals one. The evaporation front of the NAPL, denoted as s(t), initially stays right at the land surface, i.e., z = s(t) = 0 where z is the vertical axis, and moves downward with increasing time. Figure 1(b) shows the contamination scenario in which the gas phase VOC begins to evaporate to atmosphere and the NAPL starts to vaporize to gas phase after the tank being removed. Assume that the NAPL phase of VOC evaporates fully above the front and the front migrates instantaneously when the evaporation occurs. In other words, the VOC presents only in gas, liquid, and adsorbed phases and the NAPL saturation, SR, equals zero between the land surface and the front. Below the front, the residual VOC still exists. At the front, volatilization occurs from NAPL phase to gas phase instantaneously and the content

of NAPL transferring to gas phase follows the conservation of mass.

The total concentration of VOC is the sum of the concentration of each phase; that is:

R R b S L L G G

T C C C

C = θ + θ + ρ +ρ θ (1)

where CG, CL, and Cs represent the chemical concentrations in the gas, liquid, and adsorbed phases, respectively, and ρb and ρR are the soil bulk density and density of NAPL,

respectively. The symbols θG, θL, and θR are the volumetric contents of gas, liquid, and NAPL phases, respectively, and denote as

G

G φS

θ = , θ_L =φS_L, and θ_R =φS_R (2)

where φ is the soil porosity and SG and SL are gas and liquid saturations, respectively.

The equilibrium relationships between the gas and liquid phases as well as the liquid and adsorbed phases may be expressed, respectively, as

L H

G K C

C = and C_S =K_DC_L =K_ocf_ocC_L (3)

where KH is Henry’s Law constant, KD is the distribution coefficient, Koc is the organic carbon partition coefficient, and foc is the soil organic carbon fraction. The equilibrium relationships given in Equation (3) are linear and reversible and their coefficients are dependent on

different chemical characteristics and soil properties.

The equation of mass conservation for those four phases in the unsaturated soil is (Corapcioglu and Baehr, 1987)

[

^{+ + + +∇ =0} tortuosity factor accounting for the reduced flow area and increased path length of diffusing molecules in soil can be related to the soil porosity and the fluid content in the soil.

Millongton and Quirk (1961) defined that

⎟⎟

is the air-gas diffusion coefficient.

Consider that VOC has multiple components. Substituting Equations (2) and (3) into Equation (4), the mass-conservation equation in three dimension for those four phases in the unsaturated soil becomes (Zaidel and Zazovsky, 1999)

( )

where and are equilibrium gas and liquid concentrations of pure component, respectively, f is the mass fraction of organic compound in NAPL, i is the number of each component, and N

p

CG C_L^p

C is the total number of organic compounds. The saturated gas phase

concentration for pure component can be estimated from the idea gas law as

T M P

C_G^P = ⁰ /ℜ (9)

where P^P0 is the saturated vapor pressure of the VOC, M is the molecular weight of the VOC, is the ideal gas constant, and T is the absolute temperature.

ℜ

To describe the behavior of volatilization of VOC from NAPL phase to gas phase, the evaporation front is introduced as a moving boundary in the unsaturated soil. The problem domain for VOC transport with a moving boundary in the soil can be divided into three regions. The VOC transport equations in these three regions are described below.

2.1.1 Below the evaporation front

In this region, NAPL phase exists and the saturation of NAPL phase is greater than zero, i.e., S_R >0. Equation (8) can then be expressed as

with representing mole fractions of organic compounds in the NAPL phase, i.e., . For a one-dimensional system, Equation (6) becomes:

ui based on Equation (11), the mass-conservation equation for each component can be written

P Gi R i ρ φ/C σ =

as: NAPL saturation S0 and mole fraction u0 of component one as:

)

Substituting Equations (10) and (14) into Equation (12) yields

( )

2.1.2 Above the evaporation front

In this region the NAPL phase fully evaporates. i.e., S_R =0. For a two-component VOC, Equation (11) can be written for each component as follows:

2

where u1 and u2 represent normalized concentrations above the evaporation front,

2.1.3 At the evaporation front

Assume that the evaporation of NAPL phases occurs instantaneously and the sum of mole fraction equals one, i.e., u1+u2 = 1. Combining Equations (12) and (13), the equation for mass conservation at the front can be expressed as

t

Taking the limits of and , Equation (18) describing the front can then be written as

where the superscripts t+ and t– of the mole fraction in Equation (20) denote the mole fraction at the time slightly after and before the volatilization, respectively, and z+ and z– represent the mole fraction at the region slightly below and above the front, respectively. Equation (20) can be further simplified as

(

t t

) (

t t

)

R

2.1.4 Boundary and initial conditions

For simulating the ambient environment near the surface, the gas phase VOC is

considered to diffuse across a stagnant air boundary layer with thickness d and the gas phase concentration is assumed equal to zero at the top of boundary layer. The flux diffusing to the atmosphere can then be expressed as (Jury et al., 1983)

G T G G

R hC z

D C =

∂

∂ , at z=0 (22)

where . If the thickness of air boundary layer near the surface is very small and negligible, one may assume that d = 0. Equation (22) can then be expressed as

d D h= _G^air /

=0

CT , at z=0 (23)

For a multicomponent VOC, the upper boundary conditions can be written from Equation (23)

as

( )

^0, 2

( )

^{0, 0}

1 t =u t =

u (24)

The conditions of the mole fractions at the lower boundary are denoted as:

( )

L t u^o

u₁ , = , u₂

( )

L,t = 1−u^o (25) where L is the depth of the lower boundary.

The initial NAPL saturation and the mole fractions of each component are assumed

spatially uniform; they are

( )

z,0 S₀

S_R = (26)

( )

^{z uo}

u₁ ,0 = , ^u₂

( )

^z,0 ^=1−uo (27)

The evaporation front is initially located at the land surface and thus denoted as

( )

⁰ =⁰

s (28)

2.2 The numerical method in solving the model

This section presents the numerical method used to solve the two-component model.

2.2.1 Finite difference approximation

The equations of describing VOC transport in the three regions are solved separately by the finite difference method. An interpolative moving grid approach (Javierre et al. 2006) is adopted to handle the moving boundary problem. The total number of nodes within the problem domain is equal to N and r is the nodal number assigned at the evaporation front beginning from the land surface. Accordingly, the number of grids from the land surface to the front is r-1, the number of grids below the front is N-r, and the initial grid size dz is equal to L/(N-1). The grid sizes above and below the front defined as dzr and dzN-r, respectively, need to be re-estimated after each move of the front. To avoid introducing large truncation error, the grid sizes dzr and dzN-r should be close to dz. If the front move to a location between the nodal numbers initially assigned as j and j+1, the new grid sizes of dzr and dzN-r

are then estimated by s/(r-1) and (L-s)/(N-r), respectively, where r = j+1. The backward difference relative to time and central difference relative to space are used to approximate the VOC transport equations. Therefore, the difference equation for the mole fraction in the

region below the front obtained from Equation (16) is

where n is the number of time step, and dt is the time interval. The difference equations for the mole fractions of the two components in the region above the front obtained from Equations (17) and (18) are, respectively:

2

Finally, the difference equations for the mole fractions of the two components at the front also derived from Equations (17) and (18) with different grid sizes above and below the front are, respectively, expressed as

r

Substituting Equation (10) into Equation (21), the difference equation for the front location can then be obtained as

( ) (

n

)

R

2.2.2 The solution procedure of the model

Because the location of the evaporation front s(t) is unknown, a trial and error procedure is therefore taken to find the front location. The steps of the procedure are listed below and the related flowchart is shown in Figure 2:

1. Give the initial front location (Equation (28)) and nodal values of mole fraction based on the boundary conditions (Equations (24) and (25)) and initial conditions (Equations (26) and (27)).

2. Guess front location after the start of evaporation.

3. Determine the nodal number for the front and the grid sizes based on the front location.

4. Solve the solutions for the region below the front (Equation (29)), for the region above the front (Equations (30) and (31)), and at the front (Equations (32) and (33)).

5. Compute the front location based on Equation (34).

6. Proceed to next time step if the location of the front converges, i.e., the difference between two succeeding estimations of the front location is less than a very small value, e.g., 10^-10 m;

otherwise, go back to step 3.

2.3 Simplified model: Single-component case

This section presents a single-component model with the governing equation simplified from the two-component model. The domain of the single-component model is also divided into three regions based on the front location. In a one-dimensional homogeneous and unsaturated soil system, Equation (6) representing the VOC transport can be written as

2 0 2.3.1 Below the evaporation front

In this region, the NAPL has not vaporized to gas phase yet. The VOC concentrations in each phase are the initial saturated concentrations and the total concentration can be expressed as

is the saturated VOC concentrations for pure component in adsorbed phases, CT0

is the initial total concentration, and θG0, θL0, and θR0

are the initial volumetric contents, respectively.

2.3.2 Above the evaporation front

In this region, NAPL phase of the VOC completely vaporizes to gas phase; therefore, the VOC presents only in gas, liquid, and adsorbed phases, i.e., CT = CGθG+ CLθL+C^Sρb. Jury et al. (1983) used a ratio R to represent each phase in relation to the total concentration.

Accordingly, one may introduce the following equation based on Equations (1) and (3):

With Equation (37), Equation (35) can be written as

2 0 where denotes as effective diffusion coefficient. Equation (38) describes the VOC transport in gas and liquid phases between the land surface and evaporation front. In reality, the volatilization of NAPL occurs right at the evaporation front which will be discussed in the next section.

G G

E D R

D = /

2.3.3 At the evaporation front

The evaporation front moves downward while the NAPL vaporizes to the gas phase.

Assume that the transformation of the contents between these two phases occurs instantaneously and the liquid volumetric content is unchanged, i.e., θL = θL0

. In addition, the VOC concentrations in each phase are still saturated,i.e., CG = CGP

, CL = CLP

, and CS = CSP

. The total concentration at the front can therefore be expressed as

0 T

T C

C = (39)

The gas, liquid, and NAPL phases remain in the unsaturated soil pores, the sum of θG, θL, and

θR equals soil porosity φ, i.e., θ_G+θ_L+θ_R =φ. Assumes φ does not change with time,

the transformation of volumetric content with respect to time among each phase is conserved.

Thus,

following relationship from Equation (1)

( )

the front can then be written as

( )

where the superscripts + and - denote the volumetric content at the time slightly after and before the volatilization, respectively, and the concentrations at the region slightly below and above the front, respectively. Consider that the volatilization occurs instantaneously, therefore equals and equals zero after evaporation. The VOC

concentrations in each phase are the initial saturated concentrations below the front and the concentration gradient of gas phase below the front is naturally equal to zero, i.e.,

+

θR θ_R⁰ θ_R⁻

=0

∂

∂C_G⁺ z . Since the liquid density is higher than the gas phase concentration about three

orders (Falta et al., 1989, Tables 1 and 4), the term related to CGP

on the left-hand side of Equation (42) is thus negligible. Accordingly, Equation (43) representing the front z = s(t) can be expressed as

z

2.3.4 The analytical solution of single-component model

Consider that the VOC is saturated or in equilibrium state in different phases and uniformly distributed in the unsaturated soil initially. The mathematical model describing the single-component VOC transport in the soil consists of Equation (38) as the governing equation, Equations (39) and (43) as the lower boundary conditions, and Equation (23) as the upper boundary condition.

Based on Boltzmann’s transformation, a new variable is defined as ξ = z 2 D_Et . Equation (38) can then be transformed to an ordinary differential equation as

0

The solution of Equation (44) can be obtained as (Carslaw and Jaeger, 1959)

( )

B erf

A

C_T(ξ)= ⋅ ξ + (45)

where erf(ξ) is the error function and A and B are unknown coefficients. Substituting Equation (45) into Equation (23), the result for the concentration distribution is

( )

_⎟^⎟

Substituting Equation (39) into Equation (46), the evaporation front s(t) and coefficient A can then be obtained, respectively, as

t

where α is an unknown constant depending upon the soil parameters and contaminant characteristics.

The time of vanish of NAPL can be solved by Newton’s method (Yeh, 1987) from Equation (47) when the front reaches a target location below the land surface designated by the environmental or legal requirement. In addition, the moving speed of the evaporation front can also be obtained after taking the derivative of Equation (47) with respect to time and the result is

U_f t 2

= α (49)

Substituting Equations (46), (47) and (48) into Equation (43) yields

⎟⎟

The unknown constant α can then be easily determined from Equation (50) by Newton’s method. Note that the normalized total concentration is defined as CT(z,t)/CT0

, representing the mole fraction in the single-component model.

CHAPTER 3 RESULTS AND DISCUSSION

Leaks of petroleum fuels from the underground storage tanks are common problems for soil contamination. The petroleum spills are often associated with aromatic hydrocarbons such as benzene, toluene, ethyl benzene, and various xylene isomers (BTEX). In this section, the hydrocarbons of benzene and toluene are chosen to simulate their transport and mole fraction distributions in unsaturated soils using the two-component model. In the past, Carbon tetrachloride was commonly used as coolant in industry or produced as the fire

extinguishers. Carbon tetrachloride is highly toxic; a small amount of this chemical residing in the soil may pose severe environmental problems. The carbon tetrachloride in the

unsaturated soil is considered as a target VOC and analyzed using the single-component model.

Six cases are considered to address the issues in regard to the evaporation rate,

evaporation front movement, mole fraction, and concentration distributions of VOC for the present models. Case 1 is to compare the mole fractions of toluene at various evaporation times predicted by single-component and two-component models. Case 2 examines the effect of initial mole fraction on the evaporation and the changes of the mole fraction distributions of benzene and toluene. Case 3 investigates the effect of soil porosity on vaporization of carbon tetrachloride from NAPL phase to gas phase while case 4 addresses the issue of carbon tetrachloride transport based on the model with and without considering

the presence of NAPL phase. Case 5 studies the migrations of evaporation front for different contaminants, namely carbon tetrachloride and toluene. Case 6 discusses the effect of

effective diffusion coefficient of carbon tetrachloride on the moving speed of evaporation front. The values of the soil chemical properties are listed in Table 1 and the properties of benzene, toluene, and carbon tetrachloride are given in Table 2 for these six cases. Note that the depth of the lower boundary L is chosen as 5 m, the total number of nodes N is 10000, and the time interval dt is 0.1 sec in the case study when adopting the finite different

approximation for the two-component model.

3.1 Case 1: Different evaporation times in two models

This case uses the same assumptions for both single-component and two-component models and considers that toluene is the only VOC found in the soil, i.e., u0 = 1. Figure 3 shows the mole fraction distributions of toluene versus depth predicted by the

single-component and two-component models at various evaporation times. The dashed line denotes the solution of single-component model while the solid line represents the results predicted by the two-component model. Moreover, the symbols of rhombus, triangle, and circle represent the mole fractions at times 1, 10, and 100 day, respectively. This figure shows the front locations at various evaporation times and at the front the mole fraction equals its initial value for the single-component VOC. The figure also shows that the curves

predicted by both models are fairly close, implying that the results predicted by the

two-component model with the present numerical approach match well with those estimated based on the analytical solution of the single-component model. The moving speeds of the front Uf estimated by equation (49) are 8.296×10^-2, 2.624×10^-2, and 8.296×10^-3 m/day at times 1, 10, and 100 day, respectively, indicating that the moving speed decreases rapidly at early time and then slowly as time increases

3.2 Case 2: Initial mole fraction

In this case, benzene is considered to be component one and toluene is component two in the two-component model. Figures 4(a) - 4(c) show the mole fraction distributions of

benzene and toluene versus depth when the initial mole fractions of component one are 0.2, 0.5, and 0.8, respectively, at 100 day. The evaporation front of the NAPL with u0 = 0.2, 0.5, and 0.8 reaches 0.860, 0.931, and 1.002 m below the surface, respectively. In addition, at the front u1 = 0.123 and u2 = 0.877 when u0 = 0.2, u1 = 0.310 and u2 = 0.690 when u0 = 0.5, and u1 = 0.498 and u2 = 0.502 when u0 = 0.8. The figures show that the depth of the front increases with the initial mole fraction of benzene, representing the moving speed of the front depends on the initial mole fraction. The mole fraction of benzene increases with depth until reaching u1 = u0; on the other hand, the mole fraction of toluene increases above the front but decreases below the front until reaching u₂ =1 u− ₀. These results indicate that at the front, the mole fraction of benzene decreases as time increases while that of toluene increases with time. Moreover, the mole fractions of both components are not equal to their initial values at

在文檔中 揮發性有機物在未飽和土壤中自然揮發之研究 (頁 18-0)