# Emergency Evacuation of Hazardous Chemical Accidents Based on Diffusion Simulation

Authors: Jiang-Hua Zhang; Hai-Yue Liu; Rui Zhu; Yang Liu
Journal: Complexity
 (2017)

---

## Abstract

The recent rapid development of information technology, such as sensing technology, communications technology, and database, allows us to use simulation experiments for analyzing serious accidents caused by hazardous chemicals. Due to the toxicity and diffusion of hazardous chemicals, these accidents often lead to not only severe consequences and economic losses, but also traffic jams at the same time. Emergency evacuation after hazardous chemical accidents is an effective means to reduce the loss of life and property and to smoothly resume the transport network as soon as possible. This paper considers the dynamic changes of the hazardous chemicals’ concentration after their leakage and simulates the diffusion process. Based on the characteristics of emergency evacuation of hazardous chemical accidents, we build a mixed-integer programming model and design a heuristic algorithm using network optimization and diffusion simulation (hereafter NODS). We then verify the validity and feasibility of the algorithm using Jinan, China, as a computational example. In the end, we compare the results from different scenarios to explore the key factors affecting the effectiveness of the evacuation process.


---

## Body

## 1. Introduction

With the development of the world economy, hazardous chemicals have become important materials in industrial and agricultural production, national defense construction, and the daily life of people. However, in recent years, serious accidents caused by hazardous chemicals have occurred with increasing frequency in the world. For example, on February 7, 2010, an acrylate leakage accident occurred at a chemical company in Panyu, Guangzhou, and thousands of nearby villagers were evacuated urgently. On April 17, 2013, a fertilizer plant in the American region of western Texas exploded, setting at least 10 buildings on fire, causing 35 casualties and injuring more than 160 people. On July 19, 2014, an especially serious road traffic and deflagration accident occurred on the Shanghai-Kunming expressway in Shaoyang, Hunan, causing 54 casualties, injuring 6 people, and creating a direct economic loss of 53 million yuan. On November 15, 2014, a methanethiol leakage accident occurred at an American DuPont factory, which was located in the southeastern Houston area, causing 4 casualties, injuring 1 person, and affecting many surrounding residents. On August 12, 2015, huge explosions caused by improper storage of chemicals occurred in Tianjin and claimed the lives of 165 people and caused a directed economic loss of 6.8 billion yuan. Due to the toxicity and diffusion of hazardous chemicals, these accidents not only caused casualties and economic losses but also polluted the environment, causing serious social panic and traffic jams. Emergency evacuation after hazardous chemical accidents is an effective means to reduce the harm. In addition, the information technology is developing rapidly in recent years, so we can use these technologies to do some simulation experiences which can have a positive impact on emergency evacuation. Therefore, using advanced technologies to study the emergency evacuation following hazardous chemical accidents is very necessary.Research on the emergency evacuation problem first appeared in the 1960s. Evacuation models can be mainly divided into two categories: mathematical models based on network optimization and simulation models. Hamacher and Tjandra (2002) [1] reviewed the mathematical models based on network optimization. Murray-Tuite and Wolshon (2013) [2] presented a review of highway-based evacuation modeling and simulation and its evolution over the past decade. In recent years, many scholars have studied the emergency evacuation problem using different methods. Cova and Johnson (2003) [3], Xie and Turnquist (2011) [4], and Fang et al. (2013) [5] studied the lane-based evacuation network optimization problem. Lim et al. (2012) [6], Hadas and Laor (2013) [7], and Sheu and Pan (2014) [8] focused on the development of models for the design of an optimal network for emergency evacuation. Duanmu et al. (2012) [9], Fang et al. (2013) [5], Uchida (2014) [10], and Zhang and Chang (2014) [11] established evacuation models considering travel demand distribution problems and risk-averse drivers’ route choices. In addition, some scholars considered modeling evacuation networks based on the vulnerability of path networks, for example, Snelder et al. (2012) [12], Chen et al. (2012) [13], and Bell et al. (2014) [14]. And GIS-based DSS was developed to bring the emergency evacuation planning into implementation by some scholars, for example, Zhang et al. (2016) [15].Regarding the study of simulation models, Yin et al. (2014) [16] presented an agent-based travel demand model for simulating hurricane evacuations, which was capable of generating comprehensive household activity and travel plans. The system considered six typical evacuation decisions: evacuate/stay, accommodation type choice, evacuation destination choice, mode choice, vehicle usage choice, and departure time choice. Zhang et al. (2012) [17] introduced a data simulation model to express the topological and geometrical relationships of a building’s inner structures to determine the optimal routes for rescue and evacuation in a complex building. Nishino et al. (2012) [18] introduced a proposed simulation method that calculated the probability that the ratio of burned-down buildings or fire fatalities in a district would exceed a threshold within a given time period after an earthquake. Alvear et al. (2014) [19] used computer modeling and simulation to explore the impact that crew procedures had on evacuating two high-speed trains under different fire scenarios. Wu and Huang (2015) [20] simulated the dynamics of the evacuees and derived the evacuation times of a high-rise building by using the control volume model. C. Wang and J. Wang (2017) [21] presented an extended floor-field (FF) model combined with risk factors for emergency evacuation and they validated model through numerical simulations with specific room structures.Many scholars have conducted very fruitful studies on the building evacuation problem from different angles. Pursals and Garzón (2009) [22], Zhang et al. (2013) [23], and Ronchi et al. (2014) [24] formulated the building evacuation problem and generalized the evacuation function to incorporate evacuation routes using Monte Carlo and other stochastic methods. Wu and Chen (2012) [25] used a 3D geometric network model (GNM) and the Dijkstra algorithm to consider smoke movement during different times of a building fire. Bhushan and Sarda (2013) [26] introduced ladder exits in the building evacuation problem that were created by placing ladders at possible openings in a building. They modified the existing network flow based on an evacuation planning model to incorporate ladders and presented three different optimization formulations to model various scenarios.Siddiqui et al. (2012) [27] and Verma et al. (2012) [28] studied the evacuation problem of hazardous materials, in which Siddiqui et al. (2012) [27] proposed a computational fluid dynamics- (CFD-) based model for indoor risk assessment that considered accidental release of a sustained, small, undetected leak of a dense toxic gas (chlorine) in an industrial indoor environment. Li et al. (2012) [29], Najafi et al. (2013) [30], and Bish and Sherali (2013) [31] studied large-scale regional evacuations caused by hurricanes, wildfires, or earthquakes. They studied the demand-based strategies of aggregate-level staging and routing to structure the evacuation demand, both with and without congestion and then designed two heuristics. Shen et al. (2015) [32] established a conceptual decision support system for toxic leakage in a chemical industrial park and reported a method of determining the scope of evacuation following a leakage accident.As seen from above, the existing literature on emergency evacuation mostly focuses on disasters such as fires and hurricanes. However, there are few studies of emergency evacuations following hazardous chemical accidents. As everyone knows, disasters and accidents of different types have different characteristics. Thus, emergency evacuation plans should be made according to the different characteristics of specific issues. Hazardous chemical accidents are different from other accidents because the toxicity and the diffusion of chemicals make the hazardous chemicals’ concentration change with the passage of time, weather conditions, and wind speed. Therefore, it is very necessary to research the emergency evacuation of hazardous chemical accidents.In this paper, we focus on the emergency evacuation following hazardous chemical accidents. We build a mixed-integer programming model and design a heuristic algorithm based on network optimization and diffusion simulation (hereafter NODS). The distinguishing contributions of the paper are as follows. First, according to the assumptions of the hazardous chemical accident, we simulate the dynamic change of the hazardous chemicals’ concentration and the range of influence of the hazardous chemical accident with a MATLAB program. Second, considering the toxicity and the diffusion of hazardous chemicals, we build a mixed-integer programming model that aims to minimize the concentration of hazardous chemicals to which people are exposed in the entire evacuation process. Then, a heuristic algorithm is designed to solve the model. Finally, a computational example based on an incident in Jinan shows the validity and feasibility of our algorithm, and we explore the key factors affecting the evacuation process through different scenarios.The remainder of the paper is organized as follows. Section2 introduces the physicochemical property of chlorine and the results of simulation. In Section 3, we introduce the mixed-integer programming model and the NODS algorithm. The computational example is carried out in Section 4. Finally, we conclude the work with some recommended future directions in Section 5.


## 2. Diffusion Simulation

### 2.1. Physicochemical Property of Chlorine

