# A Robust Optimization Approach to Emergency Vehicle Scheduling

Authors: Xiao-Xia Rong; Yi Lu; Rui-Rui Yin; Jiang-Hua Zhang
Journal: Mathematical Problems in Engineering
 (2013)

---

## Abstract

The emergency vehicle scheduling problem is studied under the objective function to minimize the total transportation time with uncertain road travel time. Firstly, we build a stochastic programming model considering the constrained chance. Then, we analyze the model based on robust optimization method and get its equivalent set of uncertainty constraint, which has good mathematical properties with consideration of the robustness of solutions. Finally, we implement a numerical example to compare the results of robust optimization method and that of the particle swarm optimization algorithm. The case study shows that the proposed method achieves better performance on computational complexity and stability.


---

## Body

## 1. Introduction

Vehicle scheduling problem is a key problem in logistics management, which refers to organizing the appropriate route to make vehicles ordered. Under the satisfaction with certain constraints (such as goods demand and supply, delivery time, distance restriction, vehicle capacity constraint, etc.), it achieves a certain goal (such as minimum cost, minimum transportation time, minimum vehicle number, etc.). Emergency vehicle scheduling problem is a kind of special logistics activity to cope with unexpected disasters and accidents by vehicle scheduling for emergency supplies, personnel, and so on. From the 1970s, there are some methods proposed to solve the problem of disaster relief for its importance. Chaiken and Larson [1] provided a description of emergency service systems and designed a strategy to solve the corresponding operational problems. Based on their model, several improvements have been achieved for specific requirements and conditions in the disasters or incidents. Bozorgi-Amiri et al. [2] proposed the multiobjective stochastic programming model of emergency vehicle scheduling problem under uncertainty and constructed the three-layer model from the supplier to the distribution center according to the needs so that the final cost and satisfaction are optimal. Vitoriano et al. [3] developed a multicriteria optimization model for humanitarian aid distribution. This model is the core of a decision support system under development to assist organizations in charge of the distribution of humanitarian aid and is illustrated with a case study based on the 2010 Haiti catastrophic earthquake. Campbell et al. [4] specially focus on two alternative objective functions for the TSP and VRP: one that minimizes the maximum arrival time (minmax) and one that minimizes the average arrival time (minavg). Yuan and Wang [5] proposed two vehicle scheduling models of the emergency logistics management. The first model used the modified Dijkstra algorithm to solve related problems in order to minimize total travel time along a path. The second model used the ant colony algorithm to minimize the total travel time and the path complexity. Zhang et al. [6] proposed an emergency vehicle scheduling problem which considered the multiple resource constraints and possible secondary disasters and established a heuristic algorithm to solve this problem effectively.However, the above literatures rarely considered the scenario of the uncertainty of the vehicle travel time. In fact, when the unexpected disasters and accidents happened, the communication facilities could not be used as usual due to destruction of the nature force. As a result, we could not get the real road information. Besides, there are various sudden factors that cause the travel time on the road to be uncertain. At this time, the decision-makers just depend on the limited data to predict extent of the damage and the disaster range and make the emergency rescues in the shortest possible time.The robust optimization method is first proposed by Soyster in 1973 [7], which uses math set theory to describe the uncertain information of the parameters. It uses the large probability to avoid the decision deviation in bad condition, get the robust domain, and avoid the huge losses. This method can effectively avoid the instability of other algorithms, and it also has great application. For example, Ben-Tal et al. [8] applied the robust optimization method into the dynamic management of supply contract. Bertsimas and Thiele [9] applied the robust optimization method into the inventory-pricing model and inventory management. Bertsimas et al. [10] surveyed the primary research, both theoretical and applied, in the area of robust optimization. They highlighted applications of RO across a wide spectrum of domains, including finance, statistics, learning, and various areas of engineering. Therefore, this paper will use the robust optimization method to analyze the question of the emergency vehicle scheduling problem under the uncertainty of the travel time. And then, we will propose a new optimization algorithm to solve this problem.The structure of this paper is as follows. In Section2, we introduce the model of emergency vehicle routing problem in uncertain travel time. In Section 3, we give the robust counterpart of the uncertainty set and formulate a new equivalent determined set. Section 4 involves an example and the comparison between two algorithms. Section 5 contains a summary and conclusions.


## 2. Problem Description

