# Identifiability and Identification of Trace Continuous Pollutant Source

Authors: Hongquan Qu; Shouwen Liu; Liping Pang; Tao Hu
Journal: The Scientific World Journal
 (2014)

---

## Abstract

Accidental pollution events often threaten people’s health and lives, and a pollutant source is very necessary so that prompt remedial actions can be taken. In this paper, a trace continuous pollutant source identification method is developed to identify a sudden continuous emission pollutant source in an enclosed space. The location probability model is set up firstly, and then the identification method is realized by searching a global optimal objective value of the location probability. In order to discuss the identifiability performance of the presented method, a conception of a synergy degree of velocity fields is presented in order to quantitatively analyze the impact of velocity field on the identification performance. Based on this conception, some simulation cases were conducted. The application conditions of this method are obtained according to the simulation studies. In order to verify the presented method, we designed an experiment and identified an unknown source appearing in the experimental space. The result showed that the method can identify a sudden trace continuous source when the studied situation satisfies the application conditions.


---

## Body

## 1. Introduction

A pollution incident has the potential to create significant harm on numbers of people, especially when it happens in an enclosed ventilation space such as a cabin in a manned spacecraft, submarine, or aircraft. It is urgent to develop related technology to efficiently locate a sudden pollutant source so that control actions can be taken rapidly.Pollutant source identification is a process that searches a source term reversely by using limited information. It usually consists of detecting the initial emission time, estimating emission strength and finding the source location. However, inverse identification for a source is an extremely difficult challenge in the research field of inverse problems, namely, for the complexity in the coupling of source position, strength and emission time, and the incompleteness and uncertainty of observed data. With the improvement of computer processing capacity, the identification of an accidental pollutant source is becoming a more popular topic in some fields, such as research of pollution in atmospheric environment, water environment, enclosed spaces, and porous media. Many researchers attempt to use some potential method to detect contaminant source [1, 2].Most of the methods used retrospective methods. These methods use the sensor observed data to estimate the unknown pollutant source. These methods can be further classified into analytical, optimization, probabilistic, and backward computational fluid dynamics (CFD) methods.The analytical method can identify a pollutant source based on an analytical solution of a velocity field and a measured concentration distribution. It is often applied in an atmospheric environment or on surface water where a steady flow field is present [3–5].The optimization method can search for the optimal value by comparing the measured pollutant concentration data with the corresponding calculated data. The linear or nonlinear optimization and the maximum likelihood method are widely adopted. This method requires a large amount of direct calculation and this may lead to a heavy computation load. This method is successfully used to identify the source in an atmospheric environment and in underground water [6–10]. Although this method is quick and efficient for simple situations, it is difficult to obtain the analytical solution for a more complicated flow field. In order to improve the robustness of this method, some researchers attempted to use sensitivity analysis in order to reduce the uncertainty [11, 12].The probabilistic method uses probability theory to calculate the possibility that a pollution source appears at a certain position [13, 14]. By using this method, the ill-posed problem of source identification may be transformed into a well-posed problem in the extended statistical space [15, 16]. Many researchers have tried to apply this method to identify an indoor pollutant source and to optimize the placement of sensors [17–23]. Though there are some shortcomings of this method, such as the complexity of the solution of the adjoint state equations, and the fact that the uncertainty of the model cannot be overcome completely, it is indeed a promising method for identifying a source.The backward CFD method can inversely evaluate the source information by using a negative time step and a reverse flow field. The inverse CFD model is ill-posed, so regularization technique and stabilization technique should be used to improve its solution stability [24]. Some researchers carried out several studies to improve this stability [25–28]. Currently, the combination of the probabilistic method and the backward CFD method is a significant advance in the identification of a pollutant source in an enclosed space [18–20]. They used backward PDF combined with the quasireversibility method and pseudoreversibility method, respectively, to identify an instantaneous source by using the accurate observed data of a single sensor. Although this research requires an accurate velocity field and exact sensor measured data, it has proved to be promising.The above endeavors have been promoting the development of source identification. Actual pollutant source emission sometimes is a continuous and not instantaneous process, and the measured data is influenced by sensor noise. So, we develop a method to identify a sudden continuous emission source in a steady velocity field that remains unchanged or little changed during a source emission process by using single sensor information with noise. The application conditions for this method are presented and its assumptions are discussed in detail.


## 2. Model of Pollutant Transport Process