This paper considers common chlorine for the research of hazardous chemical accidents and emergency evacuation. At normal temperature and pressure, chloride is a type of gas with a yellow-green color and an irritating odor. It will transform into a liquid when the surrounding air pressure is more than 709 kPa. The molecular weight of chlorine is 70.91. Its melting point is −101 degrees Celsius, and its boiling point is −34.5 degrees Celsius. Its gas density is 3.21 g/L, and its relative vapor density is 2.5 (air is 1). Chloride poses dangers from combustion and explosion, and it can be used as an auxiliary combustion agent. Most general fuel can burn in the presence of chlorine gas, and flammable gas or vapor can form explosive mixtures with chlorine gas. Moreover, chloride can react with water, producing toxic hypochlorous acid and hydrochloric acid.In daily production and life, chlorine has a wide range of uses. It is used not only in the manufacture of vinyl chloride, allyl chloride, and epichlorohydrin but also as a disinfectant in the process of water treatment. On the other hand, chlorine is a highly toxic and hazardous chemical. Inhalation of high concentrations of chlorine gas can lead to death, and chlorine is one of the key hazardous chemicals in state supervision.Chlorine causes serious harm to the human body, for it can irritate the eyes, nose, throat, upper respiratory tract, and other body parts. When the mass concentration is higher than 1~6 mg/m3, it will cause significant irritation to humans. A concentration of 12 mg/m3 is unbearable for most people. A concentration of 90 mg/m3 can cause acute coughing. A concentration of 120~180 mg/m3 is very dangerous and can cause acute pulmonary edema and pneumonia in 30~60 min. A concentration of 300 mg/m3 can cause fatal injury. A concentration of 3000 mg/m3 can cause immediate respiratory failure, and a “stroke death” will occur. When the concentration is more than 30000 mg/m3, general viral gas masks lose their protective effect.Under normal temperature and pressure, liquid chlorine can be compressed and stored in liquefied form. Once it leaks into the atmosphere, it will expand and gasify under normal pressure, spreading to occupy a large space. In this paper, the liquid chlorine leakage is caused by an accident following damage to a container filling valve. After the accident, changes in the state of liquid chlorine can be divided into the following three processes.(
1
) Liquid chlorine leaks: pressurized liquid chlorine is leaked to the environment at normal temperature and pressure.(
2
) Liquid chlorine flashes: in a serious accident, liquid chlorine is quickly leaked to the surrounding atmosphere. Because of the pressure, a part of the liquid chlorine quickly converts into gas. The process of transitioning from the gas-liquid equilibrium under high pressure into the gas-liquid equilibrium under normal pressure is named the flash. The ratio of direct evaporation to liquid is also called the flash rate Fvap:(1)Fvap=CpΔTL.In the formula,Fvap is the flash rate, Cp is the mean specific heat of gas (kJ/kgK), ΔT is the temperature difference between the temperature and the boiling point in the container under the environmental pressure (K), and L is the heat of vaporization (kJ/kg). Because the mean specific heat of liquid chlorine is 0.98 kJ/kgK, the boiling point is −34.5 degrees Celsius, and the heat of vaporization is 289 kJ/kg, the flash rate of liquid chlorine is approximately 20% under normal atmospheric temperature (25 degrees Celsius). When liquid chlorine leaks, the liquid which evaporates directly will form clouds in the form of tiny smoke and mix with air, which then absorbs heat and evaporates. If the heat passed from air to the liquid smoke is not enough to cause evaporation, the smoke will condense into droplets and land on the ground, forming a liquid pool. According to experience, when F>0.2, usually it will not form a liquid pool. This paper assumes that when the percentage of liquid chlorine leakage is 20%, it will directly turn into gas.(
3
) The diffusion of gaseous chlorine: liquid chlorine diffuses in the air after the gasification.In this paper, we select the basic equation of fluid mechanics Bernoulli equation to describe the liquid chlorine leakage rate. To depict gas diffusion behavior after liquid chlorine flashes, the paper uses the Gaussian smoke model, which is the most widely used model of gas diffusion. The leakage model, the diffusion parameter, and the Gaussian smoke cloud model have been introduced in detail, respectively, in literatures [33] and [34]. Due to word constraints, they are not discussed further in this paper.


### 2.2. Simulation Calculation

The valve of a tank car which is loaded with 20 tons’ liquid chlorine has been completely damaged and the liquid chlorine leaks continuously. The valve diameter is 25 mm, leakage is 50%, and the wind speed is 5.4 meters/seconds. The key MATLAB codes for the simulation calculation are shown in Algorithms1, 2, and 3.Algorithm 1:GetOperMin(P).
(
1
)  % Calculate the liquid chlorine leak rate C in every moment, the unit is kg/s
(
2
)  %  P represents the pressure inside the liquid chlorine tank trucks
(
3
)  %  g represents the gravitational acceleration
(
4
)  %  Cd represents the liquid leakage coefficient
(
5
)  %  A represents the fissure area
(
6
)  %  P0 represents the environmental pressure
(
7
)  %  h represents the height of the valve from the ground
(
8
)  functionC=getQperMin(P)
(
9
)    P0=1.0e+005;
(
10
)   g=9.8;
(
11
)   h=0;
(
12
)   Cd=0.5;
(
13
)   A=pi∗0.025∧2;
(
14
)   rol=1.468e+003;
(
15
)   C=Cd∗A∗rol∗sqlt(2∗(P-P0)/rol+2gh);Algorithm 2:GetConcentrationPerQ(u,t,x,y,z).
(
1
)  % Calculate the concentration on some point at some moment, the unit is mg/s
(
2
)  %  dx, dy, dz represent the diffusion coefficient on the x, y, z axis respectively
(
3
)  %  x, y, z represent the coordinate values on the x, y, z axis
(
4
)  %  t represents the diffusion time, the unit is s
(
5
)  %  u represents the average wind speed, the unit is m/s
(
6
)  functionc=getConcentrationPerQ  (u,t,x,y,z)
(
7
)   dx=0.425809∗abs(x).∧0.901074;
(
8
)   dy=dx;
(
9
)   dz=0.0799904∗abs(x).∧1.12154;
(
10
)   c=exp(-(x-u∗t)∧2/(2∗dx∧2)-y∧2/(2∗dy∧2)-z∧2/(2∗dz∧2))/(sqrt(2)∗pi∧1.5∗dx∗dy∗dz);
(
11
)   c(isnan(c))=0;Algorithm 3:The simulating calculation (part).
(
1
)   % According to the hypothetical situation in this article, carries on the simulating calculation
(
2
)   %  u represents the average wind speed, the unit is m/s
(
3
)   %  Q0 represents the capacity of liquid chlorine tank trucks, Q represents the liquid chlorine leakage, the unit is kg
(
4
)   %  deltaT represents the time interval of measurement
(
5
)   %  stopTime represents the stop time of measurement
(
6
)   %  P1 represents the initial pressure in the liquid chlorine tank trucks
(
7
)   %  P represents the real-time pressure inside the liquid chlorine tank trucks
(
8
)   u=5.4;
(
9
)   Q0=2e+004;
(
10
)   P1=10e+005; P=P1;
(
11
)   deltaT=9;
(
12
)   stopTime=T;
(
13
)   deltaQ=0;
(
14
)   % Calculate the liquid chlorine leakage at every moment
(
15
)  whileQ>0
(
16
)   deltaQPerMin=getQPerMin(P)∗deltaT;
(
17
)   QPerMin=[QPerMin,deltaQPerMin];
(
18
)   deltaQ=deltaQ+deltaQPerMin;
(
19
)    P=P0∗0.2/(0.2+0.8∗deltaQ/Q0);
(
20
)    Q = Q-deltaQPerMin;
(
21
)
end
(
22
)   % Calculate the concentration on some point at some moment
(
23
)
for
t
=
0:deltaT:stopTime
(
24
)    Cxyt(:;:;1+t/deltaT)=getConcentrationPerQ(u,t,x,y,0);
(
25
)
end
(
26
)
for
t
=
d
e
l
t
a
T:deltaT:stopTime
(
27
)    Cxy=0;
(
28
)   fors=1:min(t/deltaT,length(QPerMin))
(
29
)     Cxy=QPerMin(s)∗1000000∗Cxyt(:;:;t/deltaT-s+1);
(
30
)   end
(
31
)
endAccording to the liquid chlorine concentration, the affected areas can be divided into lethal circle(3000×10-6), semilethal circle (300×10-6~3000×10-6), and injury circle (30×10-6~300×10-6) which will lead to death, disability, and injury, respectively. The concentration diffusion of the lethal circle, the semilethal circle, and the injury circle is shown in Figure 1 and Table 1.Table 1
Data table for lethal circle, semilethal circle, and injury circle.
Leakage time span (s)
Lethal circle
Semilethal circle
Injury circle
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
0
0
15
6
0
27
9
0
36
12
3
0
39
18
0
60
24
0
69
30
12
0
84
30
0
123
51
0
141
57
21
0
117
36
0
171
66
0
198
78
30
0
138
36
0
216
78
0
249
93
33
0
144
36
0
231
81
0
264
99
36
0
150
36
0
243
84
0
279
102
39
0
150
33
0
258
87
0
294
108
45
0
153
33
0
282
90
0
324
114
54
0
150
33
0
315
93
0
366
123
63
0
147
33
0
348
96
0
405
129
69
0
147
33
0
366
99
0
429
135
72
0
144
33
0
375
99
0
441
135
75
0
144
33
0
384
99
0
453
138
78
0
144
33
0
393
96
0
465
138
81
0
144
33
0
399
96
0
477
141
84
0
141
33
0
408
96
0
489
141
87
0
141
33
0
414
96
0
498
141
90
0
141
33
0
420
96
0
510
144
93
0
141
33
0
426
96
0
519
144
96
0
138
30
0
432
96
0
531
144
99
0
138
30
0
438
93
0
540
144
102
0
138
30
0
441
93
0
549
144
105
0
138
30
0
447
93
0
561
144
108
0
138
30
0
450
93
0
570
144
111
0
135
30
0
453
93
0
579
144
114
0
135
30
0
453
93
0
588
144
117
0
135
30
0
456
93
0
597
144
120
0
135
30
0
456
90
0
603
144
123
0
135
30
0
456
90
0
612
144
126
0
135
30
0
456
90
0
621
144
129
0
132
30
0
456
90
0
627
141
132
0
132
30
0
453
90
0
636
141
240
0
120
27
0
393
81
0
651
126
303
30
114
27
27
375
78
24
609
120
306
45
114
24
36
372
78
33
606
120
309
66
108
15
48
372
78
45
606
120
312
0
0
0
60
372
78
54
603
120
315
0
0
0
72
369
75
66
603
117
339
0
0
0
216
348
45
177
591
111
342
0
0
0
246
336
33
195
588
111
345
0
0
0
0
0
0
210
585
108
348
0
0
0
0
0
0
228
582
105
372
0
0
0
0
0
0
426
522
33Figure 1
The concentration diffusion with different leakage time.
(a)
The concentration diffusion figure when the leakage has lasted for 12 seconds
(b)
The concentration diffusion figure when the leakage has lasted for 45 seconds
(c)
The concentration diffusion figure when the leakage has lasted for 270 seconds
(d)
The concentration diffusion figure when the leakage has lasted for 315 secondsTable1 shows that when the valve of the tank car loaded with 20 tons of liquid chlorine has been completely damaged and liquid chlorine leaks continuously, the valve diameter of the tank is 25 mm, the leakage is 50%, and the wind speed is 5.4 meters/second; the lateral distances of the lethal circle, the semilethal circle, and the injury circle are 36 m, 99 m, and 144 m, respectively. The maximum distance at which people will be injured in a downwind direction is 651 m. Thus, we should evacuate people in the affected area to locations 651 m from the accident in the downwind direction and 144 m from the accident in the longitudinal distance to minimize casualties. We obtain the influenced area by simulation, as shown in Figure 2. The area of influence is 252,210 m2.Figure 2
The influenced area of liquid chlorine leakage.


