ABSTRACT

This paper presents the fundamentals of direct inverse modeling using CFD simulations to detect air pollution sources in urban areas. Generally, there are four techniques used for detecting pollution sources: the analytical technique, the optimization technique, the probabilistic technique, and the direct technique. The study discusses the potentialities and limits of each technique, where the direct inverse technique is focused. Two examples of applying the direct inverse technique in detecting pollution source are introduced. The difficulties of applying the direct inverse technique are investigated. The study reveals that the direct technique is a promising tool for detecting air pollution source in urban environments. However, more efforts are still needed to overcome the difficulties explained in the study.

1. Introduction

Nowadays, inverse modeling seems to be a promising topic in terms of environmental research. The importance of reverse modeling arises from the increased numbers of pollution sources due to continuing industrial expansion coupled with population growth, especially in large cities. In addition, terrorist attacks are becoming more frequent and deadlier, such as the sarin gas attacks on the Tokyo subway in 1995, which resulted in 12 deaths and over 6000 injured [1,2], and the terrorist attacks on New York City and Washington DC on September 11, 2001 and the following anthrax dispersion by mail. Moreover, in the US, fire accidents occur in structures at the rate of one every 61 s, and in particular residential fires occur every 79 s. Nationwide, there was a civilian fire-related death every 156 minutes and a civilian fire-related injury every 28 minutes [3]. All of these events have confirmed that the terroristic attacks are no longer hypothesis but a reality. From this perspective, the ability to predict pollution sources characteristics: location, strength, and release time has become a necessity in order to create a complete picture of the air quality conditions within the release domain and to ensure the public’s safety. Inverse modeling using Computational Fluid Dynamics (CFD) simulations is one efficient and promising tool in this regard.

Traditionally, CFD is used to explore the cause-effect relationship that is called forward-time modeling. Inverse modeling is to find the causal characteristics such as the pollutant source location and strength from the finite effectual information like the distributions of airflow and pollutant concentrations. Generally, in modeling process, there are three different types of problems [4]:

• The forward-time problem: given input and system parameters, find out output of the model.

• The reconstruction problem: given system parameters and output, find out which input has led to this output.

• The identification problem: given input and output, determine the system parameters that agree with the relationship between input and output.

Type 1 problem is an oriented cause-effect sequence. In this sense, type 2 and 3 problems are inverse problems because they are the problems of finding unknown causes of known consequences.

In order to determine the pollution source within a domain through the use of inverse modeling, some measurements such as concentration distributions are needed to infer the values of the inputs or parameters that characterize the system [5]. In reality, the concentration field is obtained from data measured by concentration sensors. The conventional CFD codes predict airflow based on given inputs (boundary conditions) and system parameters (building and system characteristics), which is a forward-time problem. Then, the reverse simulation is carried out through the solution of the reversed transport equation using the measured concentration field and the solved air flow field.

By carrying out inverse CFD simulations, one can detect air pollution source (or release) characteristics. Identification of pollution sources immediately after release is a matter of urgency in order to ensure the public’s safety. Since, inverse modeling is one of the efficient and promising tools in this regard; this paper sheds light on the existing inverse techniques and explores the feasibility of using these techniques for urban environment study. Generally, there are four techniques used for detecting pollution sources the analytical technique, the optimizition technique, the probabilistic technique, and the direct inverse technique. First, the present study discusses the potentialities and limits of each technique. Then, the direct inverse modeling technique, as a promising tool to detect air pollution sources, is focused.

2. Different Techniques of Inverse Modeling

In the present study, the existing inverse modeling techniques that have been employed to detect pollutant sources in urban areas are reviewed. Figure 1 shows classification of inverse CFD modeling techniques. Each technique has its own advantage and disadvantage. In the following sections, the four techniques will be discussed in details.

2.1. The Analytical Approach

The analytical approach requires analytical solution of the distributions of airflow and pollutant concentrations. The causal characteristics are then inversely solved. The analytical approach has been successfully applied to multi-dimensional heat conduction problems [6]. In groundwater transport, the analytical approach has been used to solve contaminant transport in one-dimensional flow [7] or in 2-D uniform flow [8]. The analytical approach has also been used to solve an inverse atmospheric transport problem in three-dimensional uniform