Assume that a source emits pollutant in a closed space with a constant temperature condition and a small viscous stress and then a gaseous pollutant transport process by air is governed by the ADE and momentum equations [29]:
(1)
∂
C
∂
t
+
∇
(
u
⃑
C
)
=
∇
(
d
·
grad
C
)
+
S
f
,
(2)
∂
(
ρ
u
x
)
∂
t
+
∇
(
ρ
u
x
u
⃑
)
=
-
∂
p
∂
x
+
F
x
,
∂
(
ρ
u
y
)
∂
t
+
∇
(
ρ
u
y
u
⃑
)
=
-
∂
p
∂
y
+
F
y
,
∂
(
ρ
u
z
)
∂
t
+
∇
(
ρ
u
z
u
⃑
)
=
-
∂
p
∂
z
+
F
z
,
where C is the pollutant mass concentration, mg/m3; t is the forward time, s; u
⃑ is the velocity vector, m/s; u
x, u
y, and u
z are the velocity components in x, y, and z directions of u
⃑; d is the diffusion coefficient, m2/s; S
f is the pollutant source which is a function of (p, S, t
e); p is the source position; S is the emission strength, mg/(m3
·s); ρ is the fluid density, kg/m3; p is the fluid pressure on the micro element; and F
x, F
y and F
z are the body forces in x, y, and z directions, N.We can obtain an accurate solution of the concentration field by solving the ideal gas equation and (1)~(2). However, this simultaneous solution is very complex and time-consuming. It will become even more difficult to identify a pollutant source reversely because it should solve u
⃑ and S together by using the partially known C.In order to avoid this complex and time-consuming calculation process, we assume the studied situation satisfies the following two conditions before and after a source appears: (1)u
⃑ is large enough; (2) ρ is approximately uniform; and (3) they remain unchanged or little changed.With these two assumptions, a pollutant source identification process can be transformed into a relatively simple inverse problem because (2) can be neglected and only (1) needs to be considered. We can estimate the three parameters (p, S, and t
e) of the source by solving (1) reversely when u
⃑ and a part of C are known.Equation (1) can be solved numerically. When the pollutant concentration distribution at time t is known, the concentration distribution sequent at time t
+
1 can be obtained as follows:
(3)
A
C
t
+
1
=
a
C
t
+
S
f
,
where C
t
+
1 and C
t are the column vectors of the pollutant concentration for each grid node at time t
+
1 and t, respectively; a is the coefficient matrix of C
t; and A is a diagonally dominant banded sparse matrix; hence, it must have the inverse matrix A
-
1:
(4)
C
t
+
1
=
A
-
1
a
C
t
+
A
-
1
S
f
.A can be uniquely determined when the computational zone, the mesh, and the discretization scheme are known. We chose the power-law scheme to solve the discretized equations of pollutant ADEs in this paper.


## 3. Source Identification Method