## 2.1. Physicochemical Property of Chlorine

This paper considers common chlorine for the research of hazardous chemical accidents and emergency evacuation. At normal temperature and pressure, chloride is a type of gas with a yellow-green color and an irritating odor. It will transform into a liquid when the surrounding air pressure is more than 709 kPa. The molecular weight of chlorine is 70.91. Its melting point is −101 degrees Celsius, and its boiling point is −34.5 degrees Celsius. Its gas density is 3.21 g/L, and its relative vapor density is 2.5 (air is 1). Chloride poses dangers from combustion and explosion, and it can be used as an auxiliary combustion agent. Most general fuel can burn in the presence of chlorine gas, and flammable gas or vapor can form explosive mixtures with chlorine gas. Moreover, chloride can react with water, producing toxic hypochlorous acid and hydrochloric acid.In daily production and life, chlorine has a wide range of uses. It is used not only in the manufacture of vinyl chloride, allyl chloride, and epichlorohydrin but also as a disinfectant in the process of water treatment. On the other hand, chlorine is a highly toxic and hazardous chemical. Inhalation of high concentrations of chlorine gas can lead to death, and chlorine is one of the key hazardous chemicals in state supervision.Chlorine causes serious harm to the human body, for it can irritate the eyes, nose, throat, upper respiratory tract, and other body parts. When the mass concentration is higher than 1~6 mg/m3, it will cause significant irritation to humans. A concentration of 12 mg/m3 is unbearable for most people. A concentration of 90 mg/m3 can cause acute coughing. A concentration of 120~180 mg/m3 is very dangerous and can cause acute pulmonary edema and pneumonia in 30~60 min. A concentration of 300 mg/m3 can cause fatal injury. A concentration of 3000 mg/m3 can cause immediate respiratory failure, and a “stroke death” will occur. When the concentration is more than 30000 mg/m3, general viral gas masks lose their protective effect.Under normal temperature and pressure, liquid chlorine can be compressed and stored in liquefied form. Once it leaks into the atmosphere, it will expand and gasify under normal pressure, spreading to occupy a large space. In this paper, the liquid chlorine leakage is caused by an accident following damage to a container filling valve. After the accident, changes in the state of liquid chlorine can be divided into the following three processes.(
1
) Liquid chlorine leaks: pressurized liquid chlorine is leaked to the environment at normal temperature and pressure.(
2
) Liquid chlorine flashes: in a serious accident, liquid chlorine is quickly leaked to the surrounding atmosphere. Because of the pressure, a part of the liquid chlorine quickly converts into gas. The process of transitioning from the gas-liquid equilibrium under high pressure into the gas-liquid equilibrium under normal pressure is named the flash. The ratio of direct evaporation to liquid is also called the flash rate Fvap:(1)Fvap=CpΔTL.In the formula,Fvap is the flash rate, Cp is the mean specific heat of gas (kJ/kgK), ΔT is the temperature difference between the temperature and the boiling point in the container under the environmental pressure (K), and L is the heat of vaporization (kJ/kg). Because the mean specific heat of liquid chlorine is 0.98 kJ/kgK, the boiling point is −34.5 degrees Celsius, and the heat of vaporization is 289 kJ/kg, the flash rate of liquid chlorine is approximately 20% under normal atmospheric temperature (25 degrees Celsius). When liquid chlorine leaks, the liquid which evaporates directly will form clouds in the form of tiny smoke and mix with air, which then absorbs heat and evaporates. If the heat passed from air to the liquid smoke is not enough to cause evaporation, the smoke will condense into droplets and land on the ground, forming a liquid pool. According to experience, when F>0.2, usually it will not form a liquid pool. This paper assumes that when the percentage of liquid chlorine leakage is 20%, it will directly turn into gas.(
3
) The diffusion of gaseous chlorine: liquid chlorine diffuses in the air after the gasification.In this paper, we select the basic equation of fluid mechanics Bernoulli equation to describe the liquid chlorine leakage rate. To depict gas diffusion behavior after liquid chlorine flashes, the paper uses the Gaussian smoke model, which is the most widely used model of gas diffusion. The leakage model, the diffusion parameter, and the Gaussian smoke cloud model have been introduced in detail, respectively, in literatures [33] and [34]. Due to word constraints, they are not discussed further in this paper.


## 2.2. Simulation Calculation

The valve of a tank car which is loaded with 20 tons’ liquid chlorine has been completely damaged and the liquid chlorine leaks continuously. The valve diameter is 25 mm, leakage is 50%, and the wind speed is 5.4 meters/seconds. The key MATLAB codes for the simulation calculation are shown in Algorithms1, 2, and 3.Algorithm 1:GetOperMin(P).
(
1
)  % Calculate the liquid chlorine leak rate C in every moment, the unit is kg/s
(
2
)  %  P represents the pressure inside the liquid chlorine tank trucks
(
3
)  %  g represents the gravitational acceleration
(
4
)  %  Cd represents the liquid leakage coefficient
(
5
)  %  A represents the fissure area
(
6
)  %  P0 represents the environmental pressure
(
7
)  %  h represents the height of the valve from the ground
(
8
)  functionC=getQperMin(P)
(
9
)    P0=1.0e+005;
(
10
)   g=9.8;
(
11
)   h=0;
(
12
)   Cd=0.5;
(
13
)   A=pi∗0.025∧2;
(
14
)   rol=1.468e+003;
(
15
)   C=Cd∗A∗rol∗sqlt(2∗(P-P0)/rol+2gh);Algorithm 2:GetConcentrationPerQ(u,t,x,y,z).
(
1
)  % Calculate the concentration on some point at some moment, the unit is mg/s
(
2
)  %  dx, dy, dz represent the diffusion coefficient on the x, y, z axis respectively
(
3
)  %  x, y, z represent the coordinate values on the x, y, z axis
(
4
)  %  t represents the diffusion time, the unit is s
(
5
)  %  u represents the average wind speed, the unit is m/s
(
6
)  functionc=getConcentrationPerQ  (u,t,x,y,z)
(
7
)   dx=0.425809∗abs(x).∧0.901074;
(
8
)   dy=dx;
(
9
)   dz=0.0799904∗abs(x).∧1.12154;
(
10
)   c=exp(-(x-u∗t)∧2/(2∗dx∧2)-y∧2/(2∗dy∧2)-z∧2/(2∗dz∧2))/(sqrt(2)∗pi∧1.5∗dx∗dy∗dz);
(
11
)   c(isnan(c))=0;Algorithm 3:The simulating calculation (part).
(
1
)   % According to the hypothetical situation in this article, carries on the simulating calculation
(
2
)   %  u represents the average wind speed, the unit is m/s
(
3
)   %  Q0 represents the capacity of liquid chlorine tank trucks, Q represents the liquid chlorine leakage, the unit is kg
(
4
)   %  deltaT represents the time interval of measurement
(
5
)   %  stopTime represents the stop time of measurement
(
6
)   %  P1 represents the initial pressure in the liquid chlorine tank trucks
(
7
)   %  P represents the real-time pressure inside the liquid chlorine tank trucks
(
8
)   u=5.4;
(
9
)   Q0=2e+004;
(
10
)   P1=10e+005; P=P1;
(
11
)   deltaT=9;
(
12
)   stopTime=T;
(
13
)   deltaQ=0;
(
14
)   % Calculate the liquid chlorine leakage at every moment
(
15
)  whileQ>0
(
16
)   deltaQPerMin=getQPerMin(P)∗deltaT;
(
17
)   QPerMin=[QPerMin,deltaQPerMin];
(
18
)   deltaQ=deltaQ+deltaQPerMin;
(
19
)    P=P0∗0.2/(0.2+0.8∗deltaQ/Q0);
(
20
)    Q = Q-deltaQPerMin;
(
21
)
end
(
22
)   % Calculate the concentration on some point at some moment
(
23
)
for
t
=
0:deltaT:stopTime
(
24
)    Cxyt(:;:;1+t/deltaT)=getConcentrationPerQ(u,t,x,y,0);
(
25
)
end
(
26
)
for
t
=
d
e
l
t
a
T:deltaT:stopTime
(
27
)    Cxy=0;
(
28
)   fors=1:min(t/deltaT,length(QPerMin))
(
29
)     Cxy=QPerMin(s)∗1000000∗Cxyt(:;:;t/deltaT-s+1);
(
30
)   end
(
31
)
endAccording to the liquid chlorine concentration, the affected areas can be divided into lethal circle(3000×10-6), semilethal circle (300×10-6~3000×10-6), and injury circle (30×10-6~300×10-6) which will lead to death, disability, and injury, respectively. The concentration diffusion of the lethal circle, the semilethal circle, and the injury circle is shown in Figure 1 and Table 1.Table 1
Data table for lethal circle, semilethal circle, and injury circle.
Leakage time span (s)
Lethal circle
Semilethal circle
Injury circle
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
The left endpoint distance away from the scene
The right endpoint distance away from the scene
Longitudinal diffusion distance
0
0
15
6
0
27
9
0
36
12
3
0
39
18
0
60
24
0
69
30
12
0
84
30
0
123
51
0
141
57
21
0
117
36
0
171
66
0
198
78
30
0
138
36
0
216
78
0
249
93
33
0
144
36
0
231
81
0
264
99
36
0
150
36
0
243
84
0
279
102
39
0
150
33
0
258
87
0
294
108
45
0
153
33
0
282
90
0
324
114
54
0
150
33
0
315
93
0
366
123
63
0
147
33
0
348
96
0
405
129
69
0
147
33
0
366
99
0
429
135
72
0
144
33
0
375
99
0
441
135
75
0
144
33
0
384
99
0
453
138
78
0
144
33
0
393
96
0
465
138
81
0
144
33
0
399
96
0
477
141
84
0
141
33
0
408
96
0
489
141
87
0
141
33
0
414
96
0
498
141
90
0
141
33
0
420
96
0
510
144
93
0
141
33
0
426
96
0
519
144
96
0
138
30
0
432
96
0
531
144
99
0
138
30
0
438
93
0
540
144
102
0
138
30
0
441
93
0
549
144
105
0
138
30
0
447
93
0
561
144
108
0
138
30
0
450
93
0
570
144
111
0
135
30
0
453
93
0
579
144
114
0
135
30
0
453
93
0
588
144
117
0
135
30
0
456
93
0
597
144
120
0
135
30
0
456
90
0
603
144
123
0
135
30
0
456
90
0
612
144
126
0
135
30
0
456
90
0
621
144
129
0
132
30
0
456
90
0
627
141
132
0
132
30
0
453
90
0
636
141
240
0
120
27
0
393
81
0
651
126
303
30
114
27
27
375
78
24
609
120
306
45
114
24
36
372
78
33
606
120
309
66
108
15
48
372
78
45
606
120
312
0
0
0
60
372
78
54
603
120
315
0
0
0
72
369
75
66
603
117
339
0
0
0
216
348
45
177
591
111
342
0
0
0
246
336
33
195
588
111
345
0
0
0
0
0
0
210
585
108
348
0
0
0
0
0
0
228
582
105
372
0
0
0
0
0
0
426
522
33Figure 1
The concentration diffusion with different leakage time.
(a)
The concentration diffusion figure when the leakage has lasted for 12 seconds
(b)
The concentration diffusion figure when the leakage has lasted for 45 seconds
(c)
The concentration diffusion figure when the leakage has lasted for 270 seconds
(d)
The concentration diffusion figure when the leakage has lasted for 315 secondsTable1 shows that when the valve of the tank car loaded with 20 tons of liquid chlorine has been completely damaged and liquid chlorine leaks continuously, the valve diameter of the tank is 25 mm, the leakage is 50%, and the wind speed is 5.4 meters/second; the lateral distances of the lethal circle, the semilethal circle, and the injury circle are 36 m, 99 m, and 144 m, respectively. The maximum distance at which people will be injured in a downwind direction is 651 m. Thus, we should evacuate people in the affected area to locations 651 m from the accident in the downwind direction and 144 m from the accident in the longitudinal distance to minimize casualties. We obtain the influenced area by simulation, as shown in Figure 2. The area of influence is 252,210 m2.Figure 2
The influenced area of liquid chlorine leakage.


