# Kinetic-Monte-Carlo-Based Parallel Evolution  Simulation Algorithm of Dust Particles

Authors: Xiaomei Hu; Zhifeng Xu; Hongxia Cai; Junjun Hu
Journal: Journal of Applied Mathematics
 (2014)

---

## Abstract

The evolution simulation of dust particles provides an important way to analyze the impact of dust on the environment. KMC-based parallel algorithm is proposed to simulate the evolution of dust particles. In the parallel evolution simulation algorithm of dust particles, data distribution way and communication optimizing strategy are raised to balance the load of every process and reduce the communication expense among processes. The experimental results show that the simulation of diffusion, sediment, and resuspension of dust particles in virtual campus is realized and the simulation time is shortened by parallel algorithm, which makes up for the shortage of serial computing and makes the simulation of large-scale virtual environment possible.


---

## Body

## 1. Introduction

Many ecological environmental problems have emerged during the process of urbanization. Sedimentation of lots of surface dust in cities caused by transportation is one of them. Dust has close relation with particulates in the atmosphere [1]. Ecological system can be hurt by dust covertly over a long period, while it is human body that dust can do most direct and greatest harm to. Dust (especially suspended particulate matter with aerodynamic diameters less than 10 μm) has been one of the most serious air pollutants in China for many years. The surface dust can be resuspension under some certain circumstances and the contaminant will make bad influence on our body. It is shown by some researches that Pb can retard children’s intellectual development and weaken their intelligence, while these acknowledged prisoners including Cu, Cd, Cr, Zn, As, and Hg are in a position to change human beings’ nervous and respiratory system. A lot of problems are caused by the pollution of dust, such as laze, respiratory disease, and lung cancer [2]. As a result, research on the surface dust is not only a crucial aspect of environment evaluation, but has also great significance on human health when the urbanization is accelerating.The evolution process of surface dust includes sediment, diffusion, and resuspension. Nevertheless, it is quite difficult to study the evolution process of surface dust by experiment equipment currently because of the extreme complication of the evolution process. Therefore, simulation has become one of the most important means to research the dust evolution process. The current simulation research includes the diffusion of dust around buildings in cities and the relationship between dust deposition and wind power [3–7]. Dust diffusion is one of the hot topics of the dust evolution. Dust diffusion is a complex dynamic process. It represents a series of particles diffusion process in the wind, such as particles deposition, resuspension, and transit [8]. A number of investigations have been carried out for predicting turbulent diffusion around buildings using wind tunnel tests [9, 10]. Numerical methods for simulating flow and diffusion fields have developed rapidly. Particle deposition in vertical square ventilation duct flows by two different numerical models has been studied [11, 12]. The effects of particle diameter, dimensionless relaxation time, flow direction, and air speed in vertical upward and downward square duct flow on particle deposition velocities are discussed; lift and gravity have also been taken into account. Ali has investigated a time-dependent partial differential equation governing the transport of heavy dust into the atmosphere. Dust concentration is expressed in terms of a series of confluent hyper geometric functions [13]. Sun et al. investigated the microparticle deposition and distribution by employing the Eulerian approach with Reynolds stress turbulent model and a Lagrangian trajectory method [14]. With the rapid development in computer technology, computational fluid dynamics (CFD) method has matured to simulate the ventilation performance and contaminant dispersion and transfer in buildings. Saha et al. have used CFD to assess the effect of wind tunnel sizes on air velocity and concentration boundary layers and on ammonia emission estimation [15]. Tominaga et al. have used CFD to predict the air diffusion around a construction [16, 17]. Besides these, Schneider et al. proposed a semiempirical two-compartment, constant parameter model [8]. Roney and White compared the near-surface wind-tunnel fugitive dust concentration profiles arising from soil surfaces beds with those calculated numerically [18]. However, most of the studies paid attention to simulate the diffusion process of particles. The whole dust evolution process is seldom traced.A lot of methods can be used in the simulation of dynamic process, such as first principal method (FP), molecular dynamic method (MD), Monte Carlo method (MC), and finite element method (FEM) [19–22]. Among these methods, MC method has been well used in simulation of air pollution dispersion and river water pollution. Three decades ago, Monte Carlo (MC) method was used to study the dynamic process. The basic idea of the MC method is that the solution of the problem is equivalent to hypothetical statistical model parameters using random number, and the parameters are estimated by the statistical model of a sample [23]. Okin proposed MC method to simulate airborne dust diffusion model [24]. Yao et al. used MC method in simulation of emission height effects on building [25]. Tian presented the MC simulation of complex terrain effect on dust diffusion [26]. Dejun GU puts forward a model for the convective boundary layer of the line source diffusion along MC method. Xu et al. simulate gas dispersion based on Monte Carlo which could satisfy the requirement of the long-distance pipeline disasters emergency decision making [27]. The object researched by kinetic Monte Carlo (KMC) method is unbalanced or a relaxation process. Time evolution correct or not is the key factor in the simulation process; simulation time step must be the system real time step, so KMC method is an effective method of studying on the kinetic behavior. Sun et al. applied KMC analysis in environmental risk assessment of a chlorine release accident [28]. Peng and Yuan simulated gas-solid flow behavior in desulphurization tower based on KMC [29]. Because dust have a similar dynamic property with water and air pollutant, KMC method has widely been used in water flow simulation and atmospheric diffusion field. The kinetic Monte Carlo (KMC) related to the MC has the advantage of simulating in a long period. In addition, KMC is a stochastic process, which makes it fit for evolution simulation of the surface dust particles. In the current conditions, it is hard to simulate the KMC evolution of dust particles in a large-scale urban environment by one computer. Therefore, the parallel computing method can be used to shorten the simulation time [30–33].Virtual reality (VR) has the property of intuition and interactive, which makes up the defeat of computer simulation and makes the results of simulation take on effect of 3D stereo display [34–38]. This paper will adopt visual reality technology (VRT) to visualize the results of dust evolution. A virtual Shanghai University campus is rendered by OpenGL, and a parallel simulation system of dust evolution in the virtual campus is created in order to study the impact of different parameters (such as nonrain period and wind) on dust evolution. By comparison with experimental data, the validity of the model is verified. The research aims at providing theory and quantitative reference for dust particles evolution.