Based on the above assumptions aboutu
⃑ and ρ, we will develop a method to identify a sudden trace pollutant source by using a single sensor with noise.A sensor can detect a sudden pollution when its observed concentration rises upon a given threshold value,Cth. We define two time, t
e and t
p. t
e is the initial emission time when a source begins to emit a pollutant. t
p is the detecting time when a sensor takes to detect the existence of the pollutant. There is a delay time t
d and t
d
=
t
p
-
t
e.An identification process can be conducted aftert
pand C
sensor are known. We use the idea of multiple parameter multiple hypothesis optimization to identify the unknown three source parameters (p, S, and t
e) of source. The measured concentration sequence of sensor in the time interval [
t
p
,
t
check
] is specially taken out from C
sensor and defined by C
senCheck for the identification. There are four required steps to implement this process, as shown in Figure 1.Figure 1
Implementation steps for a source identification.(A) Multiple Hypothesis of Source Position. A set of hypotheses of source position, (p
1
,
p
2
,
…
,
p
i
,
…
p
n) is established firstly. It includes all possible source positions. The second hypothesis assumes that a virtual pollutant source with a constant emission strength, s, will appear on every possible hypothetical position p
i, respectively. For a certain p
i, a set of hypotheses, (
t
e
,
1
i
,
t
e
,
2
i
,
…
,
t
e
,
j
i
,
…
,
t
e
,
k
i
), will be established thirdly. If the backward searching time index is j, then t
e
s
,
j
i
=
t
p
-
j
×
Δ
t.(B) Optimization Only for
p
i
. The optimal three parameters of source only for p
i can be searched with the objective function of minimum characteristic distance. It includes the following steps.(a)
Whenu
⃑ and d are known, the hypothetical measured concentration sequence, (
c
1
i
,
c
2
i
,
…
,
c
j
i
,
…
,
c
k
i
), at sensor position for p
i can be calculated by forwardly solving (4) in [
t
p
,
t
check
] with the assumption of (p
i
,
s
,
t
e
,
j
i).
(b)
Assume thatu
⃑ remains unchanged or little changed during the pollutant emission process and then the steady-state observed concentration is linear with respect to emission rate [23]; hence, the approximate relative emission strength S
j
i can be obtained as follows:
(5)
S
j
i
=
S
s
=
∫
t
p
t
check
(
C
senCheck
-
C
0
)
d
t
∫
t
p
t
check
(
c
j
i
-
C
0
)
d
t
,
where s is the constant emission strength of the virtual pollutant source, mg/(m3·s); S is the estimated emission strength, S
=
S
j
i
·
s, mg/(m3·s); S
j
i is the estimated relative emission strength with the assumption of (
p
i
,
s
,
t
e
s
,
j
i
); and C
0 is the initial concentration, mg/m3.
(c)
After(
S
1
i
,
S
2
i
,
…
,
S
j
i
,
…
,
S
k
i
) are calculated by using (5), the concentration at sensor position in the time interval [
t
p
,
t
check
] can be calculated forwardly by using (4) and is denoted as the hypothetical measured concentration sequence of sensor C
j
i:
(6)
C
j
i
=
S
j
i
(
c
j
i
-
C
0
)
+
C
0
and (
C
1
i
,
C
2
i
,
…
,
C
j
i
,
…
,
C
k
i
) can be obtained with (6).
(d)
We defined a characteristic distance to describe the similarity degree betweenC
j
i and C
senCheck in [
t
p
,
t
check
]:
(7)
D
j
i
=
∑
t
=
t
p
t
check
(
C
j
,
t
i
-
C
senCheck
,
t
)
2
,
where D
j
i is the characteristic distance between C
j
i and C
senCheck.
(e)
After the set of(
D
1
i
,
D
2
i
,
…
,
D
j
i
,
…
,
D
k
i
) has been calculated with (7), we define the minimum characteristic distance D
min
⁡
i
=
min
⁡
j
=
1,2
,
…
,
k
(
D
j
i
). D
min
⁡
i will be an objective function to search the optimal three parameters of source only for the position p
i. Let “opt” express the optimal one in the set of (
j
=
1,2
,
…
,
k
) and then t
e
,
opt
i
=
t
p
-
j
opt
i
×
Δ
t and S
opt
i can be optimized accordingly to D
min
⁡
i. Now, the optimal result for the hypothetical position p
i is (
p
i
,
t
e
,
opt
i
,
S
opt
i
).(C) Minimum Characteristic Distance Set. Similarly, the set of minimum characteristic distance, (
D
min
⁡
1
,
D
min
⁡
2
,
…
,
D
min
⁡
i
,
…
,
D
min
⁡
n
), can be finally obtained for the other hypothetical positions in the set of (
p
1
,
p
2
,
…
,
p
n
).(D) Source Identification Based on Location Probability.A location probability is established to clearly point out the most probable source positions after {
D
min
⁡
1
,
D
min
⁡
2
,
D
min
⁡
3
,
…
,
D
min
⁡
n
} has been obtained. We use the location probability to show the probability of a source appearing at the ith assumed position [14]:
(8)
P
r
i
=
e
-
(
D
min
⁡
i
/
n
σ
)
2
×
(
∑
i
=
1
n
e
-
(
D
min
⁡
i
/
n
σ
)
2
)
-
1
,
where σ is the measurement standard deviation of sensor in [
t
p
,
t
check
] and n is the number of the hypothetical source position.A distribution map of location probability can be obtained after the set of location probabilities,(
P
r
1
,
P
r
2
,
…
,
P
r
i
,
…
,
P
r
n
), has been calculated by using (8). This map is directly able to show where the source appears with max location probability. We denoted max
⁡
(
P
r
i
) as the global optimal objective function and let P
r
q
=
max
⁡
(
P
r
i
), where the superscript q is the node corresponding to the max
⁡
(
P
r
i
).AfterP
r
q has been searched, the source can be located at p
q; that is, (8) is used for the source location. Hence, the global optimal source parameter is (
p
q
,
t
e
,
opt
q
,
S
opt
q
) and this is the final identification result for all hypothetical positions, (
p
1
,
p
2
,
…
,
p
i
,
…
,
p
n
).


## 4. Identifiability

The presented identification method is relatively effective because only (1) is used for estimating a source reversely. But its performance of identification could be acceptable only when the assumptions about u
⃑ and ρ are satisfied very well. If not, then an identification error will be generated. In this section, the impact of u
⃑ and ρ on the performance of identification will be discussed in detail in order to obtain the application conditions.