flow [9]. Using analytical and graphical techniques, Islam [10,11] worked out a method for determining an unknown emission source in urban areas based on the well known Gaussian Plume Model (GPM). He assumed that the GPM describes the dispersion process fully.

The above discussion shows that the analytical approach can be accurate and efficient. However, it is only used for simple problems. Accordingly, the application of the analytical approach in problems of complex geometry is very limited.

2.2. The Optimization Approach

In the optimization approach, forward-time modeling is used to obtain the effectual data (such as distributions of pollutant concentrations) based on all possible causal characteristics. Then the approach optimizes a solution that is best-fitted with the corresponding measured data. This approach has been widely applied in identifying groundwater pollution source as linear optimization method [12], maximum likelihood method [13], and nonlinear optimization method [14]. As an example, Arvelo et al. (2002) [6] have used a modified multi-zone model to study the optimal placement of chemical/biological warfare agent sensors in a building with nine offices and a hallway. The optimal sensor locations should be in the spaces where most possible sources of chemical/biological warfare agent could be located so the sensors can detect the agent in least amount of time. Different from the previous researchers, they used the genetic algorithm to interpret the computed data to locate the sources. Because plausible combinations of possible causal characteristics are huge, this approach involves a large amount of modeling.

2.3. The Probabilistic Approach

The probabilistic approach also does forward-time modeling. The approach uses probability to express a possible causal aspect. In groundwater transport, Bagtzoglou et al. (1992) [15] calculated the possibility of a contaminant source in groundwater by reversing only the convective contaminant transport. Other researchers [16] used Bayesian theory to interpret the possibility of each contaminant source. In enclosed environment, Sohn et al. (2002) [17] also used Bayesian probability model to identify the contaminant source in a five-room building. They used a multi-zonal model to calculate the airflow and contaminant transport. Since the multi-zone model can only provide some macroscopic information about the contaminant transport, it is necessary to run CFD simulations if more accurate and detailed information is needed. Accordingly, both the probability and optimization approaches need a huge amount of forward-time modeling.

2.4. The Direct Inverse Approach

The direct inverse approach solves inverse problems by reversing directly the governing equations that describe cause-effect relations. Kato et al. (1992) applied the reverse tracking of flow field over time, to obtain the residual lifetime of air at a point. In addition, Kato et al. [18-20] used the same technique to assess local pollution from upwind regions with backward trajectory analysis of the flows in an atmospheric environment. They reversed only the pollutant transport by convection and neglected the effect of diffusion. The method has also been used in groundwater contaminant transport [15,21], where convection transport of pollutants is solved with reversed velocity field and the diffusion was left unchanged.

In the following sections, the direct inverse modeling technique is investigated in details. The fundamentals of the technique are introduced. Then, the difficulties in applying the technique are discussed. Finally, two examples of applying the direct inverse technique in detecting pollution source in urban area are introduced.

3. Fundamentals of Direct Inverse Modeling

Pollutant concentrations are predicted based on the convection-diffusion equation including meteorological data, transport diffusion, and the relevant emissions (S). It can be written as follows [22]:

(1)

where c is the pollutant concentration (kg/kg)t is the time (s)u_i is the Cartesian components of the velocity (m/s)x_i is the Cartesian coordinates (m)K is the concentration diffusivity coefficient (kg/m/s)ρ is the air density (kg/m³), and S is the pollutant source strength (kg/m³/s).

In order to explain the principles of inverse modeling, Equation (1) can be simply integrated based on one dimensional flow as shown in Figure 2 with using upwind scheme for space:

(2)

Assuming K and ρ to be constant, the integration result becomes:

(3)

Rearranging Equation (3), then:

(4)

To perform the integration with respect to time, the concentration-time relation is needed. Usually, the following equation is assumed [23]:

(5)

where f is a weighting factor used to combine the concentration value at the new time step with both of new and old concentration values.

By applying fully implicit scheme (i.e. f = 1), Equation (5) can be rearranged to the form:

(6)

Clearly, in order to calculate the concentration at previous time, one can simply apply a negative time step in Equation (1). However, the negative time step makes the coefficient of c_P(t) greater than one, which amplifies the calculation errors as time elapsed in negative direction. Thus, Equation (1) is unstable in inverse modeling. Accordingly, in order to solve the reverse transport equation with numerical stability, another approach has to be used. One of these approaches is to use a slightly different equation instead of the original equation [24-26].