## 2. Kinetic Monte Carlo Simulation of Surface Dust Evolution

### 2.1. The 3D Modeling of Virtual Campus

Shanghai University campus is simplified and the 3D model is shown in Figure1.Figure 1
3D model of Shanghai University campus.The campus is 120 meter in width and 470 meter in length. Since most of surface dust works in the height of 0 to 5 meters, the virtual campus is in the scale of 120 m * 470 m * 6 m. To verify the validity of the KMC simulation method, five positions are chosen as the sampling point and signed as 1, 2, 3, 4, and 5, shown in Figure1. Location 1 is greatly influenced by the buildings, while the transportation has less impact on it. Location 2 is located on the campus road which is near the greensward. The flow of people has greater impact on it than the buildings. Location 3 and Location 4 are moderately influenced by the buildings and the flow of people. Location 5 is on a road where traffic is quite large and slightly impacted by buildings. By collecting and analyzing the weight of dust at the 5 positions, and comparing them with the KMC simulation result, the initial parameters of dust evolution simulation are defined to suit the real situation.


### 2.2. KMC Simulation Modeling of Surface Dust Particle in Virtual Campus

The surface dust modeling is dynamic, and the dust particles in the simulation system have three periods, including “emerging,” “moving,” and “vanishing.” With the passage of time, some of the existing particles are vanishing and the new ones are emerging. The survival particles move in the simulation process randomly, and three kinds of events occur in the process of movement: sediment, diffusion, and resuspension. Therefore, the KMC simulation modeling of dust particles consists of five events.(1)
Dust particles merge and join in the simulation system with attributes provided.
(2)
Dust particles diffuse in the air and update their attributes.
(3)
Certain dust particles sediment on the ground and update their attributes.
(4)
Some particles on the ground are resuspension and update their attributes.
(5)
Dust particles exceeding the life cycle are deleted from the simulation system.


#### 2.2.1. The Initialization of Dust Particles

It is assumed that 108 dust particles are released into the virtual campus environment. Every dust particle has six initial attributes: (1) initial position; (2) initial wind speed; (3) initial size; (4) initial color; (5) shape; (6) life cycle.The initial positions of dust particles are uniformly distributed. These particles have the same size, shape, and color. Their density is 1800 kg/m3. The initial speed of each particle is the wind speed of its position which can be obtained by numerical simulation of wind field in the virtual campus. Common Computing Fluid Dynamics (CFD) software, such as Fluent, CFX, and Phonics Star-CD, can simulate the wind field.Assuming air is incompressible viscous fluid, the type of flow is turbulent; density is regarded as constant. Control equation includes continuity and momentum equation. Related equations are as follows [39]:
(1)
∂
u
i
∂
x
i
=
0
,
∂
u
i
∂
t
+
∂
(
u
i
u
j
)
∂
x
j
=
-
1
ρ
∂
P
∂
x
i
+
∂
∂
x
j
(
v
∂
u
i
∂
x
j
)
;
x
i and x
j
(
i
,
j
=
1,2
,
3
) present the distance in each coordinate (x
,
y
,
z), u
i and u
j present the velocity component in each coordinate (x
,
y
,
z), P is air pressure, ρ is air density, and t is time.To solve the strongly swirling flow problem in the numerical simulation of wind field in the virtual campus,k-ε model is used to produce certain distortion:
(2)
k
=
3
2
(
u
ref
-
T
i
)
2
;
ε
=
C
μ
3
/
4
k
3
/
2
ℓ
;
ℓ
=
0.07
L
.
u
ref
- denotes inflow point velocity; L denotes equivalent length.Besides, in order to obtain the accurate simulation results of wind field, some improvements are made.(1) There are two types of grid during Shanghai University virtual campus space grid division, including triangle and hybrid grid. Triangle grid has strong boundary adaptability and hybrid grid can save compute time. In order to simply model, campus local area and plane grid is used [40].The simulation area contains two parts: building interval area and external wind field area, shown in Figure2. The grid has been divided into unstructured grid and structured grid separately. In order to improve the accuracy of calculation, the building interval area is encrypted. Especially, the external wind field uses quadrilateral structure grid and the internal wind field uses triangle unstructured grid.Figure 2
Grid division in virtual campus.(2) There are three main boundaries in the simulation of wind field: flow inlet boundary, flow outlet boundary, and solid boundary. In order to simplify and optimize the boundary conditions, specific setting is listed in Table1.Table 1
Boundary conditions.
Boundary
Condition
Formula
Sky
No slip wall boundary
u
=
v
=
0
Inlet condition
Velocity inlet boundary
ν
=
3.2 m/s
Ground, building
Free-slip boundary
v
=
0,   ∂
u
∂
y
=
0
Outlet condition
Pressure-type boundary
∂
P
∂
x
=
0,  Δ
ρ
=
235.56 pa(3) Before the simulation, it is necessary to make sure of its validity and judge its convergence. The experimental result shows that when iteration is around 1800 times, most kinds of iterative curves are close to our setting number under the residuals:1
e
-
06. As shown in Figure 3, the wind field simulation result in this condition has good convergency and fidelity.Figure 3
Residual curve simulated by Fluent.According to the meteorological record, the average wind speed 3.2 m/s is set and wind direction is northeast. The wind speed in virtual campus is shown in Figure4.Figure 4
Profile of wind speed.CFD postprocessing software, such as Ensight, Tecplot, and FieldView, can deal with the 3D grid data. Meanwhile, the spatial coordinates of dust particles are written on the text in a grid data form. With the help of CFD, wind speed at each grid node is available and the wind speed anywhere on campus can be drawn by interpolation.


#### 2.2.2. Dust Particles’ KMC Movement