## 3. Model and Algorithms

In this paper, we specifically study emergency evacuation problems following hazardous chemical accidents. General emergency evacuation problems consider mostly how to evacuate all people to safe areas with the fastest speed or in the shortest time. Nevertheless, the biggest difference between emergency evacuation problems following hazardous chemical accidents and other evacuation problems is that the goal of emergency evacuation problems following hazardous chemical accidents is to minimize damage to people cause by the hazardous chemicals. This paper considers the emergency evacuation problem between a specific evacuation source point and a specific end point after a hazardous chemical accident in which the concentration of hazardous chemicals at each point changes at discrete points in time.


### 3.1. Symbols and Decision Variables

G
=
(
V
,
E
): traffic network diagram
V: point set; points represent path intersections, including the evacuation source point and end point; V=N, in which the evacuation source point is 1 and the end point is N
V
S: evacuation source point set
V
M: middle point set
V
E: evacuation end point set
E
=
(
e
i
j
): the set of an arc, where eij represents the arc between point i and j
E
S: starting arc set
E
E: ending arc set
I
C
S: the initial number of people who are waiting to be evacuated at the evacuation source point
M
C
N
i: the maximum capacity of point i
M
C
i
j: the maximum capacity of arc eij
C
N
i
(
t
): the capacity of point i at time
C
i
j
(
t
): the capacity of arc eij at time t
t
i
j: the length of time when travelling through arc eij
P
P
M
i
j
(
t
): the cumulative concentration on the arc of point i to point j at time t
C
o
n
i
t: the concentration on point i at time t
T
0: the beginning time of evacuation from the evacuation source point
P
k: the kth path of the evacuation source point
f
k: the flow of Pk, which is the kth path of the evacuation source point
x
i
j
(
t
): the number of people who enter arc eij at time t and arrive at point j at time t+tij
y
i
(
t
): the number of people on point i at time t
z
i
j
(
t
): the number of all people on the arc eij at time t.


### 3.2. Modeling

The mixed-integer programming problem is written as follows, selecting point 1 as the evacuation source point and pointN as the end point:(2)min∑eij∈Exijt·PPMijt+CN1t·Con1t(3)s.t.y10=ICS(4)yi0=0∀i∈V∖1(5)xij0=0∀i,j∈V(6)yit≤MCNi∀i∈V(7)xijt≤MCij∀i,j∈V(8)zijt=∑tz∈Tijxijtz∀i,j∈V(9)∑ni∈Vyit+∑eij∈Ezijt=ICS(10)yit-yit-1=∑k∈predixkit-tki-∑j∈succixijt∀i∈VM(11)yit=yit-1∑eij∈ESxijt∀i∈VS(12)yit=yit-1+∑eij∈EExijt-tij∀i∈VE(13)xijt≥0,yit≥0,zijt≥0,where predi=j∣j,i∈E, succi=j∣i,j∈E, and Tij=te∣t-tij+1≤te≤t,(i,j)∈E. The objective function (2) states that the goal of the problem is to minimize the concentration of hazardous chemicals to which all people who are waiting to be evacuated are exposed during the entire evacuation process. Constraint functions (3), (4), and (5) represent the initial state of the whole evacuation network. Equations (6) and (7) represent the capacity limit constraints of the point i and arc eij at time t. Equation (8) represents the conservation of flow on the arc eij. Equation (9) represents the conservation of the evacuation crowd in the whole evacuation network at time t. Equation (10) represents the conservation of flow at the middle point i. Equation (11) represents the conservation of flow at evacuation source point. Equation (12) represents the conservation of flow at the end point. Equation (13) represents the fact that the decision variables must be nonnegative.


### 3.3. Algorithm

Since evacuation problems are NP-hard problems and there is no polynomial algorithm, we design the NODS algorithm that can calculate evacuation plans following hazardous chemical accidents based on network optimization and diffusion simulation satisfying the constraints and objectives of the above mathematical model. The algorithm introduces a concept ofK shortest paths, which refer to the 1st,2nd,…,Kth shortest paths from evacuation source point to the end point. In this article, K shortest paths refer to the path corresponding to the 1st,2nd,…,Kth lowest cumulative concentration of chemicals from evacuation source point to end point. When designing the evacuation route, we only consider evacuating on the paths of the 1st,2nd,…,Kth lowest cumulative concentration of chemicals instead of evacuating on all the routes.


#### 3.3.1. Definitions and Characters

Definition 1.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l and the arrival time at point j is tj=ti+tij.Definition 2.
The label of pointj is (i,tj,m), and i represents the point on the path prior to point j, tj represents the arrival time at point j, and m represents the current stage.Definition 3.
L
m represents the point set that is contained in the stage of m. Pm-1(j) represents the former point set of stage m-1 of point j∀j∈Lm and that is Pm-1(j)=Lm-1∩P(j).Definition 4.
Assume that entering arceij from point i at time t, time tij is taken to arrive at point j and the concentration of hazardous chemicals at point i is ppm(i,t) at time t, and the concentration of hazardous chemicals at point j is to be ppm(j,t+tij) at time t+tij. Therefore, the concentration of hazardous chemicals that evacuees are exposed to when entering arc eij at the time of t is ppm(i,t)+ppm(j,t+tij)/2×tij, denoted by ppm(i,j,t).Definition 5.
The initial concentration of dangerous chemicals at the starting point at timet is ppm(i,t)+ppm(i,t-1).Definition 6.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l, the cumulative concentration of hazardous chemicals that the evacuees are exposed to is PPM(j,t)=PPM(i,ti)+ppm(i,j,ti) when arriving at point j at time t, where t=ti+tij.Definition 7.
The lowest cumulative concentration of hazardous chemicals at pointj in Lm when arriving at stage m at the time of tj is PPMm(j,tj); that is, (14)PPMmj,tj=mini∈Ti,tPPMm-1i,ti+ppmi,j,ti,where T(i,tj)={i∣tj=ti+tij,yi(t)≤MCNi,xij(t)≤MCij,t∈[0,tj],i∈Pm-1(j)}.Definition 8.
P
P
M
∗
(
i
) represents the lowest cumulative concentration of hazardous chemicals from the evacuation source point to point i.Property 9.
One has(15)PPM∗i=minm∈Mmin0≤ti≤TPPMmi,ti.Property 10.
The subpart of the lowest cumulative concentration of hazardous chemicals from the evacuation source point to pointi may not be the path of the lowest cumulative concentration of hazardous chemicals.