In general, assuming that the situation of all roads is completely known, the network includes a material centeroandnaffected areas. We need to deliver supplies from material centeroto thenaffected nodes. The vehiclesmstart from the centeroand return to theoafter task is finished. Each node can only be visited once. The total transportation time is required to be minimum. Then the model is given as follows ((2a), (2b), (2c), (2d), and (2e)):
(1)minxijk∑k=1m∑j=0n∑i=0ntijkxijkhhhhhhhhhhhhhhhhhhhhhhhh(2a)s.t.∑i=1nqiyik≤Qk,k=1,…,m(2b)ff∑k=1mysk=1,s=1,…,n(2c)ffy0k=1,k=1,…,m(2d)ff∑i=1nxisk+∑j=1nxsjk=2ysk,k=1,…,m,s=1,…,n(2e)ffxijk={0,1},ysk={0,1},where (2a) is the load constraint for vehicle; (2b) represents the total number of all vehicles to access the nodes; each affected node should be visited one time. The formula of (2c) ensures that each vehicle starts from the materials center. (2d) shows that the vehicle goes into the node and then leaves from it, which can ensure the right order of the access nodes. Of which,qiis the demand at the nodei;Qkis the maximum weight load of thekth vehicle;tijkis the travel time on the arc(i,j)of thekth vehicle.xijkis 0-1 integer variables, where 1 means through the arc(i,j)and 0 is otherwise.yskis also 0-1 integer decision variables and 1 means thekth vehicle visits the nodesand 0 is otherwise.However, in the reality, when unexpected disasters occur, such as earthquake and floods, the road information becomes unknown and all kinds of unexpected factors make the travel time along the road in be determinable. At this time, we can only give the probability distribution of the travel time through the limited information and historical data. Thus, it creates a special kind of emergency vehicle scheduling problem. The road information is uncertain in the transport network. The vehicle starts from a rescue center and traverses all affected locations and then back to the center, which requires the minimum total expected travel time. In the model (2a), (2b), (2c), (2d), and (2e), the travel time in which the vehiclektravels on the arc(i,j)cannot be determined, so the travel time of thekth vehicle is uncertain, which is denoted byt~ijk. The objective function is adjusted to the mean oft~ijk. Adding constraints (4f),θis the effective time of the rescue andβis the small probability. This constraint means that ensuring the relief supplies arrives in time in greater probability. The stochastic programming model is as follows ((4a), (4b), (4c), (4d), (4e), and (4f)):
(3)minE(∑k=1m∑j=0n∑i=0nt~ijkxijk)hhhhhggggggghhhhhhh(4a)s.t.∑i=1nqiyik≤Qk,k=1,…,m(4b)fff∑k=1mysk=1,s=1,…,n(4c)fffy0k=1,k=1,…,m(4d)fff∑i=1nxisk+∑j=1nxsjk=2ysk,k=1,…,m,s=1,…,n(4e)fffxijk={0,1},ysk={0,1}(4f)fffP{∑i=0n∑j=0nt~ijkxijk≤θ}≥1-β,k=1,…,m.In this model, the chance constraint with uncertainties (4f) is the difficulty of solving the model. In order to describe it conveniently, we will use∑i=0n∑j=0nt~ijkxijk≤θ’s vector formt~x-θ≤0.


## 3. Robust Optimization Analysis

### 3.1. Related Knowledge

We present uncertainties of the datat~as affinely dependent on a set of independent random variablest~i,i=1,…,l, as
(5)t~=t0+∑i=1lΔtiz~i,z~i∈W,W=(-z_,z-),
wheret0is the nominal value of the data andΔtiis a direction of data perturbation. We callz~ithe primitive uncertainty, which has mean zero and support in(-z_i,z-i).In this paper, foru∈Rnis a vector, we denote the definition of vector norm and dual norm according to [11, 12]. We restrict the vector norm to be considered in an uncertainty set as follows:∥
u
∥
=
∥
|
u
|
∥and∥u∥≤∥u∥2∀u, where|u|is the vector with thejcomponent equal to|uj|for eachjfrom 1 tol. We call this an absolute norm. The dual norm∥·∥*is defined as∥u∥*=max∥x∥≤1u′x.We next show some basic properties, which are from [11, 12].Proposition 1.
For the absolute norm,(a)
one has(6)∥w∥*=∥|w|∥*;(b)
for allv,w, such that|v|≤|w|,∥v∥*≤∥w∥*;(c)
for allv,w, such that|v|≤|w|,∥v∥≤∥w∥;(d)∥t∥*≥∥t∥2,∀t
.Proposition 2.
For the programming as follows:(7)z*=maxa′v+b′w,s.t.∥v+w∥≤δv,w≥0.
The optimal value isz*=δ∥t∥*, wheretj=max{aj,bj,0},j=1,…,l.When random variables are incorporated in optimization models, operations are often cumbersome and computationally intractable. Instead of using complete distributional information, our aim is to obtain nontrivial probability bounds against constraint violation. Here, we denote the set of values associated with forward and backward deviations of random variablez~. Letz~be a random variable and letMz~(s)=E(exp(sz~))be its moment-generating function. We denote the set of values associated with forward deviations ofz~as follows:
(8)F(z~)={α:α≥0,Mz~-E(z~)(ϕα)≤exp(ϕ22),∀ϕ≥0}.Likewise, for backward deviations, the following set is defined:(9)B(z~)={α:α≥0,Mz~-E(z~)(-ϕα)≤exp(ϕ22),∀ϕ≥0}.Whenz~is symmetrically distributed around its mean, then we haveF(z~)=B(z~).Proposition 3.
Letx~andy~be two independent random variables with zero means, such thatpx~∈F(x~),qx~∈B(x~),py~∈F(y~), andqy~∈B(y~).(a)
Ifz~=x~+y~, then(pz~,qz~)=(px~2+py~2,qx~2+qy~2)satisfiespz~∈F(z~),qz~∈B(z~).(b)
P(x~>Ωpx~)≤exp(-Ω2/2)andP(x~<-Ωqx~)≤exp(-Ω2/2).We will analyze the primitive uncertaintyz~iin two cases. The first is only about the norm uncertainty set, which is computed in the short time. In the second case, we will discuss the constraint ofWin (5) at the base of the norm uncertainty set, which has the high complexity.