Once the particles emerge, they begin to move when they get the initial attributes. Their move attributes can be derived by the initial ones, considering sediment, diffusion, and resuspension incidents, respectively.(
1) The Diffusion of Surface Dust Particles. Considering drag force, gravitational setting, Saffman lift force, and turbulent diffusions in the process of computation, motion equations of the particle can be written as [12]
(3)
d
u
p
d
t
=
F
D
(
u
-
u
p
)
+
g
x
(
ρ
p
-
ρ
)
ρ
p
.
F
D
(
u
-
u
p
) is the drag force per unit particle mass
(4)
F
D
=
18
μ
ρ
p
d
p
2
C
D
Re
24
.
u represents wind speed, u
p is particle speed, ρ is air density, ρ
p is the density of particles, d
p is particle diameter, μ is the molecular viscosity of the fluid, C
D is drag coefficient, and Re is Reynolds number.This event mainly deals with dust particles on airflow field, namely, particle moves on the MC lattice mentioned above. If particles move with no memory and have equal probability to each direction, the particle’s motion can be taken as Monte-Carlo motion [26].The trajectory of particles is shown as follows:(5)
X
i
=
v
t
t
+
X
0
i
(
i
=
1,2
,
3
)
.
i is the coordinate direction of x, y, z.The motion of a particle is defined as(6)
u
(
t
+
Δ
t
)
=
u
(
t
)
+
u
′
(
t
)
,
v
(
t
+
Δ
t
)
=
v
(
t
)
+
v
′
(
t
)
,
w
(
t
+
Δ
t
)
=
w
(
t
)
+
w
′
(
t
)
;
u
(
t
), v
(
t
), and w
(
t
) are the mean values in this time step.u
′
(
t
), v
′
(
t
), and w
′
(
t
) are compounded by two parts: relevant part and stochastic part
(7)
u
′
(
t
)
=
F
D
x
(
u
-
u
p
)
c
,
v
′
(
t
)
=
F
D
y
(
u
-
u
p
)
c
+
g
y
(
ρ
p
-
ρ
)
ρ
p
,
w
′
(
t
)
=
F
D
z
(
u
-
u
p
)
c
.
F
D
x
(
u
-
u
p
), F
D
y
(
u
-
u
p
), and F
D
z
(
u
-
u
p
) are the effectiveness of wind to the particles; c is the probability to any of the direction, c
=
v
/
v
max
⁡.So the particles diffuse in the virtual campus environment according to the formulas above.(
2) The Sediment of Surface Dust Particles. When the particles fulfill the following condition:
(8)
0
<
z
j
(
t
)
<
5
mm
,
v
′
(
t
)
<
0
they can be taken as sediment on the surface of the earth in this model. At this time, v
(
t
+
Δ
t
)
=
0.(
3) The Resuspension of Surface Dust Particles. Many researches based on Reynolds stress have taken that the particles may be resuspended when the wind speed is increasing to critical friction velocity. The equation can be shown as follows:
(9)
u
f
=
A
(
ρ
s
-
ρ
)
d
g
ρ
.
u
f is the critical friction velocity, ρ
s is particle density, d is the diameter of the particle, g is the acceleration of gravity, ρ is the density of wind, and A is a random coefficient with values between 0.16 and 0.21.When the particle is resuspended,(10)
z
j
(
t
)
=
5
mm
.And the particle velocity is equal to the wind speed.


#### 2.2.3. The Vanish of Surface Dust Particles

The particles have life cycle once they emerge in the virtual environment. They can be deemed as vanished when they move out of the boundary of the virtual campus or they move into the interior of the buildings.


## 2.1. The 3D Modeling of Virtual Campus

Shanghai University campus is simplified and the 3D model is shown in Figure1.Figure 1
3D model of Shanghai University campus.The campus is 120 meter in width and 470 meter in length. Since most of surface dust works in the height of 0 to 5 meters, the virtual campus is in the scale of 120 m * 470 m * 6 m. To verify the validity of the KMC simulation method, five positions are chosen as the sampling point and signed as 1, 2, 3, 4, and 5, shown in Figure1. Location 1 is greatly influenced by the buildings, while the transportation has less impact on it. Location 2 is located on the campus road which is near the greensward. The flow of people has greater impact on it than the buildings. Location 3 and Location 4 are moderately influenced by the buildings and the flow of people. Location 5 is on a road where traffic is quite large and slightly impacted by buildings. By collecting and analyzing the weight of dust at the 5 positions, and comparing them with the KMC simulation result, the initial parameters of dust evolution simulation are defined to suit the real situation.


## 2.2. KMC Simulation Modeling of Surface Dust Particle in Virtual Campus

The surface dust modeling is dynamic, and the dust particles in the simulation system have three periods, including “emerging,” “moving,” and “vanishing.” With the passage of time, some of the existing particles are vanishing and the new ones are emerging. The survival particles move in the simulation process randomly, and three kinds of events occur in the process of movement: sediment, diffusion, and resuspension. Therefore, the KMC simulation modeling of dust particles consists of five events.(1)
Dust particles merge and join in the simulation system with attributes provided.
(2)
Dust particles diffuse in the air and update their attributes.
(3)
Certain dust particles sediment on the ground and update their attributes.
(4)
Some particles on the ground are resuspension and update their attributes.
(5)
Dust particles exceeding the life cycle are deleted from the simulation system.


### 2.2.1. The Initialization of Dust Particles