#### 3.3.2. The Step of the NODS Algorithm

The basic idea of the NODS algorithm is to initially, at each evacuation timet, calculate the K shortest paths and their flows according to the concentration of chemicals at each point at each different moment and then evacuate through the K shortest paths that have been calculated, until all of the people at the evacuation source point have been evacuated.Step 1 (initialization).
LetICS be the number of persons that will be evacuated at the evacuation source point, MCNi be the maximum capacity of point i, MCij be the maximum capacity of arc eij, and PPM(vi,t) be the concentration of hazardous chemicals of point i at t.
Let the capacity of pointi at t be CNi(t)=MCNi, and the capacity of arc eij at t be Cij(t)=MCij;
Lett0=0.Step 2 (divide phases).
Divide each stage according to the number of the arc paths on the evacuation network. Then, we can get the number of the stagesM. Let Lm be the setoff points that belong to stage m and let m=0.Step 3.
The start time of the evacuation at the source point ist=t0.Step 4.
Search theK shortest paths from source point 1 to point N; let k=1.Step 5.
Compute the cumulative concentration and the label for every point in stagem.
Calculate the minimum cumulative concentration of pointj for which the arrival time is tj. The minimum cumulative concentration is PPMm(j,tj)=mini∈T(i,tj)⁡{PPMm-1(i,ti)+PPM(i,j,ti)} and label (i,tj,m) to point j.Step 6 (determine whether the stagem is the end).
LetLm=Lm\{vj}. If Lm=φ, then all the points in stage m have been calculated. If m=M, then turn to Step 8. Otherwise, let m=m+1 and return to Step 5.Step 7 (compute the minimum cumulative concentration).
The minimum cumulative concentration isPPM∗N=minm⁡PPMm(N,tN).Step 8 (backtrack the shortest path).
Backtrack the point according to the label that has the minimum objective value. Then, we can get thekth shortest path for evacuees from source point 1 at t0. We denote the kth shortest path from point 1 to point N by Pk=1,n1k,n2k,…,N and calculate the flow fk and the total time tk(N) of Pk.Step 9 (update the evacuation network).
Ifk=K, then return to Step 10; else, calculate MCij=MCij-fak,  MCNi=MCNi-fak. If MCij=0, then delete arc eij from the network; if MCNi=0, then delete all the former arcs adjacent to point i. If the updated network is not connected, return to Step 10. Otherwise, let k=k+1, and turn to the Step 4.Step 10 (determine whether the evacuation is finished).
Check the evacuation source point. If the victims are evacuated, then the algorithm ends. Otherwise, lett0=t0+1, and return to Step 3.


## 3.1. Symbols and Decision Variables

G
=
(
V
,
E
): traffic network diagram
V: point set; points represent path intersections, including the evacuation source point and end point; V=N, in which the evacuation source point is 1 and the end point is N
V
S: evacuation source point set
V
M: middle point set
V
E: evacuation end point set
E
=
(
e
i
j
): the set of an arc, where eij represents the arc between point i and j
E
S: starting arc set
E
E: ending arc set
I
C
S: the initial number of people who are waiting to be evacuated at the evacuation source point
M
C
N
i: the maximum capacity of point i
M
C
i
j: the maximum capacity of arc eij
C
N
i
(
t
): the capacity of point i at time
C
i
j
(
t
): the capacity of arc eij at time t
t
i
j: the length of time when travelling through arc eij
P
P
M
i
j
(
t
): the cumulative concentration on the arc of point i to point j at time t
C
o
n
i
t: the concentration on point i at time t
T
0: the beginning time of evacuation from the evacuation source point
P
k: the kth path of the evacuation source point
f
k: the flow of Pk, which is the kth path of the evacuation source point
x
i
j
(
t
): the number of people who enter arc eij at time t and arrive at point j at time t+tij
y
i
(
t
): the number of people on point i at time t
z
i
j
(
t
): the number of all people on the arc eij at time t.


## 3.2. Modeling

The mixed-integer programming problem is written as follows, selecting point 1 as the evacuation source point and pointN as the end point:(2)min∑eij∈Exijt·PPMijt+CN1t·Con1t(3)s.t.y10=ICS(4)yi0=0∀i∈V∖1(5)xij0=0∀i,j∈V(6)yit≤MCNi∀i∈V(7)xijt≤MCij∀i,j∈V(8)zijt=∑tz∈Tijxijtz∀i,j∈V(9)∑ni∈Vyit+∑eij∈Ezijt=ICS(10)yit-yit-1=∑k∈predixkit-tki-∑j∈succixijt∀i∈VM(11)yit=yit-1∑eij∈ESxijt∀i∈VS(12)yit=yit-1+∑eij∈EExijt-tij∀i∈VE(13)xijt≥0,yit≥0,zijt≥0,where predi=j∣j,i∈E, succi=j∣i,j∈E, and Tij=te∣t-tij+1≤te≤t,(i,j)∈E. The objective function (2) states that the goal of the problem is to minimize the concentration of hazardous chemicals to which all people who are waiting to be evacuated are exposed during the entire evacuation process. Constraint functions (3), (4), and (5) represent the initial state of the whole evacuation network. Equations (6) and (7) represent the capacity limit constraints of the point i and arc eij at time t. Equation (8) represents the conservation of flow on the arc eij. Equation (9) represents the conservation of the evacuation crowd in the whole evacuation network at time t. Equation (10) represents the conservation of flow at the middle point i. Equation (11) represents the conservation of flow at evacuation source point. Equation (12) represents the conservation of flow at the end point. Equation (13) represents the fact that the decision variables must be nonnegative.


## 3.3. Algorithm

Since evacuation problems are NP-hard problems and there is no polynomial algorithm, we design the NODS algorithm that can calculate evacuation plans following hazardous chemical accidents based on network optimization and diffusion simulation satisfying the constraints and objectives of the above mathematical model. The algorithm introduces a concept ofK shortest paths, which refer to the 1st,2nd,…,Kth shortest paths from evacuation source point to the end point. In this article, K shortest paths refer to the path corresponding to the 1st,2nd,…,Kth lowest cumulative concentration of chemicals from evacuation source point to end point. When designing the evacuation route, we only consider evacuating on the paths of the 1st,2nd,…,Kth lowest cumulative concentration of chemicals instead of evacuating on all the routes.


### 3.3.1. Definitions and Characters

Definition 1.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l and the arrival time at point j is tj=ti+tij.Definition 2.
The label of pointj is (i,tj,m), and i represents the point on the path prior to point j, tj represents the arrival time at point j, and m represents the current stage.Definition 3.
L
m represents the point set that is contained in the stage of m. Pm-1(j) represents the former point set of stage m-1 of point j∀j∈Lm and that is Pm-1(j)=Lm-1∩P(j).Definition 4.
Assume that entering arceij from point i at time t, time tij is taken to arrive at point j and the concentration of hazardous chemicals at point i is ppm(i,t) at time t, and the concentration of hazardous chemicals at point j is to be ppm(j,t+tij) at time t+tij. Therefore, the concentration of hazardous chemicals that evacuees are exposed to when entering arc eij at the time of t is ppm(i,t)+ppm(j,t+tij)/2×tij, denoted by ppm(i,j,t).Definition 5.
The initial concentration of dangerous chemicals at the starting point at timet is ppm(i,t)+ppm(i,t-1).Definition 6.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l, the cumulative concentration of hazardous chemicals that the evacuees are exposed to is PPM(j,t)=PPM(i,ti)+ppm(i,j,ti) when arriving at point j at time t, where t=ti+tij.Definition 7.
The lowest cumulative concentration of hazardous chemicals at pointj in Lm when arriving at stage m at the time of tj is PPMm(j,tj); that is, (14)PPMmj,tj=mini∈Ti,tPPMm-1i,ti+ppmi,j,ti,where T(i,tj)={i∣tj=ti+tij,yi(t)≤MCNi,xij(t)≤MCij,t∈[0,tj],i∈Pm-1(j)}.Definition 8.
P
P
M
∗
(
i
) represents the lowest cumulative concentration of hazardous chemicals from the evacuation source point to point i.Property 9.
One has(15)PPM∗i=minm∈Mmin0≤ti≤TPPMmi,ti.Property 10.
The subpart of the lowest cumulative concentration of hazardous chemicals from the evacuation source point to pointi may not be the path of the lowest cumulative concentration of hazardous chemicals.


### 3.3.2. The Step of the NODS Algorithm