### 3.2. The Robust Analysis Based on Norm Uncertainty Set

Ifz~is symmetrically distributed, we use the symmetrically norm uncertainty set as follows:A
1
=
{
z
~
:
∥
z
~
∥
≤
Ω
}, whereΩ∈Ris a fixed number. Generally, ifz~is asymmetrically distributed, we may use the asymmetrically norm uncertainty set as follows:A
2
=
{
z
~
:
z
~
=
v
~
-
w
~
,
∥
P
-
1
v
~
+
Q
-
1
w
~
∥
≤
Ω
,
v
~
,
w
~
≥
0
,
v
~
,
w
~
∈
R
l
}, whereP=diag(p1,…,pl)andQ=diag(q1,…,ql),pi,qi>0. Particularly, whenP,Qare identity matrix, the setA2is equivalent to theA1. Then, (5) can be written as
(10)t~=t0+∑i=1lΔti(v~i-w~i).
And the constraintt~x-θ≤0is equivalent to
(11)(t0+∑i=1lΔti(v~i-w~i))x≤θ,with∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0.Theorem 4.
ForΩ=2ln(1/β),∃P=diag(p1,…,pl),Q=diag(q1,…,ql),u∈Rl, ifxsatisfies the robust counterpart as follows:
(12)t0′+Ω∥u∥*≤θuj≥pj(Δtj′x)uj≥-qj(Δtj′x)Ffj=1,…,l,
thenP{t~x>θ}≤βholds.Proof.
We first express how to get the set (12). Let
(13)φ=(φ1,φ2,…,φl)=(xΔt1,…,xΔtl).
Thus, (11) can be rewritten as
(14)t0x+∑i=1lφi(v~i-w~i)≤θ,∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0⟺t0x+max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{φ(v~-w~)}≤θ.
According to Proposition 2, we have
(15)t0x+Ω∥u∥*≤θ,
where
(16)u=(u1,u2,…,ul),uj=max{pjφj,-qjφj,0}=max{pjφj,-qjφj}.
So a new set which is called robust counterpart is got:
(17)t0+Ω∥u∥*≤θuj≥pj(Δtjx)uj≥-qj(Δtjx)u∈Rl,uj≥0,j=(1,…,l).
Next, we will prove thatxsatisfies (12) if and only ifP(t~x>θ)<βholds.
From (16), we have
(18)P{t~x>θ}=P(t0x+∑i=1lz~(Δtix)>θ)=P(φz~>θ-t0x)≤P(φz~>Ω∥u∥*).
Thus,
(19)P(∑i=0n∑j=0nt~ijxij>θ)≤P(φz~>Ω∥u∥*)≤P(φz~>Ω∥u∥2).
Let the elements ofP,Qform the setsP(z~),Q(z~), according to Proposition 3; then∥u∥2∈P(φ′z~). At the same time, whenΩ=2ln(1/β)holds, we have
(20)P(∑i=0n∑j=0nt~ijxij>θ)≤exp(-Ω22)=β.In this case,A1,A2are discussed only considering the absolute normΩ. The distribution of random variable is not computed, such as the constraint ofz~i∈W,W=(-z_,z-)in (5), so there is some deviation between the counterpart set and the original set. Further, we consider the distributionWof the random variable based on norm uncertainty set in (6).


### 3.3. The Robust Analysis Based on Norm Uncertainty Set andW