It is assumed that 108 dust particles are released into the virtual campus environment. Every dust particle has six initial attributes: (1) initial position; (2) initial wind speed; (3) initial size; (4) initial color; (5) shape; (6) life cycle.The initial positions of dust particles are uniformly distributed. These particles have the same size, shape, and color. Their density is 1800 kg/m3. The initial speed of each particle is the wind speed of its position which can be obtained by numerical simulation of wind field in the virtual campus. Common Computing Fluid Dynamics (CFD) software, such as Fluent, CFX, and Phonics Star-CD, can simulate the wind field.Assuming air is incompressible viscous fluid, the type of flow is turbulent; density is regarded as constant. Control equation includes continuity and momentum equation. Related equations are as follows [39]:
(1)
∂
u
i
∂
x
i
=
0
,
∂
u
i
∂
t
+
∂
(
u
i
u
j
)
∂
x
j
=
-
1
ρ
∂
P
∂
x
i
+
∂
∂
x
j
(
v
∂
u
i
∂
x
j
)
;
x
i and x
j
(
i
,
j
=
1,2
,
3
) present the distance in each coordinate (x
,
y
,
z), u
i and u
j present the velocity component in each coordinate (x
,
y
,
z), P is air pressure, ρ is air density, and t is time.To solve the strongly swirling flow problem in the numerical simulation of wind field in the virtual campus,k-ε model is used to produce certain distortion:
(2)
k
=
3
2
(
u
ref
-
T
i
)
2
;
ε
=
C
μ
3
/
4
k
3
/
2
ℓ
;
ℓ
=
0.07
L
.
u
ref
- denotes inflow point velocity; L denotes equivalent length.Besides, in order to obtain the accurate simulation results of wind field, some improvements are made.(1) There are two types of grid during Shanghai University virtual campus space grid division, including triangle and hybrid grid. Triangle grid has strong boundary adaptability and hybrid grid can save compute time. In order to simply model, campus local area and plane grid is used [40].The simulation area contains two parts: building interval area and external wind field area, shown in Figure2. The grid has been divided into unstructured grid and structured grid separately. In order to improve the accuracy of calculation, the building interval area is encrypted. Especially, the external wind field uses quadrilateral structure grid and the internal wind field uses triangle unstructured grid.Figure 2
Grid division in virtual campus.(2) There are three main boundaries in the simulation of wind field: flow inlet boundary, flow outlet boundary, and solid boundary. In order to simplify and optimize the boundary conditions, specific setting is listed in Table1.Table 1
Boundary conditions.
Boundary
Condition
Formula
Sky
No slip wall boundary
u
=
v
=
0
Inlet condition
Velocity inlet boundary
ν
=
3.2 m/s
Ground, building
Free-slip boundary
v
=
0,   ∂
u
∂
y
=
0
Outlet condition
Pressure-type boundary
∂
P
∂
x
=
0,  Δ
ρ
=
235.56 pa(3) Before the simulation, it is necessary to make sure of its validity and judge its convergence. The experimental result shows that when iteration is around 1800 times, most kinds of iterative curves are close to our setting number under the residuals:1
e
-
06. As shown in Figure 3, the wind field simulation result in this condition has good convergency and fidelity.Figure 3
Residual curve simulated by Fluent.According to the meteorological record, the average wind speed 3.2 m/s is set and wind direction is northeast. The wind speed in virtual campus is shown in Figure4.Figure 4
Profile of wind speed.CFD postprocessing software, such as Ensight, Tecplot, and FieldView, can deal with the 3D grid data. Meanwhile, the spatial coordinates of dust particles are written on the text in a grid data form. With the help of CFD, wind speed at each grid node is available and the wind speed anywhere on campus can be drawn by interpolation.


### 2.2.2. Dust Particles’ KMC Movement

Once the particles emerge, they begin to move when they get the initial attributes. Their move attributes can be derived by the initial ones, considering sediment, diffusion, and resuspension incidents, respectively.(
1) The Diffusion of Surface Dust Particles. Considering drag force, gravitational setting, Saffman lift force, and turbulent diffusions in the process of computation, motion equations of the particle can be written as [12]
(3)
d
u
p
d
t
=
F
D
(
u
-
u
p
)
+
g
x
(
ρ
p
-
ρ
)
ρ
p
.
F
D
(
u
-
u
p
) is the drag force per unit particle mass
(4)
F
D
=
18
μ
ρ
p
d
p
2
C
D
Re
24
.
u represents wind speed, u
p is particle speed, ρ is air density, ρ
p is the density of particles, d
p is particle diameter, μ is the molecular viscosity of the fluid, C
D is drag coefficient, and Re is Reynolds number.This event mainly deals with dust particles on airflow field, namely, particle moves on the MC lattice mentioned above. If particles move with no memory and have equal probability to each direction, the particle’s motion can be taken as Monte-Carlo motion [26].The trajectory of particles is shown as follows:(5)
X
i
=
v
t
t
+
X
0
i
(
i
=
1,2
,
3
)
.
i is the coordinate direction of x, y, z.The motion of a particle is defined as(6)
u
(
t
+
Δ
t
)
=
u
(
t
)
+
u
′
(
t
)
,
v
(
t
+
Δ
t
)
=
v
(
t
)
+
v
′
(
t
)
,
w
(
t
+
Δ
t
)
=
w
(
t
)
+
w
′
(
t
)
;
u
(
t
), v
(
t
), and w
(
t
) are the mean values in this time step.u
′
(
t
), v
′
(
t
), and w
′
(
t
) are compounded by two parts: relevant part and stochastic part
(7)
u
′
(
t
)
=
F
D
x
(
u
-
u
p
)
c
,
v
′
(
t
)
=
F
D
y
(
u
-
u
p
)
c
+
g
y
(
ρ
p
-
ρ
)
ρ
p
,
w
′
(
t
)
=
F
D
z
(
u
-
u
p
)
c
.
F
D
x
(
u
-
u
p
), F
D
y
(
u
-
u
p
), and F
D
z
(
u
-
u
p
) are the effectiveness of wind to the particles; c is the probability to any of the direction, c
=
v
/
v
max
⁡.So the particles diffuse in the virtual campus environment according to the formulas above.(
2) The Sediment of Surface Dust Particles. When the particles fulfill the following condition:
(8)
0
<
z
j
(
t
)
<
5
mm
,
v
′
(
t
)
<
0
they can be taken as sediment on the surface of the earth in this model. At this time, v
(
t
+
Δ
t
)
=
0.(
3) The Resuspension of Surface Dust Particles. Many researches based on Reynolds stress have taken that the particles may be resuspended when the wind speed is increasing to critical friction velocity. The equation can be shown as follows:
(9)
u
f
=
A
(
ρ
s
-
ρ
)
d
g
ρ
.
u
f is the critical friction velocity, ρ
s is particle density, d is the diameter of the particle, g is the acceleration of gravity, ρ is the density of wind, and A is a random coefficient with values between 0.16 and 0.21.When the particle is resuspended,(10)
z
j
(
t
)
=
5
mm
.And the particle velocity is equal to the wind speed.