### 4.1. Velocity Field

The velocity field,u
⃑, has an important impact on the identification of pollutant source. In our analysis, u
⃑ is assumed to be large enough and remains unchanged or little changed in the whole process, so the source emission or other factors should not change velocity field greatly.In order to test the change of velocity fields before and after a source happening, we define a synergy degree of velocity fields to quantitatively discuss how to satisfy the above velocity field assumption according to the field synergy principle [30].For two scalar fields of the same variableφ in the same space, denoted by the field A and B, we define a nonnegative real number Φ as the synergy degree of the field A and the field B for the same variable. Φ is expressed as a percentage:
(9)
Φ
=
max
⁡
{
φ
B
(
x
,
y
,
z
)
-
φ
A
(
x
,
y
,
z
)
max
⁡
{
|
φ
A
|
}
}
,
where φ
B
(
x
,
y
,
z
) is the value of variable φat (
x
,
y
,
z
) position in the field B; φ
A
(
x
,
y
,
z
) is the value of φat (
x
,
y
,
z
) in the field A; and max
⁡
{
|
φ
A
|
} expresses the maximum absolute value of the scalar φ in the field A.Ifφ is a vector, take its components in every axis direction and obtain their corresponding synergy degrees. The largest component of synergy degrees will be the synergy degree of vector φ between the field A and the field B.u
⃑
A and u
⃑
B are the velocity fields before and after a sudden source appears. We can obtain the synergy degrees of velocity fields with the definition of (9). If Φ
≈
0, then u
⃑
A
≈
u
⃑
B, which represents the velocity field, is not changed greatly; If Φ
>
0, then u
⃑
A
≠
u
⃑
B, which represents the velocity field, is changed. The greater the value of Φ is, the greater the velocity field will be changed.One of the prerequisites to using our identification method isu
⃑
A
≈
u
⃑
B; hence, Φ should be required to be as small as possible. In order to discuss the impact of v
in or u
⃑ on Φ, we use a 2D simulation case, as shown in Figure 2, for the analysis.Figure 2
Size of the simulation space (unit: mm).It is a rectangular enclosed space, 560 mm long and 360 mm wide, with a 20 mm wide inlet and a 20 mm wide outlet. There is a pollutant source atP
0, and it is a 10 mm × 10 mm square pollution source. There is a pollutant sampling point at S
0. The direction of gravity is the direction of -
y.In our study, CO2 is the aim pollutant. Its molecular weight is larger than air. The source begins to release the CO2 pollutant with constant emission strength at 0 second after the velocity field has become steady. The source emits a trace amount of the pollutant continuously.Four cases with differentv
in and S were analyzed. According to the above definition, we can obtain the synergy degrees of velocity fields before (t
=
0) and after a sudden source appears in the gravity field space for the case in Figure 2 with different inlet velocity and source conditions. They are listed in Table 1.Table 1
Values ofΦ in different simulation cases.
Case
v
in (m/s)
S kg/(m2·s)
Φ (%)
Case 1
1.0
0.5
1.47
Case 2
0.5
0.5
10.43
Case 3
0.3
0.5
13.12
Case 4
0.3
0.25
9.34From Table1, we can see that(1)
the synergy degree of velocity fieldsΦ will increase with the decrease of v
in. For example, when S
=
0.5 kg/(m2·s), Φ is changed from 1.47% to 13.12% corresponding to v
in changing from 1 m/s to 0.3 m/s.
(2)
If we assume that the application condition for the source identification isΦ
<
11% (in fact, Φ should be as small one as possible). Φ
=
10.43% in Case 2, in which v
in
=
0.5 m/s and S
=
0.5 kg/(m2·s), and it satisfies the assumption condition of Φ
<
11%. While for the same S and a smaller inlet velocity (v
in
=
0.3 m/s) in Case 3, its Φ
=
13.12% and it does not satisfy Φ
<
11%. But if S becomes smaller, for example, S
=
0.25 kg/(m2·s) in Case 4 with v
in
=
0.3 m/s, then its Φ
=
9.34% and it satisfies Φ
<
11% again. So, there is a close relationship between v
in and S. Therefore, a larger inlet velocity is more ideal to satisfy the assumption condition of application than a smaller one and will ensure that the source could be identified.


### 4.2. Momentum Equation