In this case, we consider all probable values oft~, including the worst-caseW=(-z_,z-). The constraint ofz~in norm uncertainty setA1may be written as
(21)B1={ξ~:∥ξ~∥≤Ω,-z_≤ξ~≤z-}.Corresponding toA2, we have
(22)B2={ξ~:ξ~=v~i-w~i,∥P-1v~+Q-1w~∥≤Ω,g-z_≤v~i-w~i≤z-,v~,w~≥0∥P-1v~+Q-1w~∥},
whereP=diag(p1,…,pl)and likewiseQ=diag(q1,…,ql)withpi,qi>0,i=1,…,l.Then,(23)t~=t0+∑i=1lΔti(v~i-w~i).So,t~x-θ≤0can be expressed as follows:
(24)(t0+∑i=1lΔti(v~i-w~i))x≤θ,with∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0,-z_≤v~-w~≤z-.
Similarly, we have the following conclusions.Theorem 5.
ForP=diag(p1,…,pl),Q=diag(q1,…,ql),u∈Rl, andΩ=2ln(1/β), ifxsatisfies the robust counterpart as follows:
(25)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,uj≥pj(Δtjx-rj+sj),uj≥-qj(Δtjx-rj+sj),rj,sj≥0,∀j=1,2,…,l,
thenP{t~x>θ}≤βholds.Proof.
We first express how to get the set (25).
Let(26)φ=(φ1,φ2,…,φl)=(xΔt1,xΔt2,…,xΔtl).
Thus, (24) can be rewritten as(27)t0x+∑i=1lΔti(v~i-w~i)′φi≤θ,∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0,-z_≤v~-w~≤z-⟺t0x+max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,-z_≤v~-w~≤z-,v~,w~≥0}{(v~-w~)y}≤θ⟺t0x+minr,s>0{max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{(v~-w~)φ+r(z--v~+w~)+s(z_+v~-w~)}}≤θ⟺t0x+minr,s>0{max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{(φ-r+s)v~-(φ-r+s)w~+rz-+sz_}}≤θ.
According to Proposition2, we have
(28)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,
where
(29)uj(r,s)=max{pj(φj-rj+sj),-qj(φj-rj+sj),0}=max{pj(φj-rj+sj),-qj(φj-rj+sj)}.
So, a new set which is called robust counterpart is got:(30)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,uj≥pj(Δtjx-rj+sj),uj≥-qj(Δtjx-rj+sj),rj,sj≥0,∀j=1,2,…,l.
Next, we will prove thatΩ=2ln(1/β)if and only ifP(t~x>θ)<β.
From (28),
(31)P{t~x>θ}=P(φξ~>θ-t0x)≤P(φξ~>minr,s≥0{Ω∥u∥*+rz-+sz_}).
According to Proposition 1,
(32)P(∑i=0n∑j=0nt~ijxij>θ)≤P(φξ~>minr,s≥0{Ω∥u∥*+rz-+sz_})≤P(φξ~>minr,s≥0{Ω∥u∥2+rz-+sz_}).
Letu*=u(r*,s*)be the optimal solution ofminr,s≥0{Ω∥u∥2+rz-+sz_}, and, considering the inequalityξ~∈(-z_,z-), then
(33)P(φξ~>minr,s≥0{Ω∥u∥2+rz-+sz_})=P(φξ~>Ω∥u*∥2+r*z-+s*z_)≤P((φ-r*+s*)ξ~>Ω∥u*∥2).
Let the elements ofP,Qform the setsP(z~),Q(z~), according to Proposition 3; we have∥u*∥2∈P((φ-r*+s*)ξ~). So forΩ=2ln(1/β), the following inequality holds
(34)P(∑i=0n∑j=0nt~ijxij>θ)≤exp(-Ω22)=β.


## 3.1. Related Knowledge

We present uncertainties of the datat~as affinely dependent on a set of independent random variablest~i,i=1,…,l, as
(5)t~=t0+∑i=1lΔtiz~i,z~i∈W,W=(-z_,z-),
wheret0is the nominal value of the data andΔtiis a direction of data perturbation. We callz~ithe primitive uncertainty, which has mean zero and support in(-z_i,z-i).In this paper, foru∈Rnis a vector, we denote the definition of vector norm and dual norm according to [11, 12]. We restrict the vector norm to be considered in an uncertainty set as follows:∥
u
∥
=
∥
|
u
|
∥and∥u∥≤∥u∥2∀u, where|u|is the vector with thejcomponent equal to|uj|for eachjfrom 1 tol. We call this an absolute norm. The dual norm∥·∥*is defined as∥u∥*=max∥x∥≤1u′x.We next show some basic properties, which are from [11, 12].Proposition 1.
For the absolute norm,(a)
one has(6)∥w∥*=∥|w|∥*;(b)
for allv,w, such that|v|≤|w|,∥v∥*≤∥w∥*;(c)
for allv,w, such that|v|≤|w|,∥v∥≤∥w∥;(d)∥t∥*≥∥t∥2,∀t
.Proposition 2.
For the programming as follows:(7)z*=maxa′v+b′w,s.t.∥v+w∥≤δv,w≥0.
The optimal value isz*=δ∥t∥*, wheretj=max{aj,bj,0},j=1,…,l.When random variables are incorporated in optimization models, operations are often cumbersome and computationally intractable. Instead of using complete distributional information, our aim is to obtain nontrivial probability bounds against constraint violation. Here, we denote the set of values associated with forward and backward deviations of random variablez~. Letz~be a random variable and letMz~(s)=E(exp(sz~))be its moment-generating function. We denote the set of values associated with forward deviations ofz~as follows:
(8)F(z~)={α:α≥0,Mz~-E(z~)(ϕα)≤exp(ϕ22),∀ϕ≥0}.Likewise, for backward deviations, the following set is defined:(9)B(z~)={α:α≥0,Mz~-E(z~)(-ϕα)≤exp(ϕ22),∀ϕ≥0}.Whenz~is symmetrically distributed around its mean, then we haveF(z~)=B(z~).Proposition 3.
Letx~andy~be two independent random variables with zero means, such thatpx~∈F(x~),qx~∈B(x~),py~∈F(y~), andqy~∈B(y~).(a)
Ifz~=x~+y~, then(pz~,qz~)=(px~2+py~2,qx~2+qy~2)satisfiespz~∈F(z~),qz~∈B(z~).(b)
P(x~>Ωpx~)≤exp(-Ω2/2)andP(x~<-Ωqx~)≤exp(-Ω2/2).We will analyze the primitive uncertaintyz~iin two cases. The first is only about the norm uncertainty set, which is computed in the short time. In the second case, we will discuss the constraint ofWin (5) at the base of the norm uncertainty set, which has the high complexity.