The basic idea of the NODS algorithm is to initially, at each evacuation timet, calculate the K shortest paths and their flows according to the concentration of chemicals at each point at each different moment and then evacuate through the K shortest paths that have been calculated, until all of the people at the evacuation source point have been evacuated.Step 1 (initialization).
LetICS be the number of persons that will be evacuated at the evacuation source point, MCNi be the maximum capacity of point i, MCij be the maximum capacity of arc eij, and PPM(vi,t) be the concentration of hazardous chemicals of point i at t.
Let the capacity of pointi at t be CNi(t)=MCNi, and the capacity of arc eij at t be Cij(t)=MCij;
Lett0=0.Step 2 (divide phases).
Divide each stage according to the number of the arc paths on the evacuation network. Then, we can get the number of the stagesM. Let Lm be the setoff points that belong to stage m and let m=0.Step 3.
The start time of the evacuation at the source point ist=t0.Step 4.
Search theK shortest paths from source point 1 to point N; let k=1.Step 5.
Compute the cumulative concentration and the label for every point in stagem.
Calculate the minimum cumulative concentration of pointj for which the arrival time is tj. The minimum cumulative concentration is PPMm(j,tj)=mini∈T(i,tj)⁡{PPMm-1(i,ti)+PPM(i,j,ti)} and label (i,tj,m) to point j.Step 6 (determine whether the stagem is the end).
LetLm=Lm\{vj}. If Lm=φ, then all the points in stage m have been calculated. If m=M, then turn to Step 8. Otherwise, let m=m+1 and return to Step 5.Step 7 (compute the minimum cumulative concentration).
The minimum cumulative concentration isPPM∗N=minm⁡PPMm(N,tN).Step 8 (backtrack the shortest path).
Backtrack the point according to the label that has the minimum objective value. Then, we can get thekth shortest path for evacuees from source point 1 at t0. We denote the kth shortest path from point 1 to point N by Pk=1,n1k,n2k,…,N and calculate the flow fk and the total time tk(N) of Pk.Step 9 (update the evacuation network).
Ifk=K, then return to Step 10; else, calculate MCij=MCij-fak,  MCNi=MCNi-fak. If MCij=0, then delete arc eij from the network; if MCNi=0, then delete all the former arcs adjacent to point i. If the updated network is not connected, return to Step 10. Otherwise, let k=k+1, and turn to the Step 4.Step 10 (determine whether the evacuation is finished).
Check the evacuation source point. If the victims are evacuated, then the algorithm ends. Otherwise, lett0=t0+1, and return to Step 3.


## 3.3.1. Definitions and Characters

Definition 1.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l and the arrival time at point j is tj=ti+tij.Definition 2.
The label of pointj is (i,tj,m), and i represents the point on the path prior to point j, tj represents the arrival time at point j, and m represents the current stage.Definition 3.
L
m represents the point set that is contained in the stage of m. Pm-1(j) represents the former point set of stage m-1 of point j∀j∈Lm and that is Pm-1(j)=Lm-1∩P(j).Definition 4.
Assume that entering arceij from point i at time t, time tij is taken to arrive at point j and the concentration of hazardous chemicals at point i is ppm(i,t) at time t, and the concentration of hazardous chemicals at point j is to be ppm(j,t+tij) at time t+tij. Therefore, the concentration of hazardous chemicals that evacuees are exposed to when entering arc eij at the time of t is ppm(i,t)+ppm(j,t+tij)/2×tij, denoted by ppm(i,j,t).Definition 5.
The initial concentration of dangerous chemicals at the starting point at timet is ppm(i,t)+ppm(i,t-1).Definition 6.
Assume thatP=(1,…,i,j,…,l) is a path from the evacuation source point to point l, the cumulative concentration of hazardous chemicals that the evacuees are exposed to is PPM(j,t)=PPM(i,ti)+ppm(i,j,ti) when arriving at point j at time t, where t=ti+tij.Definition 7.
The lowest cumulative concentration of hazardous chemicals at pointj in Lm when arriving at stage m at the time of tj is PPMm(j,tj); that is, (14)PPMmj,tj=mini∈Ti,tPPMm-1i,ti+ppmi,j,ti,where T(i,tj)={i∣tj=ti+tij,yi(t)≤MCNi,xij(t)≤MCij,t∈[0,tj],i∈Pm-1(j)}.Definition 8.
P
P
M
∗
(
i
) represents the lowest cumulative concentration of hazardous chemicals from the evacuation source point to point i.Property 9.
One has(15)PPM∗i=minm∈Mmin0≤ti≤TPPMmi,ti.Property 10.
The subpart of the lowest cumulative concentration of hazardous chemicals from the evacuation source point to pointi may not be the path of the lowest cumulative concentration of hazardous chemicals.


## 3.3.2. The Step of the NODS Algorithm

The basic idea of the NODS algorithm is to initially, at each evacuation timet, calculate the K shortest paths and their flows according to the concentration of chemicals at each point at each different moment and then evacuate through the K shortest paths that have been calculated, until all of the people at the evacuation source point have been evacuated.Step 1 (initialization).
LetICS be the number of persons that will be evacuated at the evacuation source point, MCNi be the maximum capacity of point i, MCij be the maximum capacity of arc eij, and PPM(vi,t) be the concentration of hazardous chemicals of point i at t.
Let the capacity of pointi at t be CNi(t)=MCNi, and the capacity of arc eij at t be Cij(t)=MCij;
Lett0=0.Step 2 (divide phases).
Divide each stage according to the number of the arc paths on the evacuation network. Then, we can get the number of the stagesM. Let Lm be the setoff points that belong to stage m and let m=0.Step 3.
The start time of the evacuation at the source point ist=t0.Step 4.
Search theK shortest paths from source point 1 to point N; let k=1.Step 5.
Compute the cumulative concentration and the label for every point in stagem.
Calculate the minimum cumulative concentration of pointj for which the arrival time is tj. The minimum cumulative concentration is PPMm(j,tj)=mini∈T(i,tj)⁡{PPMm-1(i,ti)+PPM(i,j,ti)} and label (i,tj,m) to point j.Step 6 (determine whether the stagem is the end).
LetLm=Lm\{vj}. If Lm=φ, then all the points in stage m have been calculated. If m=M, then turn to Step 8. Otherwise, let m=m+1 and return to Step 5.Step 7 (compute the minimum cumulative concentration).
The minimum cumulative concentration isPPM∗N=minm⁡PPMm(N,tN).Step 8 (backtrack the shortest path).
Backtrack the point according to the label that has the minimum objective value. Then, we can get thekth shortest path for evacuees from source point 1 at t0. We denote the kth shortest path from point 1 to point N by Pk=1,n1k,n2k,…,N and calculate the flow fk and the total time tk(N) of Pk.Step 9 (update the evacuation network).
Ifk=K, then return to Step 10; else, calculate MCij=MCij-fak,  MCNi=MCNi-fak. If MCij=0, then delete arc eij from the network; if MCNi=0, then delete all the former arcs adjacent to point i. If the updated network is not connected, return to Step 10. Otherwise, let k=k+1, and turn to the Step 4.Step 10 (determine whether the evacuation is finished).
Check the evacuation source point. If the victims are evacuated, then the algorithm ends. Otherwise, lett0=t0+1, and return to Step 3.


## 4. Case Study

This paper uses Jinan as an example to verify the validity and feasibility of the NODS algorithm. For simplicity, we get the number of evacuees from the influenced area and the average influenced population. That is, the number of evacuees equals the influenced area times the average influenced population.Suppose the tank car, which has a valve diameter of 25 mm and is loaded with 20 tons of liquid chlorine, had a traffic accident around the Jinbaolai sea cucumber store on Shanda Road at night that damaged the filling valve, causing liquid chlorine to continuously leak from the tank. Assuming the wind is from the west and the wind speed is 5.4 meters/second, we can obtain the influenced area in Lixia district, which is 252,210 m2. The average population of Lixia district is 7466.34 per/square kilometer. Thus, we can calculate the total number of people who need to be evacuated as 1,884. In addition, this paper assumes that the transverse distance from the evacuation source point (marked as the triangle in Figure 3) to the leakage point is 20 meters, the longitudinal distance is 20 meters, and all people who need to be evacuated are gathered at the evacuation source point. We assume that rescue vehicles are used for evacuation, the departure interval time of each vehicle is 9 seconds, and every 50 evacuees are one unit. Therefore, there are 38 units of victims that need to be evacuated to the emergency end point (marked as square in Figure 3) urgently. The evacuation map, as shown in Figure 3(a), is the road map for the evacuation of the whole evacuation zone. Figure 3(b) shows the specific locations of the evacuation source point and the leakage point.Figure 3
The evacuation map.
(a)
The road map for evacuation
(b)
The specific locations of the evacuation origin and the leakage pointThe evacuation map shown in Figure3 is represented by the network shown in Figure 4, where ① is the evacuation source point and its coordinate is (20, 20) and ㉘ is the emergency end point and its coordinate is (3010,150). Other points are arc crossing points. The first number in brackets is the maximum capacity of the corresponding arc, and the second number is the time to go through the arc. The concentration of hazardous chemicals at each point depends on time. For ease of calculation, we assume that time is discrete units. In this paper, we take 9 seconds as a unit. There are 38 units of victims at the evacuation source point ①, and we assume that the maximum capacity of all vertices is large enough to accommodate them all. We get the evacuation route by the algorithm that have been mentioned above. The result is shown in Table 2. Among them, the evacuation time is the time when evacuees leave the evacuation source point, and the evacuation completion time is the time when evacuees arrive at the emergency end point. In this section, all computational runs are made using MyEclipse 8.6 on a PC having Intel i3-3217U CPU and 4 GB memory.Table 2
Evacuation route results.
Evacuation time
Evacuation routes
Evacuation unit
Cumulative concentration(mg/m3)
Evacuation completion time
0
1-2-4-8-14-22-26-27-28
2
43.80
58
0
1-2-7-8-14-22-26-27-28
1
43.83
57
0
1-3-5-9-13-17-18-23-24-28
2
43.91
50
1
1-2-4-8-14-22-26-27-28
2
85.28
59
1
1-2-7-8-14-22-26-27-28
1
85.31
58
1
1-3-5-9-13-17-18-23-24-28
2
85.38
51
2
1-2-4-8-14-22-26-27-28
2
125.01
60
2
1-2-7-8-14-22-26-27-28
1
125.03
59
2
1-3-5-9-13-17-18-23-24-28
2
125.11
52
3
1-2-4-8-14-22-26-27-28
2
163.24
61
3
1-2-7-8-14-22-26-27-28
1
163.26
60
3
1-3-5-9-13-17-18-23-24-28
2
163.33
53
4
1-2-4-8-14-22-26-27-28
2
200.14
62
4
1-2-7-8-14-22-26-27-28
1
200.16
61
4
1-3-5-9-13-17-18-23-24-28
2
200.23
54
5
1-2-4-8-14-22-26-27-28
2
235.86
63
5
1-2-7-8-14-22-26-27-28
1
235.88
62
5
1-3-5-9-13-17-18-23-24-28
2
235.96
55
6
1-2-4-8-14-22-26-27-28
2
270.53
64
6
1-2-7-8-14-22-26-27-28
1
270.55
63
6
1-3-5-9-13-17-18-23-24-28
2
270.63
56
7
1-2-4-8-14-22-26-27-28
2
304.26
65
7
1-2-7-8-14-22-26-27-28
1
304.27
64Figure 4
The evacuation network.From Table2, we can see that when the number of people who need to be evacuated is 38 units, the shortest paths quantity is K=3, and the coordinates of the evacuation source point are (20,20). The three paths, 1-2-4-8-14-22-26-27-28, 1-2-7-8-14-22-26-27-28, and 1-3-5-9-13-17-18-23-24-28, are used for evacuation from time 0 to time 7. From Figure 4, we can see that the three paths selected in every moment are all on the edge of the evacuation network. Their longitudinal distances from the leakage point are all approximately 1000 meters, which is beyond the influence of the leakage. These three routes have the minimum cumulative concentration. Thus, it can be seen that in order to reduce the cumulative damage of evacuees in the entire evacuation process, the vertical direction of the wind should be the priority for evacuation, letting evacuees leave the area that has been affected by the leakage as soon as possible. After the evacuees have left the danger zone, they should then travel along the shortest paths to reach the emergency end point. From Table 2, we can see that it takes 7 units of time, or 63 seconds, to allow all the evacuees to evacuate out of the evacuation source point. It takes 65 units of time, or 585 seconds, for all evacuees to arrive at the emergency end point. As time goes on, the cumulative concentration of chlorine gas is increasing. The highest cumulative concentration reaches approximately 304.27 mg/m3 at time 7. This concentration will cause serious physical injury, but it is not a lethal concentration of chlorine gas. After 100 tests, the average program execution time of our algorithm is 20.2 seconds, and change in program execution time of each experiment is shown in Figure 5.Figure 5
The diagram of program execution time in 100 tests.Next, we change the evacuation scale, the shortest path quantity, and evacuation source point coordinates (the distance from the evacuation source point to the leakage point) to create different scenarios. Then, we solve the different scenarios with the algorithm that has been mentioned above and analyze the results.