In this section, we will discuss in what situation the impact of the momentum equation can be ignored.If only considering thez-axis direction (u
x
=
u
y
=
0), the one-dimensional flow velocity u
z is taken as a fixed value and μ is very small for the air and then (2) can be written as follows:
(10)
∂
ρ
∂
t
=
-
1
u
z
∂
(
ρ
u
z
)
∂
z
-
1
u
z
(
∂
p
∂
z
+
ρ
g
)
.According to the assumptions aboutu
⃑ and ρ, we know that the absolute value of -
1
/
u
z is smaller when u
z is greater; thus, the influence of the second term on the right (gravity and air buoyancy) becomes weaker. If we assume that ρ remains approximately uniform in the space, then the first term on the right of (10) can be ignored. Hence, the flow velocity and the amount of emitted pollutant are the key factors for effectively weakening the influence of the momentum equation.We denoteC
t and C
ias the true concentration and the calculated concentration at S
0 in Figure 2, respectively. C
t can be obtained by solving the simultaneous equations of continuity, advection-dispersion, momentum, and energy equations. C
i can be calculated by only solving (4). Because it is easier and more efficient to calculate C
i than to calculate C
t, we will compare C
i with C
t and determine in what situation is it more suitable to use C
i instead of C
t for the identification.The above four cases in Table1 are also used again for the comparison. (
C
i
-
C
t
) is used to show the difference of C
i and C
t. Figures 3(a)~6(a) show the comparison results of (
C
i
-
C
t
) for the four cases. We let C
t stand for the virtual measured concentration at S
0 with an observed noise. We assume that the source is unknown and estimate it by using the presented method. Figures 3(b)~6(b) show the location results for the four cases. The deep red region indicates higher location probability and the blue region indicates that zero location probability.Results for Case 1 (v
in
=
1 m/s, S
=
0.5 kg/m2·s, Φ
=
1.47%).
(a)
Concentration comparison
(b)
Location resultWe draw some conclusions from Figures3~6.(1)
The difference ofC
t and C
i is very small in most time when the inlet velocity is v
in
≥
0.5 m/s, as shown in Figures 3(a)~4(a). This implies that the influence of the momentum equation can be negligible. In this situation, we can use C
i instead of C
t for identification. The location results in Cases 1 and 2 can be acceptable because the deep red region is concentrated near the inlet and near the actual source position, as shown in Figures 3(b) and 4(b). But the difference of (
C
i
-
C
t
) and Φ are smaller in Case 1 than Case 2; then the location result in Case 1 is better than the one in Case 2.
(2)
Whenv
in
<
0.5 m/s, the difference between C
t and C
i becomes large as shown in Figure 5(a). This implies that the momentum equation has played a leading role in the CO2 transmission and can not be ignored. We can not obtain a reasonable location result, as shown in Figure 5(b). However, if the source strength becomes tracer, for example, S
=
0.25 kg/(m2·s) in Case 4, then we can obtain a more reasonable location result, as shown in Figure 6(b). This is because Φ in Case 4 is equal to 9.34% and it is less than the one in Case 3.Results for Case 2 (v
in
=
0.5 m/s, S
=
0.5 kg/m2·s, and Φ
=
10.43%).
(a)
Concentration comparison
(b)
Location resultResults for Case 3 (v
in
=
0.3 m/s, S
=
0.5 kg/m2·s, and Φ
=
13.12%).
(a)
Concentration comparison
(b)
Location resultResults for Case 4 (v
in
=
0.3 m/s, S
=
0.25 kg/m2·s, and Φ
=
9.34%).
(a)
Concentration comparison
(b)
Location result


### 4.3. Application Conditions

According to the above analyses, we know the following.(1)
The value of the inlet velocity or the velocity field is the key to satisfy the application assumptions of the presented method. An inlet velocity should be large enough to ensure that the source is effectively identified.
(2)
The inlet velocity and the source emission strength should be matching. The source emission process should not disturb the velocity field. So, a source with strong emission strength could be located in a large velocity field, but it cannot be located accurately in a small velocity field. For a small velocity field, a source needs to be one with small emission strength to be estimated accurately.Therefore, the identification method is actually suitable for estimating a trace continuous pollutant source in a steady velocity field. In the actual application, we can refer to this idea to analyze a studied situation and determine its application conditions.


## 4.1. Velocity Field