### 2.2.3. The Vanish of Surface Dust Particles

The particles have life cycle once they emerge in the virtual environment. They can be deemed as vanished when they move out of the boundary of the virtual campus or they move into the interior of the buildings.


## 2.2.1. The Initialization of Dust Particles

It is assumed that 108 dust particles are released into the virtual campus environment. Every dust particle has six initial attributes: (1) initial position; (2) initial wind speed; (3) initial size; (4) initial color; (5) shape; (6) life cycle.The initial positions of dust particles are uniformly distributed. These particles have the same size, shape, and color. Their density is 1800 kg/m3. The initial speed of each particle is the wind speed of its position which can be obtained by numerical simulation of wind field in the virtual campus. Common Computing Fluid Dynamics (CFD) software, such as Fluent, CFX, and Phonics Star-CD, can simulate the wind field.Assuming air is incompressible viscous fluid, the type of flow is turbulent; density is regarded as constant. Control equation includes continuity and momentum equation. Related equations are as follows [39]:
(1)
∂
u
i
∂
x
i
=
0
,
∂
u
i
∂
t
+
∂
(
u
i
u
j
)
∂
x
j
=
-
1
ρ
∂
P
∂
x
i
+
∂
∂
x
j
(
v
∂
u
i
∂
x
j
)
;
x
i and x
j
(
i
,
j
=
1,2
,
3
) present the distance in each coordinate (x
,
y
,
z), u
i and u
j present the velocity component in each coordinate (x
,
y
,
z), P is air pressure, ρ is air density, and t is time.To solve the strongly swirling flow problem in the numerical simulation of wind field in the virtual campus,k-ε model is used to produce certain distortion:
(2)
k
=
3
2
(
u
ref
-
T
i
)
2
;
ε
=
C
μ
3
/
4
k
3
/
2
ℓ
;
ℓ
=
0.07
L
.
u
ref
- denotes inflow point velocity; L denotes equivalent length.Besides, in order to obtain the accurate simulation results of wind field, some improvements are made.(1) There are two types of grid during Shanghai University virtual campus space grid division, including triangle and hybrid grid. Triangle grid has strong boundary adaptability and hybrid grid can save compute time. In order to simply model, campus local area and plane grid is used [40].The simulation area contains two parts: building interval area and external wind field area, shown in Figure2. The grid has been divided into unstructured grid and structured grid separately. In order to improve the accuracy of calculation, the building interval area is encrypted. Especially, the external wind field uses quadrilateral structure grid and the internal wind field uses triangle unstructured grid.Figure 2
Grid division in virtual campus.(2) There are three main boundaries in the simulation of wind field: flow inlet boundary, flow outlet boundary, and solid boundary. In order to simplify and optimize the boundary conditions, specific setting is listed in Table1.Table 1
Boundary conditions.
Boundary
Condition
Formula
Sky
No slip wall boundary
u
=
v
=
0
Inlet condition
Velocity inlet boundary
ν
=
3.2 m/s
Ground, building
Free-slip boundary
v
=
0,   ∂
u
∂
y
=
0
Outlet condition
Pressure-type boundary
∂
P
∂
x
=
0,  Δ
ρ
=
235.56 pa(3) Before the simulation, it is necessary to make sure of its validity and judge its convergence. The experimental result shows that when iteration is around 1800 times, most kinds of iterative curves are close to our setting number under the residuals:1
e
-
06. As shown in Figure 3, the wind field simulation result in this condition has good convergency and fidelity.Figure 3
Residual curve simulated by Fluent.According to the meteorological record, the average wind speed 3.2 m/s is set and wind direction is northeast. The wind speed in virtual campus is shown in Figure4.Figure 4
Profile of wind speed.CFD postprocessing software, such as Ensight, Tecplot, and FieldView, can deal with the 3D grid data. Meanwhile, the spatial coordinates of dust particles are written on the text in a grid data form. With the help of CFD, wind speed at each grid node is available and the wind speed anywhere on campus can be drawn by interpolation.


## 2.2.2. Dust Particles’ KMC Movement

Once the particles emerge, they begin to move when they get the initial attributes. Their move attributes can be derived by the initial ones, considering sediment, diffusion, and resuspension incidents, respectively.(
1) The Diffusion of Surface Dust Particles. Considering drag force, gravitational setting, Saffman lift force, and turbulent diffusions in the process of computation, motion equations of the particle can be written as [12]
(3)
d
u
p
d
t
=
F
D
(
u
-
u
p
)
+
g
x
(
ρ
p
-
ρ
)
ρ
p
.
F
D
(
u
-
u
p
) is the drag force per unit particle mass
(4)
F
D
=
18
μ
ρ
p
d
p
2
C
D
Re
24
.
u represents wind speed, u
p is particle speed, ρ is air density, ρ
p is the density of particles, d
p is particle diameter, μ is the molecular viscosity of the fluid, C
D is drag coefficient, and Re is Reynolds number.This event mainly deals with dust particles on airflow field, namely, particle moves on the MC lattice mentioned above. If particles move with no memory and have equal probability to each direction, the particle’s motion can be taken as Monte-Carlo motion [26].The trajectory of particles is shown as follows:(5)
X
i
=
v
t
t
+
X
0
i
(
i
=
1,2
,
3
)
.
i is the coordinate direction of x, y, z.The motion of a particle is defined as(6)
u
(
t
+
Δ
t
)
=
u
(
t
)
+
u
′
(
t
)
,
v
(
t
+
Δ
t
)
=
v
(
t
)
+
v
′
(
t
)
,
w
(
t
+
Δ
t
)
=
w
(
t
)
+
w
′
(
t
)
;
u
(
t
), v
(
t
), and w
(
t
) are the mean values in this time step.u
′
(
t
), v
′
(
t
), and w
′
(
t
) are compounded by two parts: relevant part and stochastic part
(7)
u
′
(
t
)
=
F
D
x
(
u
-
u
p
)
c
,
v
′
(
t
)
=
F
D
y
(
u
-
u
p
)
c
+
g
y
(
ρ
p
-
ρ
)
ρ
p
,
w
′
(
t
)
=
F
D
z
(
u
-
u
p
)
c
.
F
D
x
(
u
-
u
p
), F
D
y
(
u
-
u
p
), and F
D
z
(
u
-
u
p
) are the effectiveness of wind to the particles; c is the probability to any of the direction, c
=
v
/
v
max
⁡.So the particles diffuse in the virtual campus environment according to the formulas above.(
2) The Sediment of Surface Dust Particles. When the particles fulfill the following condition:
(8)
0
<
z
j
(
t
)
<
5
mm
,
v
′
(
t
)
<
0
they can be taken as sediment on the surface of the earth in this model. At this time, v
(
t
+
Δ
t
)
=
0.(
3) The Resuspension of Surface Dust Particles. Many researches based on Reynolds stress have taken that the particles may be resuspended when the wind speed is increasing to critical friction velocity. The equation can be shown as follows:
(9)
u
f
=
A
(
ρ
s
-
ρ
)
d
g
ρ
.
u
f is the critical friction velocity, ρ
s is particle density, d is the diameter of the particle, g is the acceleration of gravity, ρ is the density of wind, and A is a random coefficient with values between 0.16 and 0.21.When the particle is resuspended,(10)
z
j
(
t
)
=
5
mm
.And the particle velocity is equal to the wind speed.