## 3.2. The Robust Analysis Based on Norm Uncertainty Set

Ifz~is symmetrically distributed, we use the symmetrically norm uncertainty set as follows:A
1
=
{
z
~
:
∥
z
~
∥
≤
Ω
}, whereΩ∈Ris a fixed number. Generally, ifz~is asymmetrically distributed, we may use the asymmetrically norm uncertainty set as follows:A
2
=
{
z
~
:
z
~
=
v
~
-
w
~
,
∥
P
-
1
v
~
+
Q
-
1
w
~
∥
≤
Ω
,
v
~
,
w
~
≥
0
,
v
~
,
w
~
∈
R
l
}, whereP=diag(p1,…,pl)andQ=diag(q1,…,ql),pi,qi>0. Particularly, whenP,Qare identity matrix, the setA2is equivalent to theA1. Then, (5) can be written as
(10)t~=t0+∑i=1lΔti(v~i-w~i).
And the constraintt~x-θ≤0is equivalent to
(11)(t0+∑i=1lΔti(v~i-w~i))x≤θ,with∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0.Theorem 4.
ForΩ=2ln(1/β),∃P=diag(p1,…,pl),Q=diag(q1,…,ql),u∈Rl, ifxsatisfies the robust counterpart as follows:
(12)t0′+Ω∥u∥*≤θuj≥pj(Δtj′x)uj≥-qj(Δtj′x)Ffj=1,…,l,
thenP{t~x>θ}≤βholds.Proof.
We first express how to get the set (12). Let
(13)φ=(φ1,φ2,…,φl)=(xΔt1,…,xΔtl).
Thus, (11) can be rewritten as
(14)t0x+∑i=1lφi(v~i-w~i)≤θ,∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0⟺t0x+max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{φ(v~-w~)}≤θ.
According to Proposition 2, we have
(15)t0x+Ω∥u∥*≤θ,
where
(16)u=(u1,u2,…,ul),uj=max{pjφj,-qjφj,0}=max{pjφj,-qjφj}.
So a new set which is called robust counterpart is got:
(17)t0+Ω∥u∥*≤θuj≥pj(Δtjx)uj≥-qj(Δtjx)u∈Rl,uj≥0,j=(1,…,l).
Next, we will prove thatxsatisfies (12) if and only ifP(t~x>θ)<βholds.
From (16), we have
(18)P{t~x>θ}=P(t0x+∑i=1lz~(Δtix)>θ)=P(φz~>θ-t0x)≤P(φz~>Ω∥u∥*).
Thus,
(19)P(∑i=0n∑j=0nt~ijxij>θ)≤P(φz~>Ω∥u∥*)≤P(φz~>Ω∥u∥2).
Let the elements ofP,Qform the setsP(z~),Q(z~), according to Proposition 3; then∥u∥2∈P(φ′z~). At the same time, whenΩ=2ln(1/β)holds, we have
(20)P(∑i=0n∑j=0nt~ijxij>θ)≤exp(-Ω22)=β.In this case,A1,A2are discussed only considering the absolute normΩ. The distribution of random variable is not computed, such as the constraint ofz~i∈W,W=(-z_,z-)in (5), so there is some deviation between the counterpart set and the original set. Further, we consider the distributionWof the random variable based on norm uncertainty set in (6).


## 3.3. The Robust Analysis Based on Norm Uncertainty Set andW