The inverse direct technique requires less-efforts compared with the probabilistic technique and the optimization technique. However, the instability problem caused by the negative diffusion term represents the main difficulty when applying such technique. This will be investigated in the following section.

4. Solution Instability in Direct Approach

Since the reversed governing equations are unstable due to the negative diffusion term, many researchers have worked out to improve the solution stability through imposing a bound on the solution. The regularization technique (or the stabilization technique) has been used to improve the solution stability. The solution is obtained by minimizing the objective function with a regularized term. The stabilization technique introduces some stabilization terms into the reversed governing equations or solves some auxiliary equations to improve the solution stability. For example, Zhang and Chen (2007) [26] and Liu et al. (2011) [27] changed the diffusion term in the reversed transport equation as (i.e. Equation (1)):

(7)

where ε is a stabilization coefficient.

The equation is numerically stable compared with Equation (6), since the coefficient of c_P can be smaller than one, where:

(8)

The required condition is just to make ε larger than u_e∆x³/6, which is easy to be satisfied. Then the numerical scheme becomes stable and the reversed transport equation becomes dispersive.

By comparing the stabilization and the regularization methods, Skaggs and Kabala (1995) [25] concluded that the stabilization method used significantly less computational effort than the regularization method, although the stabilization method might provide slightly inferior results.

In the present study, in order to make the equation solvable with numerical stability, this study proposes the use of a filter to deal with negative concentration gradients in order to avoid unrealistic solutions. Then, Equation (1) can be rewritten as:

(9)

where the symbol “ˆ” means that this term is filtered.

The filter for the diffusion term is set in such a way that the negative concentration values are changed to the minimum positive value of the concentration within the solution domain. This can be expressed as:

(10)

Indeed, the filtration process of the diffusion term in the reversed transport equation affects the solution accuracy and decreases the method’s effectiveness in identifying the sources of pollution. This will be investigated in the following examples.

5. Examples of Applying Direct Inverse Modeling in Identifying Source Locations

5.1. Simple Laminar Flow

A 3-D numerical model for modeling a simple laminar flow is designed as shown in Figure 3. The dimensions of the computational domain are 90 m × 30 m × 30 m (L × W × H). The minimal mesh size above the ground is set as 0.18 m, and the vertical growing factor for mesh points is 1.05. The total number of cells is 594,000. A pollution source of strength 0.01 kg/m³/s and dimensions 0.5 m × 0.5 m × 0.5 m is considered at a location 12 m from the inlet face at a height of 2 m.

5.2. Turbulent Flow around a Single Building

The second example for identifying the source of pollution in urban areas using direct inverse modeling is shown in Figure 4, in which wind flows around a single building. The building dimensions are 10(x) × 10(y) × 30(z) m and the computational domain dimensions are 300(x) × 110(y) × 70(z) m with a mesh size of 444,800 cells.

Two different locations for the source are considered in order to examine the effect of pollutant source location on the prediction accuracy of the reverse technique: location (I), which is 2.5 m upwind the building, and location (II), which is 2.5 m downwind the building. In both cases, the source strength is set at 0.01 kg/m³/s with dimensions of 0.5 m × 0.5 m × 0.6 m.

5.3. Numerical Simulation

Numerical simulations are carried out using the CFD code Star-CD, based on the finite-volume discretization method. During forward-time simulation, steady-state analysis is adopted for the flow field and the Monotone Advection and Reconstruction Scheme (MARS) [28] is applied to the spatial difference. The standard k-ε model is used to simulate the turbulence effects and the pressure/velocity linkage is solved via the SIMPLE algorithm [23].

At the inflow boundary, a constant flux layer is assumed

for the turbulent energy, and turbulent intensity is assumed to be 10% of the inflow wind velocity at a representative height (z_o) of 74.6 m. Free slip condition is applied to the top and side boundaries. The generalized logarithmic law with parameter E = 9 (m) is applied to building walls and the ground surface as smooth walls.

Tables 1 and 2 summarize the parameters used in the numerical simulations for both examples, together with the applied boundary conditions.