### 4.1. Comparing the Calculation Results of Different Evacuation Scales

We take 50 evacuees as one unit of evacuees, the shortest path quantity to beK=3, and the evacuation source point coordinate to be (20,20) as a constant condition. We change the evacuation scale to create different scenarios and calculate them with the algorithm. Among them, 8, 18, and 28 units are considered small-scale, 38, 48, and 58 units are considered middle-scale, and 68, 78, and 88 units are considered large-scale. The calculated results are shown in Table 3. The evacuation time is the time when the evacuees leave the evacuation source point. All the time units are in seconds.Table 3
The calculative results of different evacuation scales.
Evacuation scale (units)
Evacuation routes
Evacuation time (s)
The highest cumulative concentration(mg/m3)
Program execution time (s)
8
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
9
85.31
5.6
18
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
27
163.26
10.1
28
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
45
235.89
14.9
38
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
63
304.27
19.8
48
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
81
369.16
24.8
58
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
99
431.09
28.3
68
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
117
490.47
32.8
78
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
135
547.61
39.2
88
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
162
602.75
43.2To analyze the influence of different evacuation scales on the result, we use the line chart to display the calculation results in Table3, as shown in Figure 6.Figure 6
The results of different evacuation scales.From Figure6, we can see that when the shortest path quantity and the evacuation source point coordinate are constant, with the increase of the evacuation scale, the evacuation time curve, the highest cumulative concentration curve, and the program execution time curve all increase with a linear growth trend. Compared to the program execution time curve, the evacuation time curve and the highest cumulative concentration curve show a relatively fast growth trend. This is because in the evacuation process with different evacuation scales, the shortest paths selected are the same in every moment, and the number of leaving evacuees is fixed in every moment. Therefore, the greater the evacuation scales, the longer the evacuation time. In every moment, the evacuees who stay in the evacuation source point will get more damaged with the passage of time, so the highest cumulative concentration will increase with the increase in evacuation time. The program execution time curve shows a relatively smooth linear growth trend, thus indicating that the algorithm is efficient and stable and demonstrates good performance. From what has been discussed above, when the shortest path quantity and the evacuation source point coordinate are constant, we should try to increase the number of evacuees on every path, such as using large rescue vehicles instead of small rescue vehicles, and clear the evacuation routes in order to improve the traffic load so that we can shorten the evacuation time and reduce the highest cumulative concentration.


### 4.2. Comparing the Calculation Results of Different Shortest Path Quantity

We take 50 evacuees to be one unit, the evacuation scale to be 38 units, and the evacuation source point coordinate to be(20,20) as a constant condition. We change the shortest path quantity K to create different scenarios and calculate them with the algorithm. Because the evacuation routes are limited and the traffic load of every path is limited at the same time, we can only find 6 shortest paths at most. Table 4 shows the calculation results when the shortest path quantity K takes values of 1, 2, 3, 4, 5, and 6. The evacuation time is the time when the evacuees leave the evacuation source point. All time units are seconds.Table 4
The calculative results of different shortest path quantity.
Shortest path quantity
Evacuation routes
Evacuation time (s)
The highest cumulative concentration(mg/m3)
Program execution time (s)
1
1-2-4-8-14-22-26-27-28
162
629.63
21.8
2
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28
108
461.06
20.3
3
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
304.27
19.8
4
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-28
54
270.54
20.3
5
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-281-3-6-7-11-15-16-19-24-28
45
235.89
22.2
6
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-281-3-6-7-11-15-16-19-24-281-3-6-10-12-16-19-24-28
36
200.47
22.3To analyze the influence of different shortest path quantity on the result, we use the line chart to display the calculation results in Table4, as shown in Figure 7.Figure 7
The results of different shortest path quantity.
(a)
The comprehensive comparison of every indicator
(b)
Program execution time curveWe can see from Figure7 that as the shortest path quantity increases, the evacuation time curve and the highest cumulative concentration curve show a downward trend. The shortest path quantity K=3 is a turning point. When K<3, the downward slope is big, and the curves are rapidly declining. When K>3, the downward slope is diminishing and the decline speed is also slowing down, which shows the shortest path quantity has an influence on the evacuation time and the highest cumulative concentration, but the influence reduces as the shortest path quantity increases. From Figure 7(b), we can see that the program execution time curve has low middle and high ends, which suggests that a larger shortest path quantity or a smaller shortest path quantity will reduce the efficiency of the algorithm. From what has been discussed above, the shortest path quantity should be taken as 3 or 4. This can ensure that the evacuation time, the highest cumulative concentration, and the program execution time all fall in an ideal range.


### 4.3. Comparing the Calculation Results of Different Evacuation Source Point Coordinates

In this case, we take every 50 evacuees to be one unit, the evacuation scale to be 38 units, and the shortest path quantity to beK=3 as a constant condition. We change the evacuation source point coordinates to create different scenarios and calculate them with the algorithm. According to the shown in Figure 1, we select (20,10), (20, 15), (20, 20), (20, 25), (20, 30), (20, 35), and (20, 40) as the coordinates of the evacuation source point. The calculation results are shown in Table 5. The evacuation time is the time when the evacuees leave the evacuation source point. All time units are in seconds.Table 5
The calculative results of different evacuation source point coordinates.
Evacuation sourcepoint coordinate
Evacuation routes
Evacuation time (s)
The highest cumulative concentration (mg/m3)
Program executiontime (s)
(20, 10)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
12816.29
19.2
(20, 15)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
2696.45
19.8
(20, 20)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
304.27
19.9
(20, 25)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
18.59
20.4
(20, 30)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.81
20.6
(20, 35)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.22
20.6
(20, 40)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.21
20.3To analyze the influence of different evacuation source point coordinates on the result, we use the line chart to display the calculated results in Table5, as shown in Figure 8.Figure 8
The results of different evacuation source point coordinates.
(a)
The comprehensive comparison of every indicator
(b)
The highest cumulative concentration of different evacuation source point abscissaIn this case, the evacuation scale and the shortest path quantity are held constant, making the evacuation source point abscissa constant. Changing the evacuation source point ordinate, from Figure8(a), we can see that the evacuation time curve shows a horizontal trend. From Table 5, we can see that changing the evacuation source point coordinates has no effect on selecting the shortest routes. It can be seen, with different evacuation source point coordinates, that the number of evacuees who leave the evacuation source point is always the same at every moment, so the evacuation time is also the same. This suggests that changing the evacuation source point coordinate has no significant effect on the evacuation time. As the evacuation source point coordinates change, the program execution curve has a certain fluctuation, but the range is just 1.4 seconds. This suggests that changing the evacuation source point coordinate has no significant effect on program execution time. As the evacuation source point abscissa increases, the highest cumulative concentration curve declines rapidly towards zero. This suggests that changing the evacuation source point coordinates has a significant effect on the highest cumulative concentration of the whole evacuation process. When the evacuation source point abscissa is a constant, the farther the ordinate and the lower the concentration of hazardous chemicals at the evacuation source point. This makes the highest cumulative concentration lower in the whole evacuation process. From Figure 8, we can see that as the evacuation source point ordinates increase, the disparity of highest cumulative concentration between different evacuation source point abscissas decreases, and the speed of decrease is faster and faster. Combined with Figure 1, we can see that the influence distance of the hazardous chemicals leakage in the downwind direction is much longer than that in the vertical direction of the wind, so we can know that the influence of the evacuation source point ordinates on the highest cumulative concentration is more significant than that of the abscissa. When we choose the evacuation source point, we should increase the ordinate to make the source point further away from the leakage point in the longitudinal direction, thereby reducing the highest cumulative concentration in the whole evacuation process.