The velocity field,u
⃑, has an important impact on the identification of pollutant source. In our analysis, u
⃑ is assumed to be large enough and remains unchanged or little changed in the whole process, so the source emission or other factors should not change velocity field greatly.In order to test the change of velocity fields before and after a source happening, we define a synergy degree of velocity fields to quantitatively discuss how to satisfy the above velocity field assumption according to the field synergy principle [30].For two scalar fields of the same variableφ in the same space, denoted by the field A and B, we define a nonnegative real number Φ as the synergy degree of the field A and the field B for the same variable. Φ is expressed as a percentage:
(9)
Φ
=
max
⁡
{
φ
B
(
x
,
y
,
z
)
-
φ
A
(
x
,
y
,
z
)
max
⁡
{
|
φ
A
|
}
}
,
where φ
B
(
x
,
y
,
z
) is the value of variable φat (
x
,
y
,
z
) position in the field B; φ
A
(
x
,
y
,
z
) is the value of φat (
x
,
y
,
z
) in the field A; and max
⁡
{
|
φ
A
|
} expresses the maximum absolute value of the scalar φ in the field A.Ifφ is a vector, take its components in every axis direction and obtain their corresponding synergy degrees. The largest component of synergy degrees will be the synergy degree of vector φ between the field A and the field B.u
⃑
A and u
⃑
B are the velocity fields before and after a sudden source appears. We can obtain the synergy degrees of velocity fields with the definition of (9). If Φ
≈
0, then u
⃑
A
≈
u
⃑
B, which represents the velocity field, is not changed greatly; If Φ
>
0, then u
⃑
A
≠
u
⃑
B, which represents the velocity field, is changed. The greater the value of Φ is, the greater the velocity field will be changed.One of the prerequisites to using our identification method isu
⃑
A
≈
u
⃑
B; hence, Φ should be required to be as small as possible. In order to discuss the impact of v
in or u
⃑ on Φ, we use a 2D simulation case, as shown in Figure 2, for the analysis.Figure 2
Size of the simulation space (unit: mm).It is a rectangular enclosed space, 560 mm long and 360 mm wide, with a 20 mm wide inlet and a 20 mm wide outlet. There is a pollutant source atP
0, and it is a 10 mm × 10 mm square pollution source. There is a pollutant sampling point at S
0. The direction of gravity is the direction of -
y.In our study, CO2 is the aim pollutant. Its molecular weight is larger than air. The source begins to release the CO2 pollutant with constant emission strength at 0 second after the velocity field has become steady. The source emits a trace amount of the pollutant continuously.Four cases with differentv
in and S were analyzed. According to the above definition, we can obtain the synergy degrees of velocity fields before (t
=
0) and after a sudden source appears in the gravity field space for the case in Figure 2 with different inlet velocity and source conditions. They are listed in Table 1.Table 1
Values ofΦ in different simulation cases.
Case
v
in (m/s)
S kg/(m2·s)
Φ (%)
Case 1
1.0
0.5
1.47
Case 2
0.5
0.5
10.43
Case 3
0.3
0.5
13.12
Case 4
0.3
0.25
9.34From Table1, we can see that(1)
the synergy degree of velocity fieldsΦ will increase with the decrease of v
in. For example, when S
=
0.5 kg/(m2·s), Φ is changed from 1.47% to 13.12% corresponding to v
in changing from 1 m/s to 0.3 m/s.
(2)
If we assume that the application condition for the source identification isΦ
<
11% (in fact, Φ should be as small one as possible). Φ
=
10.43% in Case 2, in which v
in
=
0.5 m/s and S
=
0.5 kg/(m2·s), and it satisfies the assumption condition of Φ
<
11%. While for the same S and a smaller inlet velocity (v
in
=
0.3 m/s) in Case 3, its Φ
=
13.12% and it does not satisfy Φ
<
11%. But if S becomes smaller, for example, S
=
0.25 kg/(m2·s) in Case 4 with v
in
=
0.3 m/s, then its Φ
=
9.34% and it satisfies Φ
<
11% again. So, there is a close relationship between v
in and S. Therefore, a larger inlet velocity is more ideal to satisfy the assumption condition of application than a smaller one and will ensure that the source could be identified.


## 4.2. Momentum Equation