Once the flow field is solved, it is considered steady. In fact, airflow characteristics within outdoor environments are unsteady due to fluctuations in both speed and direction. However, in the present stage of this study, fluctuations in the applied wind are not considered in the analysis. Accordingly, the wind flow is treated as steady. After solving the flow field, the concentration fields in both forward-time and reverse simulations are solved against time. The wind flow distribution together with the pollutant concentrations of the forward-time simulation are used as initial conditions for the reverse simulation. A time step of ∆t = 0.1 s is used with the implicit scheme.

In forward-time simulation, a pollutant source of strength 0.01 kg/m³/s is released in the period from t = 0 to 30 s. Then, the solution of the transport equation is continued with the absence of the source until t = 100 s

for the case of laminar flow and until t = 200 s for the case of wind flow around a single building. Figure 5 shows the characteristics of the pollutant release as a function against time.

In reverse simulation, the flow field calculated by forward-time simulation is reversed at first and the source term is set to zero. Then, the transport equation is solved from the moment of t = 100 s in the reverse time direction. The time at which the reversed simulation is stopped is not known in reality. In other words, starting from the moment of solving the transport equation in the reverse time direction, the time t = 0 at which the calculations have to be stopped is not known in advance. Since the present examples are introduced just to explain the technique of direct inverse modeling, the end time is assumed.

6. Results and Discussions

6.1. Case of Laminar Flow

Figure 6 shows the concentration fields for two different constant inflow velocities of 0.01 and 0.5 m/s. In subplots (a) and (b), the concentration fields at t = 1 s and 100 s obtained by forward-time simulation are presented, while the subplot (c) shows the distribution of pollutant concentrations obtained through reverse simulation.

With the steady-state airflow pattern, forward-time CFD simulation was used to calculate pollutant concentration distribution at t = 100 s, where a pollution source was released from t = 0 to 30 s. The distributions of steady-state airflow and transient pollutant concentration at t = 100 s were used as initial data for the inverse CFD modeling. The inverse CFD modeling calculated backwards pollutant transport from t = 100 s to 0 s is shown in Figure 6(c).

In order to determine the pollutant release location, the distribution of pollutant concentration should be in a small region around the release source as shown in Figure 6(c), and by using the maximum pollutant concentration over all locations at this instance (the peak pollutant

concentration), one could identify the pollution source. The peak concentration computed by inverse CFD modeling in Figure 6(c) clearly shows the position of the pollutant source.

By comparing the two concentration fields of u = 0.01 m/s and 0.5 m/s given in Figure 6(c), the case where u = 0.5 m/s shows a more dispersive concentration field compared with the case of u = 0.01 m/s. This can be attributed to the increased plume size with the increase in wind velocity. However, the direct reverse simulation appears to effectively identify the release sources in both cases.

6.2. Case of Turbulent Flow around a Single Building

Figure 7 shows the computed airflow pattern around a single building under steady state conditions. Figure 7(a) shows the flow field calculated by the forward-time simulation and Figure 7(b) shows the reversed flow field which was used to carry out the reverse simulation. In Figure 7(a), a symmetrical flow field is shown around the building, where two identical circulation regions are formed around the building.

It is expected that the effectiveness of inverse direct modeling in identifying the pollution source is affected by its location. So, as mentioned before, two locations for the pollutant source were examined here, namely locations (I) and (II). The first location is upwind of the building and the second location is in the wake region downwind of the building.

Source location (I):

The concentration fields around the building obtained for the source at location (I) are shown in Figure 8. At t = 1, shown in Figure 8(a), the plume starts to diffuse where the concentration field makes contact with the upwind face of the building. The plume then travels with the wind downwind of the building. The pollutant continues to be released until the moment t = 30 s. At that time, the emission source was stopped and the scalar transport equation was solved against time using inverse

direct simulation without the source. At t = 200 s, two identical regions were formed far from the building. The conditions at such time were used as the initial conditions for the inverse modeling.

Figure 8(c) shows the concentration field obtained by the inverse modeling. In such figure, the concentration field area is wider than that of the original forward-time simulation. By detecting the location of the maximum concentration overall in the domain volume, the source location is clearly determined. It is located downwind of the building along the domain centerline. However, the figure implies two problems. The first one is that the maximum concentration occupies a wide area in front of the block, which means that estimations for the location of the source will not be 100% accurate. It is thought that the wide spread of the concentration field is attributed to what is called “false diffusion”, introduced by the first order approximation of the advection term in the upwind difference scheme [29]. The second problem is that, compared to the peak concentration obtained with the forward-time simulation (see Figure 8(b)), the source strength identified by the reverse simulation (see Figure 8(c)) is more dispersive. The reason is that the reversed transport equation is not exactly the same as the governing transport equation due to the presence of a filter. Indeed, these two problems affect the prediction accuracy of the inverse modeling and this appears clearly in the wide area of the maximum concentration, which gives a wide range of possibilities for the pollution source location. This is not considered to be the ideal distribution when a gas release position is required. At the same time, the dispersive property of the reversed transport equation renders the estimated source strength inaccurate.