In this case, we consider all probable values oft~, including the worst-caseW=(-z_,z-). The constraint ofz~in norm uncertainty setA1may be written as
(21)B1={ξ~:∥ξ~∥≤Ω,-z_≤ξ~≤z-}.Corresponding toA2, we have
(22)B2={ξ~:ξ~=v~i-w~i,∥P-1v~+Q-1w~∥≤Ω,g-z_≤v~i-w~i≤z-,v~,w~≥0∥P-1v~+Q-1w~∥},
whereP=diag(p1,…,pl)and likewiseQ=diag(q1,…,ql)withpi,qi>0,i=1,…,l.Then,(23)t~=t0+∑i=1lΔti(v~i-w~i).So,t~x-θ≤0can be expressed as follows:
(24)(t0+∑i=1lΔti(v~i-w~i))x≤θ,with∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0,-z_≤v~-w~≤z-.
Similarly, we have the following conclusions.Theorem 5.
ForP=diag(p1,…,pl),Q=diag(q1,…,ql),u∈Rl, andΩ=2ln(1/β), ifxsatisfies the robust counterpart as follows:
(25)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,uj≥pj(Δtjx-rj+sj),uj≥-qj(Δtjx-rj+sj),rj,sj≥0,∀j=1,2,…,l,
thenP{t~x>θ}≤βholds.Proof.
We first express how to get the set (25).
Let(26)φ=(φ1,φ2,…,φl)=(xΔt1,xΔt2,…,xΔtl).
Thus, (24) can be rewritten as(27)t0x+∑i=1lΔti(v~i-w~i)′φi≤θ,∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0,-z_≤v~-w~≤z-⟺t0x+max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,-z_≤v~-w~≤z-,v~,w~≥0}{(v~-w~)y}≤θ⟺t0x+minr,s>0{max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{(v~-w~)φ+r(z--v~+w~)+s(z_+v~-w~)}}≤θ⟺t0x+minr,s>0{max{v~,w~:∥P-1v~+Q-1w~∥≤Ω,v~,w~≥0}{(φ-r+s)v~-(φ-r+s)w~+rz-+sz_}}≤θ.
According to Proposition2, we have
(28)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,
where
(29)uj(r,s)=max{pj(φj-rj+sj),-qj(φj-rj+sj),0}=max{pj(φj-rj+sj),-qj(φj-rj+sj)}.
So, a new set which is called robust counterpart is got:(30)t0x+minr,s≥0{Ω∥u∥*+rz-+sz_}≤θ,uj≥pj(Δtjx-rj+sj),uj≥-qj(Δtjx-rj+sj),rj,sj≥0,∀j=1,2,…,l.
Next, we will prove thatΩ=2ln(1/β)if and only ifP(t~x>θ)<β.
From (28),
(31)P{t~x>θ}=P(φξ~>θ-t0x)≤P(φξ~>minr,s≥0{Ω∥u∥*+rz-+sz_}).
According to Proposition 1,
(32)P(∑i=0n∑j=0nt~ijxij>θ)≤P(φξ~>minr,s≥0{Ω∥u∥*+rz-+sz_})≤P(φξ~>minr,s≥0{Ω∥u∥2+rz-+sz_}).
Letu*=u(r*,s*)be the optimal solution ofminr,s≥0{Ω∥u∥2+rz-+sz_}, and, considering the inequalityξ~∈(-z_,z-), then
(33)P(φξ~>minr,s≥0{Ω∥u∥2+rz-+sz_})=P(φξ~>Ω∥u*∥2+r*z-+s*z_)≤P((φ-r*+s*)ξ~>Ω∥u*∥2).
Let the elements ofP,Qform the setsP(z~),Q(z~), according to Proposition 3; we have∥u*∥2∈P((φ-r*+s*)ξ~). So forΩ=2ln(1/β), the following inequality holds
(34)P(∑i=0n∑j=0nt~ijxij>θ)≤exp(-Ω22)=β.


## 4. Numerical Experiments

Rescue centerois responsible for the rescue of eight affected areas. The vehicle began in the rescue center and back to rescue center. The problem is how to plan the vehicle running routes, making the recipient regions all visited once. In order to obtain the maximum effectiveness of the rescue, the affected people should be rescued as soon as possible. Particularly, because the disaster has just occurred, the traffic situation is not clear. We know that the largest transit time may be 1.25 times of the normal. The simulation data includes the coordinate of rescue center and eight affected areas and the demands of emergency resources in affected areas; they are shown in Table 1.Table 1
The simulation data.
Center
o
1
2
3
4
5
6
7
8
Coordinate
(18, 20)
(22, 26)
(33, 21)
(38, 20)
(18, 5)
(18, 30)
(13, 10)
(28, 19)
(17, 35)
Demand
0
2.9
2.6
3.7
4.2
3.4
4.7
3.5
3.2The rescue center can be assembled with three available vehicles. The carrying capacity is 10 units. Under normal circumstances, the average transport speed is 40 km/h and the normal travel time is given astk*. Based on historical data, it can be informed that the estimated average travel time is1.25tk*fluctuating around the average time. Thus, we get the chance constrained model as follows:
(35)minE(∑k=13∑j=08∑i=08t~ijkxijk)s.t.∑i=18qiyik≤Qk,k=1,2,3s.t.∑k=13ysk=1,s=1,2,…,8s.t.y0k=1,k=1,2,3s.t.∑i=08xisk+∑j=08xsjk=2yks,k=1,2,3s.t.P{∑i=0n∑j=0nt~ijkxijk>θ}≤β,k=1,2,3s.t.xijk={0,1},ysk={0,1},
where the demand vector isq=(2.9,2.6,3.7,4.2,3.4,4.7,3.5,and3.2)′.In this section, we will solve the model (35) using the robust analysis in Section 3.2 and hybrid particle swarm method and then carry on the comparison and analysis of results.


### 4.1. The Robust Analysis