## 2.2.3. The Vanish of Surface Dust Particles

The particles have life cycle once they emerge in the virtual environment. They can be deemed as vanished when they move out of the boundary of the virtual campus or they move into the interior of the buildings.


## 3. KMC-Based Parallel Simulation of Dust Evolution

The premise of KMC-based simulation of dust evolution is that the surface dust location can be described by a point in the virtual campus. Thus, it is hard for a computer to complete the simulation task with the increase of virtual environment or dust particles. At this time, the virtual environment can be divided into several small subspaces and particles in the subspaces are assigned to multiprocessors to simulate concurrently.


### 3.1. The Way of Data Distribution

Since the evolution of surface dust particles is a random process, the virtual campus space should be divided into continuous space and every space contains the same number of dust particles. Each process simulates the evolution of dust particles in a subspace, which can effectively ensure load balancing of each process [41–43].When virtual environment is divided into subspaces, the dust particles in regional boundary of a subspace may move to another subspace, which causes the communication between processes. When the particles are in the regional boundary, two adjacent regions need to communicate and confirm where particles are and their specific locations. In order to divide the virtual campus, block data distribution can be achieved by two ways: lattice and sheet divisions, as shown in Figure5.KMC division of virtual campus.
(a)
Lattice division
(b)
Sheet divisionThe division in Figure5(a) makes each region be required to communicate with at least three adjacent regions which have the same boundary, while the division in Figure 5(b) makes each region be required to communicate with at most two adjacent regions which have the same boundary. So the division in Figure 5(b) will reduce the traffic and it is used in the division of the virtual campus.


### 3.2. Communication Optimization Strategy

Particles’ sediment, diffusion, and resuspension should be taken into account when the dust evolution based on KMC is used. To optimize the communication strategy and reduce the traffic between processors, the data storage space is divided into 2 parts: Local Store and Neighbor Copy. Local Store keeps the data of particles in the local subspace, while Neighbor Copy keeps the data of particles in the regional boundary of other neighboring subspace.When the sediment and resuspension processes are simulated, the attributes of particles in Local Store need renewing and the particles which meet the condition sediment on the ground or re-suspend in the air. So the processors do not need to communicate with each other and reduce the communication frequency.In the diffusion process, particles in the regional boundary update their attributes both in Local Store and Neighbor Copy. However, particles which are not in the regional boundary only update their attributes in Local Store, which will reduce the communication frequency and traffic between processors.According to the description above, each process not only stores the data of particles in the local subspace, but also stores the data of particles in the neighborhood space. In Figure6, the virtual environment is divided into three processors and dashed areas contain the Neighbor Copy space of process 2 because it needs to communicate with processes 1 and 3 in the diffusion process.Figure 6
Dust data storage model in three processes.Data structure of particle in the communication is described as follows:(11)
A
particle
=
(
pos
,
vel
,
size
,
color
,
shape
,
lifecyle
)
,
where pos and vel are the particle’s location and velocity in the virtual campus at this time; size, color, shape, and lifecycle of each particle in the Neighbor Copy have the same value.


### 3.3. KMC-Based Dust Evolution Parallel Simulation Algorithm

Algorithm1 is the frame of KMC-based dust evolution parallel simulation algorithm. Three main processes of dust evolution are shown in Algorithms 2, 3, and 4.Algorithm 1:Frame of KMC-based dust evolution parallel simulation algorithm.
(1) Begin
(2) Initialization MPI, Define the number of processes and simulation timet;
(3) Master process reads the data structure of dust particles;
(4) Master process distributes the data to Local_Store of each sub-processors according to
the data distribution way of sheet division, and the simulation timeri
=
0;
(5) ifi
<
=
t
,go to (6); otherwise, go to (10);
(6) Each sub-processor executes sediment of dust particles;
(7) Each sub-processor executes re-suspension of dust particles;
(8) Each sub-processor executes diffusion of dust particles;
(9) Each sub-processor updates the value of simulation timeri, go to (5);
(10) Each sub-processor sends data structure of dust particles to the master process;
(11) The master processor collects the data sent by each sub-processor and writes to the file Output;
(12) EndAlgorithm 2:The sediment simulation algorithm of dust particles.
(1) Begin
(2) Each processor updates the position of each particle in the Local_Store;
(3) IfA
particle
(
j
)
·
pos
·
z
<
5
mm
&
&
A
particle
(
j
)
·
vel
·
z
<
0, go to (4); Otherwise, go to (5);
(4) Particlej is sediment on the ground, A
particle
(
j
)
·
pos
·
z
=
0;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) EndAlgorithm 3:The re-suspension simulation algorithm of dust particles.
(1) Begin
(2) Each processor updates the position of each particle in the Local_Store;
(3) IfA
particle
(
j
)
·
pos
·
z
=
0, go to (5); Otherwise go to (5);
(4) If the wind speed in the location of particlej is greater than the critical friction velocity,
Particlej is re-suspended and A
particle
(
j
)
·
pos
·
z
=
5
mm;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) EndAlgorithm 4:The diffusion simulation algorithm of dust particles.
(1) Begin
(2) Each processor sends the update data of particles in the border to other processors, then
receives the update data of particles from other processors and copies them in Neighbor_Copy;
(3) Calculate the diffusion probability and diffusion direction of the particlej in Local_Store;
(4) Execute the diffusion of particles, updateA
particle
(
j
)
·
pos, mark particles which enter other Neighbor_Copy;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) Communicate with other processor to update the data of particles in Neighbor_Copy;
(7) End