Source location (II):

Figure 9 shows the concentration fields around the building when the pollutant source is at location (II). Figure 9(a) shows the source location, which is 2.5 m downwind of the building’s rear wall. As with the case for location (I), the pollutant is emitted in the period from t = 0 to 30 s, and then the transport equation was solved in the absence of the source term until the moment t = 200 s. In subplot (b), two high concentration regions are formed behind the building. However, diffusion of the pollutant in this case is limited to a narrow region compared to the case for location (I).

Figure 9(c) shows the concentration field obtained using inverse modeling. The figure demonstrates that the concentration increases in the direction of the reversed wind from right to left. Also, the figure shows that the peak concentration region occupies a narrow region compared with the case for location (I). This can be attributed to the presence of the source in the wake region where a lower wind velocity exists and convection is weak. In such case, the location of the pollution source as estimated by the reverse simulation is not clear. In Figure 9(c), there are two peak concentration regions near the edges of the block. So the location of the source is not accurately identified. This indicates that the prediction

accuracy of the inverse modeling is diminished in case of weak convection.

6.3. Prediction Accuracy Improvement

From the results of the above examples, it is clear that the accuracy of the method is limited due to the wide diffusion fields around the source location caused by the filter. Accordingly, a technique for improving the prediction accuracy of the method is needed. The proposed technique used what is called a “sink” in the reversed transport equation in order to decrease the widespread of the concentration field around the source and render identification of the source more easily. The reversed transport equation is then rewritten as:

(11)

where the last term on the right-hand side is the sink term.

In order to apply the sink term technique to identify pollution source characteristics, some criteria are needed. Such criteria have to satisfy the following requirements:

• Sufficiently general to be applied to any case.

• It should show where to put the sink within the study domain.

• The criteria should give a reasonable value for the sink strength to carry out the calculations.

• Then, the following procedure is applied which couples CFD with optimization approach.

• Setting the sink at any location within the domain with any strength (S).

• Carrying out CFD simulation to solve the scalar transport equation with any release time (T)

(12)

• At time (T), we calculate the whole concentration in the domain, which is:

• Then, the optimization approach is used to find the minimum value of the concentration within the domain.

• If the initial guess doesn’t give the minimum domain concentration, another trial with new sink location and new values for S and T selected by the optimization technique is used, until the minimum concentration within the domain is obtained.

A flow chart for the above mentioned procedure is shown in Figure 10.

7. Conclusions

This paper discusses the advantages and feasibility of using inverse CFD modeling to detect air pollution sources in urban environments. The study reviews vari-

ous inverse modeling approaches and categorizes them into four techniques: the analytical technique, the optimization technique, the probabilistic technique, and the direct inverse technique. Each technique has its own advantage and limitations. The analytical approach requires analytical solution of the distributions of airflow and pollutant concentrations. However, it is only used for simple problems, which means that the applications of such technique are very limited. The optimization approach uses forward-time modeling to obtain the effectual data based on all possible causal characteristics. The limitation of this approach is that it needs a large amount of forward-time modeling. The probabilistic approach uses probability to express a possible causal aspect. In such technique all possible causal characteristics should be known before doing the modeling. The direct inverse approach solves inverse problems by reversing directly the governing equations that describe cause-effect relations. The instability problem caused by the negative diffusion term represents the main difficulty when applying such technique. The study proposed the use of a filter to account overcome the instability caused by the negative diffusion term in the scalar transport equation. In addition, in order to decrease the wide diffusion field around the source location caused by the filter, the study proposed a “sink-term” technique in the reversed transport equation to render identification of the source more easily. As a conclusion, the study reveals that the direct inverse modeling is a promising tool for detecting air pollution source in urban environments. However, more efforts are still needed to overcome the difficulties explained in the study.

REFERENCES