For only considering the time delay of the road, we present uncertainties of the datat~kas the following representation:
(36)t~k=tk0+tkz~k.Base on the experiment above,t~k∈(tk*,1.5tk*), so lettk0=1.25tk*. Here, we choose the symmetrically norm uncertaintyA1to describet~k; that is,
(37)t~k=tk0+Δtkz~k,z~k∈(-0.25tk*,0.25tk*).
Then,Δtk=tk*; we can use∥·∥*=∥·∥=|·|; thus,
(38)∥z~k∥=|z~k|≤0.25.Letβ=0.01andΩ=0.1; according to Theorem 4, we get the equivalent form of the chance constrained in (35) as follows:
(39)1.25tk*′x+0.1h≤θ,∃u,h∈R,|u|≤h,u≥(t*′x),u≥-(t*′x).
So, model (35) is transformed into the following model,
(40)minE(∑k=1m∑j=0n∑i=0nt~ijkxijk)s.t.∑i=18qiyik≤Qk,k=1,…,3s.t.∑k=13ysk=1,s=1,…,8s.t.y0k=1,k=1,2,3s.t.∑i=18xisk+∑j=18xsjk=2ysk,k=1,2,3,s=0,1,…,8s.t.1.25tk*x+0.1h≤θ,k=1,2,3s.t.h-1.25tk*x≥0,k=1,2,3s.t.xijk={0,1},ysk={0,1}.This model is a 0-1 integer programming with the variablesxijk,ysk, andh. We get the optimal value 5.3242 and the optimal rescue route 0- 1- 5 -8- 0 -6 -4 -0- 2 -3 -7- 0 by LINGO solver. The optimal rescue route is shown in Figure 1.Figure 1
The optimal rescue route in model (40).


### 4.2. Hybrid Particle Swarm Optimization Algorithm

Particle swarm optimization (PSO) algorithm is a swarm intelligence algorithm, originating from the simulation of birds flocks looking for food and has been paid attention and researched wildly. In this paper, we use global and local hybrid particle swarm algorithm to solve the stochastic optimization problem, where the selected update formula for the velocity vector is as follows:(41)vh=α*vq+(1-α)*vl,α∈(0,1).In order to obtain better experimental effect, momentum factorwin a nonlinear concave function declined.According to model (5), we select the fitness function as follows:
(42)E(∑k=1m∑j=0n∑i=0nt~ijkxijk)+M1∑k=1mmax{qiyi-Qk,0}+M2∑k=1mmax{E(∑j=0n∑i=0nt~ijkxijk)-θ,0}.Encoded mode uses real number coding way, which transforms 0-1 problem of the vehicle routing problem into the continuous problem.The process of algorithm is as follows.Step 1.
Population initialization and calculation of the fitness of initial.Step 2.
Update formula using the particle velocity and the location of the particledisxd(k+1)=xd(k)+vd(k+1).Step 3.
Calculate and record the latest fitness of the particle population.Step 4.
Determine whether accord with the termination condition or go back Step2.


### 4.3. The Analysis of the Operation Results

Figure2 shows the trajectory of fitness function of the group size of 20 particles in an iterative 1000 times. It is observed that it has better convergence for a single test.Figure 2
The convergence graph of a single test of fitness function.Figure3 is the record of 50-time tests results, and each test has 1000-time iterations. Because the random time is different, so the computation time is 2.5 hours.Figure 3
The convergence graph of 50-time tests of fitness function.Table2 is the result of the 50-time tests and group size of each test is 20. In the 50-time tests, the first 4 lines presented repeated solutions and the fifth line is the result only appearing once. So this table shows that the hybrid particle swarm algorithm results have strong randomness and stability is obviously insufficient.Table 2
Result of PSO.
The optimal solution
The results of the path
The number of times
x
1
=
0 2 5 3 0 4 6 0 1 7 8 0
6
x
2
=
0 1 3 8 0 4 6 0 2 5 7 0
4
x
3
=
0 2 3 7 0 4 6 0 1 5 8 0
13
x
4
=
0 8 7 2 0 6 4 0 5 3 1 0
5
Other
22In the result of the hybrid particle swarm algorithm, as we can see from Table2, the optimal solutionx3appears most frequently. The result is identical with robust optimization, which illustrates feasibility of the robust optimization methods, and compared with the hybrid particle swarm algorithm, the robust optimization methods have better stability. In this case, the hybrid particle swarm optimization algorithm occupies a large space, and it requires multiple iterations and several tests in order to obtain a stable solution, which is far greater than the robust method. In a word, from the perspective of the computational complexity and the stability of solution results, robust optimization method has good application value.


## 4.1. The Robust Analysis