## 3.1. The Way of Data Distribution

Since the evolution of surface dust particles is a random process, the virtual campus space should be divided into continuous space and every space contains the same number of dust particles. Each process simulates the evolution of dust particles in a subspace, which can effectively ensure load balancing of each process [41–43].When virtual environment is divided into subspaces, the dust particles in regional boundary of a subspace may move to another subspace, which causes the communication between processes. When the particles are in the regional boundary, two adjacent regions need to communicate and confirm where particles are and their specific locations. In order to divide the virtual campus, block data distribution can be achieved by two ways: lattice and sheet divisions, as shown in Figure5.KMC division of virtual campus.
(a)
Lattice division
(b)
Sheet divisionThe division in Figure5(a) makes each region be required to communicate with at least three adjacent regions which have the same boundary, while the division in Figure 5(b) makes each region be required to communicate with at most two adjacent regions which have the same boundary. So the division in Figure 5(b) will reduce the traffic and it is used in the division of the virtual campus.


## 3.2. Communication Optimization Strategy

Particles’ sediment, diffusion, and resuspension should be taken into account when the dust evolution based on KMC is used. To optimize the communication strategy and reduce the traffic between processors, the data storage space is divided into 2 parts: Local Store and Neighbor Copy. Local Store keeps the data of particles in the local subspace, while Neighbor Copy keeps the data of particles in the regional boundary of other neighboring subspace.When the sediment and resuspension processes are simulated, the attributes of particles in Local Store need renewing and the particles which meet the condition sediment on the ground or re-suspend in the air. So the processors do not need to communicate with each other and reduce the communication frequency.In the diffusion process, particles in the regional boundary update their attributes both in Local Store and Neighbor Copy. However, particles which are not in the regional boundary only update their attributes in Local Store, which will reduce the communication frequency and traffic between processors.According to the description above, each process not only stores the data of particles in the local subspace, but also stores the data of particles in the neighborhood space. In Figure6, the virtual environment is divided into three processors and dashed areas contain the Neighbor Copy space of process 2 because it needs to communicate with processes 1 and 3 in the diffusion process.Figure 6
Dust data storage model in three processes.Data structure of particle in the communication is described as follows:(11)
A
particle
=
(
pos
,
vel
,
size
,
color
,
shape
,
lifecyle
)
,
where pos and vel are the particle’s location and velocity in the virtual campus at this time; size, color, shape, and lifecycle of each particle in the Neighbor Copy have the same value.


## 3.3. KMC-Based Dust Evolution Parallel Simulation Algorithm

Algorithm1 is the frame of KMC-based dust evolution parallel simulation algorithm. Three main processes of dust evolution are shown in Algorithms 2, 3, and 4.Algorithm 1:Frame of KMC-based dust evolution parallel simulation algorithm.
(1) Begin
(2) Initialization MPI, Define the number of processes and simulation timet;
(3) Master process reads the data structure of dust particles;
(4) Master process distributes the data to Local_Store of each sub-processors according to
the data distribution way of sheet division, and the simulation timeri
=
0;
(5) ifi
<
=
t
,go to (6); otherwise, go to (10);
(6) Each sub-processor executes sediment of dust particles;
(7) Each sub-processor executes re-suspension of dust particles;
(8) Each sub-processor executes diffusion of dust particles;
(9) Each sub-processor updates the value of simulation timeri, go to (5);
(10) Each sub-processor sends data structure of dust particles to the master process;
(11) The master processor collects the data sent by each sub-processor and writes to the file Output;
(12) EndAlgorithm 2:The sediment simulation algorithm of dust particles.
(1) Begin
(2) Each processor updates the position of each particle in the Local_Store;
(3) IfA
particle
(
j
)
·
pos
·
z
<
5
mm
&
&
A
particle
(
j
)
·
vel
·
z
<
0, go to (4); Otherwise, go to (5);
(4) Particlej is sediment on the ground, A
particle
(
j
)
·
pos
·
z
=
0;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) EndAlgorithm 3:The re-suspension simulation algorithm of dust particles.
(1) Begin
(2) Each processor updates the position of each particle in the Local_Store;
(3) IfA
particle
(
j
)
·
pos
·
z
=
0, go to (5); Otherwise go to (5);
(4) If the wind speed in the location of particlej is greater than the critical friction velocity,
Particlej is re-suspended and A
particle
(
j
)
·
pos
·
z
=
5
mm;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) EndAlgorithm 4:The diffusion simulation algorithm of dust particles.
(1) Begin
(2) Each processor sends the update data of particles in the border to other processors, then
receives the update data of particles from other processors and copies them in Neighbor_Copy;
(3) Calculate the diffusion probability and diffusion direction of the particlej in Local_Store;
(4) Execute the diffusion of particles, updateA
particle
(
j
)
·
pos, mark particles which enter other Neighbor_Copy;
(5) If all particles in Local_Store are finished scanning, go to (6); Otherwise, go to (3);
(6) Communicate with other processor to update the data of particles in Neighbor_Copy;
(7) End


## 4. Visualization on Surface Dust Evolution