## 4.1. Comparing the Calculation Results of Different Evacuation Scales

We take 50 evacuees as one unit of evacuees, the shortest path quantity to beK=3, and the evacuation source point coordinate to be (20,20) as a constant condition. We change the evacuation scale to create different scenarios and calculate them with the algorithm. Among them, 8, 18, and 28 units are considered small-scale, 38, 48, and 58 units are considered middle-scale, and 68, 78, and 88 units are considered large-scale. The calculated results are shown in Table 3. The evacuation time is the time when the evacuees leave the evacuation source point. All the time units are in seconds.Table 3
The calculative results of different evacuation scales.
Evacuation scale (units)
Evacuation routes
Evacuation time (s)
The highest cumulative concentration(mg/m3)
Program execution time (s)
8
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
9
85.31
5.6
18
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
27
163.26
10.1
28
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
45
235.89
14.9
38
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
63
304.27
19.8
48
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
81
369.16
24.8
58
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
99
431.09
28.3
68
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
117
490.47
32.8
78
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
135
547.61
39.2
88
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28 1-3-5-9-13-17-18-23-24-28
162
602.75
43.2To analyze the influence of different evacuation scales on the result, we use the line chart to display the calculation results in Table3, as shown in Figure 6.Figure 6
The results of different evacuation scales.From Figure6, we can see that when the shortest path quantity and the evacuation source point coordinate are constant, with the increase of the evacuation scale, the evacuation time curve, the highest cumulative concentration curve, and the program execution time curve all increase with a linear growth trend. Compared to the program execution time curve, the evacuation time curve and the highest cumulative concentration curve show a relatively fast growth trend. This is because in the evacuation process with different evacuation scales, the shortest paths selected are the same in every moment, and the number of leaving evacuees is fixed in every moment. Therefore, the greater the evacuation scales, the longer the evacuation time. In every moment, the evacuees who stay in the evacuation source point will get more damaged with the passage of time, so the highest cumulative concentration will increase with the increase in evacuation time. The program execution time curve shows a relatively smooth linear growth trend, thus indicating that the algorithm is efficient and stable and demonstrates good performance. From what has been discussed above, when the shortest path quantity and the evacuation source point coordinate are constant, we should try to increase the number of evacuees on every path, such as using large rescue vehicles instead of small rescue vehicles, and clear the evacuation routes in order to improve the traffic load so that we can shorten the evacuation time and reduce the highest cumulative concentration.


## 4.2. Comparing the Calculation Results of Different Shortest Path Quantity

We take 50 evacuees to be one unit, the evacuation scale to be 38 units, and the evacuation source point coordinate to be(20,20) as a constant condition. We change the shortest path quantity K to create different scenarios and calculate them with the algorithm. Because the evacuation routes are limited and the traffic load of every path is limited at the same time, we can only find 6 shortest paths at most. Table 4 shows the calculation results when the shortest path quantity K takes values of 1, 2, 3, 4, 5, and 6. The evacuation time is the time when the evacuees leave the evacuation source point. All time units are seconds.Table 4
The calculative results of different shortest path quantity.
Shortest path quantity
Evacuation routes
Evacuation time (s)
The highest cumulative concentration(mg/m3)
Program execution time (s)
1
1-2-4-8-14-22-26-27-28
162
629.63
21.8
2
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-28
108
461.06
20.3
3
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
304.27
19.8
4
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-28
54
270.54
20.3
5
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-281-3-6-7-11-15-16-19-24-28
45
235.89
22.2
6
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-281-2-7-11-15-20-21-25-27-281-3-6-7-11-15-16-19-24-281-3-6-10-12-16-19-24-28
36
200.47
22.3To analyze the influence of different shortest path quantity on the result, we use the line chart to display the calculation results in Table4, as shown in Figure 7.Figure 7
The results of different shortest path quantity.
(a)
The comprehensive comparison of every indicator
(b)
Program execution time curveWe can see from Figure7 that as the shortest path quantity increases, the evacuation time curve and the highest cumulative concentration curve show a downward trend. The shortest path quantity K=3 is a turning point. When K<3, the downward slope is big, and the curves are rapidly declining. When K>3, the downward slope is diminishing and the decline speed is also slowing down, which shows the shortest path quantity has an influence on the evacuation time and the highest cumulative concentration, but the influence reduces as the shortest path quantity increases. From Figure 7(b), we can see that the program execution time curve has low middle and high ends, which suggests that a larger shortest path quantity or a smaller shortest path quantity will reduce the efficiency of the algorithm. From what has been discussed above, the shortest path quantity should be taken as 3 or 4. This can ensure that the evacuation time, the highest cumulative concentration, and the program execution time all fall in an ideal range.


## 4.3. Comparing the Calculation Results of Different Evacuation Source Point Coordinates

In this case, we take every 50 evacuees to be one unit, the evacuation scale to be 38 units, and the shortest path quantity to beK=3 as a constant condition. We change the evacuation source point coordinates to create different scenarios and calculate them with the algorithm. According to the shown in Figure 1, we select (20,10), (20, 15), (20, 20), (20, 25), (20, 30), (20, 35), and (20, 40) as the coordinates of the evacuation source point. The calculation results are shown in Table 5. The evacuation time is the time when the evacuees leave the evacuation source point. All time units are in seconds.Table 5
The calculative results of different evacuation source point coordinates.
Evacuation sourcepoint coordinate
Evacuation routes
Evacuation time (s)
The highest cumulative concentration (mg/m3)
Program executiontime (s)
(20, 10)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
12816.29
19.2
(20, 15)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
2696.45
19.8
(20, 20)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
304.27
19.9
(20, 25)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
18.59
20.4
(20, 30)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.81
20.6
(20, 35)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.22
20.6
(20, 40)
1-2-4-8-14-22-26-27-281-2-7-8-14-22-26-27-281-3-5-9-13-17-18-23-24-28
63
0.21
20.3To analyze the influence of different evacuation source point coordinates on the result, we use the line chart to display the calculated results in Table5, as shown in Figure 8.Figure 8
The results of different evacuation source point coordinates.
(a)
The comprehensive comparison of every indicator
(b)
The highest cumulative concentration of different evacuation source point abscissaIn this case, the evacuation scale and the shortest path quantity are held constant, making the evacuation source point abscissa constant. Changing the evacuation source point ordinate, from Figure8(a), we can see that the evacuation time curve shows a horizontal trend. From Table 5, we can see that changing the evacuation source point coordinates has no effect on selecting the shortest routes. It can be seen, with different evacuation source point coordinates, that the number of evacuees who leave the evacuation source point is always the same at every moment, so the evacuation time is also the same. This suggests that changing the evacuation source point coordinate has no significant effect on the evacuation time. As the evacuation source point coordinates change, the program execution curve has a certain fluctuation, but the range is just 1.4 seconds. This suggests that changing the evacuation source point coordinate has no significant effect on program execution time. As the evacuation source point abscissa increases, the highest cumulative concentration curve declines rapidly towards zero. This suggests that changing the evacuation source point coordinates has a significant effect on the highest cumulative concentration of the whole evacuation process. When the evacuation source point abscissa is a constant, the farther the ordinate and the lower the concentration of hazardous chemicals at the evacuation source point. This makes the highest cumulative concentration lower in the whole evacuation process. From Figure 8, we can see that as the evacuation source point ordinates increase, the disparity of highest cumulative concentration between different evacuation source point abscissas decreases, and the speed of decrease is faster and faster. Combined with Figure 1, we can see that the influence distance of the hazardous chemicals leakage in the downwind direction is much longer than that in the vertical direction of the wind, so we can know that the influence of the evacuation source point ordinates on the highest cumulative concentration is more significant than that of the abscissa. When we choose the evacuation source point, we should increase the ordinate to make the source point further away from the leakage point in the longitudinal direction, thereby reducing the highest cumulative concentration in the whole evacuation process.


## 5. Conclusions

Emergency evacuation is one of the key steps of emergency rescue. It is of great significance to study the evacuation problem in order to decrease the harm of emergency accidents and disasters. However, different from other disasters and accidents, hazardous chemical accidents may diffuse their damage over a wide area. In this paper, considering the dynamic changes of hazardous chemicals’ concentration, we establish a mixed-integer programming model for the emergency evacuation of hazardous chemical accidents that aims to minimize the concentration of hazardous chemicals to which people are exposed during the entire evacuation process. Then, based on network optimization, we design the NODS algorithm that can calculate evacuation plans following hazardous chemical accidents. The numerical examples illustrate the feasibility and effectiveness of NODS algorithm. Finally, we explore the key factors affecting the evacuation process through different scenarios. In future research, we will further consider hazardous chemical accidents with multi-source points.



---

*Source: 4927649-2017-12-14.xml*