In this section, we will discuss in what situation the impact of the momentum equation can be ignored.If only considering thez-axis direction (u
x
=
u
y
=
0), the one-dimensional flow velocity u
z is taken as a fixed value and μ is very small for the air and then (2) can be written as follows:
(10)
∂
ρ
∂
t
=
-
1
u
z
∂
(
ρ
u
z
)
∂
z
-
1
u
z
(
∂
p
∂
z
+
ρ
g
)
.According to the assumptions aboutu
⃑ and ρ, we know that the absolute value of -
1
/
u
z is smaller when u
z is greater; thus, the influence of the second term on the right (gravity and air buoyancy) becomes weaker. If we assume that ρ remains approximately uniform in the space, then the first term on the right of (10) can be ignored. Hence, the flow velocity and the amount of emitted pollutant are the key factors for effectively weakening the influence of the momentum equation.We denoteC
t and C
ias the true concentration and the calculated concentration at S
0 in Figure 2, respectively. C
t can be obtained by solving the simultaneous equations of continuity, advection-dispersion, momentum, and energy equations. C
i can be calculated by only solving (4). Because it is easier and more efficient to calculate C
i than to calculate C
t, we will compare C
i with C
t and determine in what situation is it more suitable to use C
i instead of C
t for the identification.The above four cases in Table1 are also used again for the comparison. (
C
i
-
C
t
) is used to show the difference of C
i and C
t. Figures 3(a)~6(a) show the comparison results of (
C
i
-
C
t
) for the four cases. We let C
t stand for the virtual measured concentration at S
0 with an observed noise. We assume that the source is unknown and estimate it by using the presented method. Figures 3(b)~6(b) show the location results for the four cases. The deep red region indicates higher location probability and the blue region indicates that zero location probability.Results for Case 1 (v
in
=
1 m/s, S
=
0.5 kg/m2·s, Φ
=
1.47%).
(a)
Concentration comparison
(b)
Location resultWe draw some conclusions from Figures3~6.(1)
The difference ofC
t and C
i is very small in most time when the inlet velocity is v
in
≥
0.5 m/s, as shown in Figures 3(a)~4(a). This implies that the influence of the momentum equation can be negligible. In this situation, we can use C
i instead of C
t for identification. The location results in Cases 1 and 2 can be acceptable because the deep red region is concentrated near the inlet and near the actual source position, as shown in Figures 3(b) and 4(b). But the difference of (
C
i
-
C
t
) and Φ are smaller in Case 1 than Case 2; then the location result in Case 1 is better than the one in Case 2.
(2)
Whenv
in
<
0.5 m/s, the difference between C
t and C
i becomes large as shown in Figure 5(a). This implies that the momentum equation has played a leading role in the CO2 transmission and can not be ignored. We can not obtain a reasonable location result, as shown in Figure 5(b). However, if the source strength becomes tracer, for example, S
=
0.25 kg/(m2·s) in Case 4, then we can obtain a more reasonable location result, as shown in Figure 6(b). This is because Φ in Case 4 is equal to 9.34% and it is less than the one in Case 3.Results for Case 2 (v
in
=
0.5 m/s, S
=
0.5 kg/m2·s, and Φ
=
10.43%).
(a)
Concentration comparison
(b)
Location resultResults for Case 3 (v
in
=
0.3 m/s, S
=
0.5 kg/m2·s, and Φ
=
13.12%).
(a)
Concentration comparison
(b)
Location resultResults for Case 4 (v
in
=
0.3 m/s, S
=
0.25 kg/m2·s, and Φ
=
9.34%).
(a)
Concentration comparison
(b)
Location result


## 4.3. Application Conditions

According to the above analyses, we know the following.(1)
The value of the inlet velocity or the velocity field is the key to satisfy the application assumptions of the presented method. An inlet velocity should be large enough to ensure that the source is effectively identified.
(2)
The inlet velocity and the source emission strength should be matching. The source emission process should not disturb the velocity field. So, a source with strong emission strength could be located in a large velocity field, but it cannot be located accurately in a small velocity field. For a small velocity field, a source needs to be one with small emission strength to be estimated accurately.Therefore, the identification method is actually suitable for estimating a trace continuous pollutant source in a steady velocity field. In the actual application, we can refer to this idea to analyze a studied situation and determine its application conditions.


## 5. Experiment Study

### 5.1. Experimental System

An experimental device was set up to test the proposed method. CO2 is used as the aim pollutant. CO2 gas injected into the enclosed space from a point is supplied by a compressed gas cylinder through a pressure relief valve and a mass flow controller. A QIC-20 mass spectrometer made by the UK company Hiden Analytical Ltd. was used to measure CO2 concentrations in real time at the sampling point.The photo of the enclosed space is shown in Figure7 and its size is given in Figure 2. It is a flat cavity topped and bottomed by two plexiglass flat panels and four-side sealing strips. Due to the thickness of the enclosed space being only 8 mm and relatively very thin, the flow velocity, which is formed by a small centrifugal fan, in the thickness direction can be considered approximately zero; thus, the flow field is close to a 2D case.Figure 7
Experimental device photo of 2D space.


### 5.2. Source Identification