The coordinates of particles obtained from the concurrent calculation are recorded in the text, and OpenGL makes the dust evolution visible. To see the dust evolution clearly, the boundaries of virtual campus are drawn as transparent.In the visualization on surface dust evolution, the dust particles are drawn pro rata because of the large amount, and the coefficient scaleK is 10−3:
(12)
N
P
a
r
t
s
(
t
)
=
S
u
m
D
N
(
t
)
*
K
.Suppose that the total particle amount in output att moment is S
u
m
D
N
(
t
). Then the visible amount of dust particles is N
P
a
r
t
s
(
t
) at the simulation of t moment.The visualization algorithm of dust particles’ KMC evolution is shown in Algorithm5. And Figure 7 shows the visualization result. It is clear to see the evolution of dust particles in the virtual campus.Algorithm 5:The visualization algorithm of dust particles’ KMC evolution.
(1) Begin
(2) Create a transparent border of virtual campus;
(3) Render a yellow ground;
(4) Render the buildings in the virtual campus;
(5) Set the attributes of dust particles, including their size, color and shape;
(6) Read the coordinates of dust particles from file Output, calculate NParts;
(7) If the particles are in the area of sampling pointsL
i,
calculate the accumulation amount of dust particlesS
(
L
i
);
(8) Render the particles according to the value of NParts;
(9) Output the accumulation amount of dust particles in five sampling points;
(10) EndFigure 7
The visualization result of dust evolution in virtual campus.


## 5. Results and Analysis

In order to evaluate the effectiveness of KMC-based parallel simulation algorithm of dust evolution in virtual campus environment, the experiment is designed as follows.Dust in five collection areas of campus is collected each nonrain day during four months. The weight of dust is gained by a delicate electronic balance and recorded. At the same time, the weather condition like wind scale, rainy day, and nonrain period is marked. According to the records, the northeast wind is the most frequent wind during the experimental period. So the following analysis is based on the condition of northeast wind. The record of dust accumulation in the northeast wind is shown in Table2.Table 2
The record of dust accumulation in different locations.
Nonrain period (day)
Location 1
Location 2
Location 3
Location 4
Location 5
1
0.406
1.095
2.095
0.745
0.171
1
0.037
1.378
0.782
0.215
0.663
1
0.098
1.790
3.964
1.846
0.285
1
0.347
2.524
1.208
2.558
0.705
1
0.166
1.251
0.920
1.974
0.520
1
0.234
2.372
1.272
2.705
0.984
1
0.336
1.128
2.449
1.612
2.056
1
0.282
0.877
0.920
0.614
0.551
1
0.112
0.654
0.346
1.115
1.554
1
0.193
0.511
0.628
1.016
2.063
1
0.205
0.960
0.922
0.838
0.924
2
0.056
1.226
2.198
0.985
0.221
2
0.296
1.865
1.088
1.084
0.547
2
0.294
1.384
0.887
2.015
0.552
2
0.259
0.597
0.611
1.224
1.866
2
0.180
0.599
0.528
1.054
1.586
2
0.181
0.497
0.314
0.215
1.470
3
0.264
1.426
1.149
0.446
0.623
3
0.226
1.495
0.808
0.791
0.753
3
0.184
1.064
0.839
0.439
3.740
4
0.045
1.804
2.656
1.575
0.219
4
0.055
1.452
0.484
0.895
1.098
4
0.214
3.577
1.203
1.388
0.315
4
0.214
0.833
0.804
1.034
1.401
4
0.195
0.804
0.510
0.964
0.854
4
0.073
1.506
1.565
0.921
1.286
4
0.156
0.691
1.476
1.366
1.481
5
0.082
1.533
1.194
0.271
0.148
5
0.323
0.602
1.958
0.978
0.694
5
0.194
0.399
1.430
1.138
2.630
6
0.049
2.416
2.467
1.561
0.453
6
0.116
0.385
0.348
0.859
1.806
7
0.080
0.804
1.070
0.999
0.575
7
0.150
0.482
0.602
0.471
1.404Figure8 compares the experimental results and simulation results based on KMC serial and parallel algorithm of dust evolution by the effect of different nonrain periods in five collection areas.Comparison of experimental and simulation results of dust fall accumulation.
(a)
Relationship of nonrain period and dust fall at Location 1
(b)
Relationship of nonrain period and dust fall at Location 2
(c)
Relationship of nonrain period and dust fall at Location 3
(d)
Relationship of nonrain period and dust fall at Location 4
(e)
Relationship of nonrain period and dust fall at Location 5Figures8(a) to 8(e) show that the dust fall accumulation became heavier and heavier with the increase of the nonrain period, which proves the effectiveness of simulation algorithm. From Figure 8, serial and parallel simulations had the same accumulation amount of dust particles in five points, which shows the accuracy of parallel simulation algorithm.To evaluate the validity of parallel simulation algorithm, its acceleration ratio and efficiency are calculated and results are shown in Table3.Table 3
Comparison study of results of parallel calculation in different processors.
The number of processes
S
P
T
P/hour
1
1
26.31
2
1.49
17.66
4
1.91
13.77
8
2.87
9.17
16
6.08
4.33The parallel acceleration ratio is defined as(13)
S
P
=
T
S
T
P
,
where T
S is the time used by serial algorithm and T
P is the time used by parallel algorithm in P processes.From Table2, the value of the acceleration ratio is small because the algorithm is related to the text operation, while the acceleration ratio increases with the number of the processes and the computation time reduces evidently. This indicates that the parallel algorithm on dust evolution can promote the efficiency although KMC evolution algorithm needs lots of boundary exams on the particles, which makes the simulation of large-scale virtual environment possible.


## 6. Conclusion

It is efficient to use the parallel algorithm to simulate the KMC evolution of surface dust particles in large-scale virtual environment. A parallel simulation algorithm of particles’ KMC evolution is proposed. It is useful to balance the load of every process and reduce the communication expense among processes with the help of data distribution way of sheet division and communication optimizing strategy. The experiment results show that simulation operation time is shortened enormously, the acceleration ratio is easy to get, and the parallel efficiency is promoted due to the reasonable process numbers in the parallel simulation algorithm which also compensates the disability of single computer. With the 3D visible simulation result, researchers can have a good understanding of the segmentation, diffusion, and resuspension of dust particles and analyze their movement discipline to lay a theoretical foundation for the dust prevention.



---

*Source: 839726-2014-01-21.xml*