For only considering the time delay of the road, we present uncertainties of the datat~kas the following representation:
(36)t~k=tk0+tkz~k.Base on the experiment above,t~k∈(tk*,1.5tk*), so lettk0=1.25tk*. Here, we choose the symmetrically norm uncertaintyA1to describet~k; that is,
(37)t~k=tk0+Δtkz~k,z~k∈(-0.25tk*,0.25tk*).
Then,Δtk=tk*; we can use∥·∥*=∥·∥=|·|; thus,
(38)∥z~k∥=|z~k|≤0.25.Letβ=0.01andΩ=0.1; according to Theorem 4, we get the equivalent form of the chance constrained in (35) as follows:
(39)1.25tk*′x+0.1h≤θ,∃u,h∈R,|u|≤h,u≥(t*′x),u≥-(t*′x).
So, model (35) is transformed into the following model,
(40)minE(∑k=1m∑j=0n∑i=0nt~ijkxijk)s.t.∑i=18qiyik≤Qk,k=1,…,3s.t.∑k=13ysk=1,s=1,…,8s.t.y0k=1,k=1,2,3s.t.∑i=18xisk+∑j=18xsjk=2ysk,k=1,2,3,s=0,1,…,8s.t.1.25tk*x+0.1h≤θ,k=1,2,3s.t.h-1.25tk*x≥0,k=1,2,3s.t.xijk={0,1},ysk={0,1}.This model is a 0-1 integer programming with the variablesxijk,ysk, andh. We get the optimal value 5.3242 and the optimal rescue route 0- 1- 5 -8- 0 -6 -4 -0- 2 -3 -7- 0 by LINGO solver. The optimal rescue route is shown in Figure 1.Figure 1
The optimal rescue route in model (40).


## 4.2. Hybrid Particle Swarm Optimization Algorithm

Particle swarm optimization (PSO) algorithm is a swarm intelligence algorithm, originating from the simulation of birds flocks looking for food and has been paid attention and researched wildly. In this paper, we use global and local hybrid particle swarm algorithm to solve the stochastic optimization problem, where the selected update formula for the velocity vector is as follows:(41)vh=α*vq+(1-α)*vl,α∈(0,1).In order to obtain better experimental effect, momentum factorwin a nonlinear concave function declined.According to model (5), we select the fitness function as follows:
(42)E(∑k=1m∑j=0n∑i=0nt~ijkxijk)+M1∑k=1mmax{qiyi-Qk,0}+M2∑k=1mmax{E(∑j=0n∑i=0nt~ijkxijk)-θ,0}.Encoded mode uses real number coding way, which transforms 0-1 problem of the vehicle routing problem into the continuous problem.The process of algorithm is as follows.Step 1.
Population initialization and calculation of the fitness of initial.Step 2.
Update formula using the particle velocity and the location of the particledisxd(k+1)=xd(k)+vd(k+1).Step 3.
Calculate and record the latest fitness of the particle population.Step 4.
Determine whether accord with the termination condition or go back Step2.


## 4.3. The Analysis of the Operation Results

Figure2 shows the trajectory of fitness function of the group size of 20 particles in an iterative 1000 times. It is observed that it has better convergence for a single test.Figure 2
The convergence graph of a single test of fitness function.Figure3 is the record of 50-time tests results, and each test has 1000-time iterations. Because the random time is different, so the computation time is 2.5 hours.Figure 3
The convergence graph of 50-time tests of fitness function.Table2 is the result of the 50-time tests and group size of each test is 20. In the 50-time tests, the first 4 lines presented repeated solutions and the fifth line is the result only appearing once. So this table shows that the hybrid particle swarm algorithm results have strong randomness and stability is obviously insufficient.Table 2
Result of PSO.
The optimal solution
The results of the path
The number of times
x
1
=
0 2 5 3 0 4 6 0 1 7 8 0
6
x
2
=
0 1 3 8 0 4 6 0 2 5 7 0
4
x
3
=
0 2 3 7 0 4 6 0 1 5 8 0
13
x
4
=
0 8 7 2 0 6 4 0 5 3 1 0
5
Other
22In the result of the hybrid particle swarm algorithm, as we can see from Table2, the optimal solutionx3appears most frequently. The result is identical with robust optimization, which illustrates feasibility of the robust optimization methods, and compared with the hybrid particle swarm algorithm, the robust optimization methods have better stability. In this case, the hybrid particle swarm optimization algorithm occupies a large space, and it requires multiple iterations and several tests in order to obtain a stable solution, which is far greater than the robust method. In a word, from the perspective of the computational complexity and the stability of solution results, robust optimization method has good application value.


## 5. Conclusions

It is a key research content of emergency management of disasters. It is also a key step of emergency rescue and assistance. This research will provide prompt and effective guidance for scientific decision making and support in disaster management. This paper studied the emergency vehicle scheduling problem under the uncertainty of travel time and proposed the stochastic programming model, and objective constraint is the minimum total transportation time. Considering the chance constraint in the model, we used the robust optimization method to analyze this model from the model’s robustness and get the equivalence set of the uncertainty constraint. Finally, this paper compared this method with PSO algorithm through a numerical example. The numerical simulation provides evidence for its effectiveness and efficiency. The results showed that the robust optimization method has better effect on computational complexity and stability. In the future, considering specific case, such as “Wenchuan” earthquake, we will study the application of this method to routing vehicles with variable route times due to traffic or other uncertainties.



---

*Source: 848312-2013-11-27.xml*