During the experiment, the inlet velocity was 1.68 m/s. CO2 concentration curve measured at S
0 is shown in Figure 8 when the CO2 was suddenly injected from P
0. The source emission strength was constant and its actual emission strength was 16.37 mg/s.Figure 8
Measured CO2 concentration data.For the experimental space in Figure7, we know that(1)
v
in
=
1.68 m/s and it is larger than 1.0 m/s of Case 1 in Table 1. At the same time, its Φ is larger than 1.47% and u
⃑ is steady. Hence the assumption of u
⃑ is satisfied very well;
(2)
S
=
16.37 mg/s and it is less than 50 mg/s in Case 1. This means that the emission amount of source is small and the assumption of ρ is satisfied.Therefore, this experimental situation satisfies the application conditions in Section4.3.We apply the identification method to identify the source reversely by using the experimental data in Figure8. The other parameters are t
p
=
3.3 s, t
check
-
t
p
=
10 s, C
0
=
655 mg/m3, and Δ
t
=
0.1 s.The estimated emission strength is 20.862 mg/s, and the final distribution map of location probability is shown in Figure9. From Figure 9, we know that the position having the largest location probability is very near to the actual pollutant source position.Figure 9
Distribution map of location probability for the experimental case.


## 5.1. Experimental System

An experimental device was set up to test the proposed method. CO2 is used as the aim pollutant. CO2 gas injected into the enclosed space from a point is supplied by a compressed gas cylinder through a pressure relief valve and a mass flow controller. A QIC-20 mass spectrometer made by the UK company Hiden Analytical Ltd. was used to measure CO2 concentrations in real time at the sampling point.The photo of the enclosed space is shown in Figure7 and its size is given in Figure 2. It is a flat cavity topped and bottomed by two plexiglass flat panels and four-side sealing strips. Due to the thickness of the enclosed space being only 8 mm and relatively very thin, the flow velocity, which is formed by a small centrifugal fan, in the thickness direction can be considered approximately zero; thus, the flow field is close to a 2D case.Figure 7
Experimental device photo of 2D space.


## 5.2. Source Identification

During the experiment, the inlet velocity was 1.68 m/s. CO2 concentration curve measured at S
0 is shown in Figure 8 when the CO2 was suddenly injected from P
0. The source emission strength was constant and its actual emission strength was 16.37 mg/s.Figure 8
Measured CO2 concentration data.For the experimental space in Figure7, we know that(1)
v
in
=
1.68 m/s and it is larger than 1.0 m/s of Case 1 in Table 1. At the same time, its Φ is larger than 1.47% and u
⃑ is steady. Hence the assumption of u
⃑ is satisfied very well;
(2)
S
=
16.37 mg/s and it is less than 50 mg/s in Case 1. This means that the emission amount of source is small and the assumption of ρ is satisfied.Therefore, this experimental situation satisfies the application conditions in Section4.3.We apply the identification method to identify the source reversely by using the experimental data in Figure8. The other parameters are t
p
=
3.3 s, t
check
-
t
p
=
10 s, C
0
=
655 mg/m3, and Δ
t
=
0.1 s.The estimated emission strength is 20.862 mg/s, and the final distribution map of location probability is shown in Figure9. From Figure 9, we know that the position having the largest location probability is very near to the actual pollutant source position.Figure 9
Distribution map of location probability for the experimental case.


## 6. Conclusions

A method to identify a trace continuous source in a steady flow field was proposed in this paper. The identification method is realized by globally searching the optimal values with the objective function of maximum location probability.In order to ensure the performance of identification, the assumptions aboutu
⃑ and ρ are discussed in detail and the application conditions for the identification method are obtained. A synergy degree of velocity fields Φ is introduced to analyze the application conditions of velocity field quantitatively. The value of Φ should be as small as possible in order to satisfy this assumption of the velocity field. Based on this premise, we draw the following application conditions from the simulation studies:(1)
the value of the velocity field (or the inlet velocity) is one key factor. A largeru
⃑ is more ideal in satisfying the application conditions;
(2)
The source emission strengthS is another key factor in addition to u
⃑. The source emission process can not disturb u
⃑ greatly and S and u
⃑ should be matching each other. Therefore, the identification method is suitable for estimating a trace continuous pollutant source in a steady velocity field.In order to investigate the practical application effect of the above studies, the experimental system was set up and subsequently a sudden pollution experiment was conducted. The presented identification method was used. The studies show that the sudden trace continuous source can be successfully located when the experimental situation satisfies the application conditions.



---

*Source: 215104-2014-04-27.xml*