# Vehicle Routing Problem with Time Windows and Simultaneous Delivery and Pick-Up Service Based on MCPSO

Authors: Xiaobing Gan; Yan Wang; Shuhai Li; Ben Niu
Journal: Mathematical Problems in Engineering
 (2012)

---

## Abstract

This paper considers two additional factors of the widely researched vehicle routing problem with time windows (VRPTW). The two factors, which are very common characteristics in realworld, are uncertain number of vehicles and simultaneous delivery and pick-up service. Using minimization of the total transport costs as the objective of the extension VRPTW, a mathematic model is constructed. To solve the problem, an efficient multiswarm cooperative particle swarm optimization (MCPSO) algorithm is applied. And a new encoding method is proposed for the extension VRPTW. Finally, comparing with genetic algorithm (GA) and particle swarm optimization (PSO) algorithm, the MCPSO algorithm performs best for solving this problem.


---

## Body

## 1. Introduction

Vehicle routing problem with time windows (VRPTW) is an important issue in logistics system which has been researched widely in recent years. The problem can be described as choosing routes for limited number of vehicles to serve a group of customers in the time windows. Each vehicle has a limited capacity. It starts from the depot and terminates at the depot. Each customer should be served exactly once. The objective of the VRPTW is to minimize the total transport costs. Many researchers have contribute to the problem. In 1981, Schrage [1] proposed vehicle routing and scheduling problem with time window constraints as an important area for progress in handling realistic complications and generalizations of the basic routing model. Feng [2] divide the VRPTW with delivery and pickup into five categories. Wang and Lang [3] proposes multiperiod vehicle routing problem with recurring dynamic time windows.However, in reality, two factors need to be considered, the uncertain number of vehicles and simultaneous delivery and pick-up. This situation is very common in transport activity, because reduction of the number of the vehicles can save costs, while simultaneous delivery and pick-up can improve transport efficient. VRPTW with simultaneous delivery and pick-up service (VRPTW-SDP) is an extension of VRPTW, where customers require simultaneous delivery and pick-up.Many researchers have made great effort in solving the VRPTW in recent years. Fisher [4] proposed k-trees to solve the VRP and VRPTW. Kolen et al. [5] described a branch-and-bound method. However, those classical approaches less efficient in solving complex problems. Thus evolutionary algorithm (EA) is proposed to deal with optimization problems. Among those algorithms, particle swarm optimization (PSO) has turned out to be an efficient algorithm in dealing with many complex optimization problems. While PSO also has some disadvantages. It sometimes immerses the local optimal value, thus the accuracy is limited. Then the multiswarm cooperative particle swarm optimization (MCPSO) was proposed as an improved PSO in [6]. It takes a multi-swarm cooperative evolutionary strategy where the master swarms change their particles based on their own knowledge and the knowledge of the particles in the slave swarms, while the slave swarms carry out PSO independently.In this paper, we study a case of the VRPTW both with uncertain number of vehicles and simultaneous delivery and pickup service. The objective of the proposed VRPTW-SDP is minimizing the transport costs. In order to solve the VRPTW-SDP, we set the objective of VRPTW-SDP as the fitness function of MCPSO. In MCPSO, each particle should contain two aspects of the customers, the order of served by which vehicle. Thus two dimensions encoding methods in MCPSO are proposed and used for the VRPTW-SDP.The rest of this paper is arranged as follows. Section2 describes the VRPTW-SDP, Section 3 describes the MCPSO for the proposed problem, Section 4 presents experiment study, and Section 5 of the paper contains the conclusion.


## 2. Description of VRPTW-SDP

Vehicle routing problem with time windows (VRPTW) can be defined as choosing routes for limited number of vehicles to serve a group of customers in the time windows. Each vehicle has a limited capacity. It starts from the depot and terminates at the depot. Each customer should be served exactly once. If vehicles arrive before the time window “opens” or after the time window “closes,” there will be waiting cost and late cost. Assume that there areK vehicles in depot 0. N customers are waiting to be served and each of the customers has a demands of gi (i=1,2,…,N) units. The distance between customer i and j is dij. Each vehicle has a capacity of qk(k=1,2,…,K) units. That is, the total demands of customers served by each vehicle cannot exceed qk units. Therefore, the vehicle has to periodically returned to the depot for reloading or a new vehicle needed to be arranged for delivery. Besides, each customer must be visited once and only once by exactly one vehicle. si represents service time needed by customer i. Therefore, the vehicle has to stay at the location of customer i for a time interval at least si (s0=0 is associated with the depot 0) for service. A time window [ETi,LTi] is considered. Therefore, if a vehicle arrives at customer i before ETi, it has to wait until the beginning of the time window to serve the customer. Thus there is a cost e for waiting. On the other hand, if a vehicle cannot arrive at i before LTi, there will be a cost f for late. The velocity of each vehicle is vk. ti represents the moment when the vehicle arrives i from the depot. And the unit freight of each vehicle is Ck. This paper considers the condition that vehicles of the depot are the same, that is, the vehicles have the same velocity v, the same capacity q and the same unit freight C.Define variable(2.1)xijk={1ifthevehiclektravelsfromitoj,0else.The goal of VRPTW is(2.2)min⁡⁡z=∑i=0N∑j=0N∑k=1KC*xijk*dij+∑i=1Nmax⁡{e*(ETi-ti);0;f*(ti-LTi)},
where
(2.3)tijk=∑xijk(ti+dijv+si)(t0=0,s0=0),s
.
t.
(2.4)∑j=1N∑k=1Kxjik=∑j=1N∑k=1Kxijk=K(i=0),(2.5)∑j=0N∑k=1Kxijk=1(i∈N),(2.6)∑i=0N∑k=1Kxijk=1(i∈N),(2.7)∑j=1Nxijk=∑j=1Nxjik=1(i=0,k∈K),(2.8)∑i=0N∑j=0Nxijk*gi≤q.In the model, formula∑i=1Nmax{e*(ETi-ti);0;f*(ti-LTi)} defines the time window constraint, where ti is the time when the vehicle arrives at customer i, ETi-ti is the waiting time of a vehicle at the customer i. From formula (2.3), si is the service time, and tijk  is the traveling time of vehicle k between customer i and j. When i=0, the vehicle is at the depot, s0=t0=0.Constraint (2.4) represents the number of the vehicles which start from the depot and go back to depot is k. Constraint (2.5) and constraint (2.6) stand for each customer can only be served by one vehicle. Constraint (2.7) means all the vehicles which start from the depot go back to the depot. Constraint (2.8) denotes the quantity of goods that each vehicle carries could not exceed the capacity q.However, two extra aspects are involved. Firstly, what if the number of the vehicles is uncertain? Zhang et al. [7] proposed that the number of the vehicles can be calculate by  ∑qi/Q, where qi is the demands of the customer i, Q is the vehicle’s capacity. This is quite simple, but it may not meet the requirement of the time window. When the number of the vehicles is uncertain, it is not necessary to use all of the k vehicles. The reduction of the vehicles may contribute to reducing the total costs. Consider the constraint (2.4), where k can be replaced by m(m≤k). Secondly, simultaneous delivery and pickup is involved to the vehicles. Under this condition, two variables are related to the customer, demands (delivery) gi and supplies (pickup)  pi.Define parameters:C  is unit freight of the vehicle, K is the total number of the vehicles, N  is the total number of the customer, dij represents the distance between customer i and j, v is the velocity of each vehicle, qi is the demands of customer i, pi is the supplies of customer i, si is the service time of customer i, [ETi,LTi] represents the time window, e is the cost for waiting, f is the cost for late, q is the specified capacity of the vehicle, ti is the time when the vehicle arrives customer i, yijk represents the capacity of the vehicle k when it travels from customer i to j, zijk represents the pick-ups of the vehicle k from customer i when it travels from customer i to j,
(2.9)xijk={1ifthevehiclektravelsfromitoj,0else.When the vehicle starts from the depot, its total cargo quantity should be equal to the sum of all the customers’ demands, who are on the route of the vehicle travels along, so(2.10)∑j=1Nyojk=∑i=0N∑j=0Nxijk*gi(k∈K,i≠j).
While, when the vehicle k returns to the depot, the total pick-ups quantities are the sum of each customer’s supplies. That is
(2.11)∑i=1Nzi0k=∑i=0N∑j=0Nxijk*pi(k∈K,i≠j).
Then the model is as follows:
(2.12)min⁡⁡z∑i=0N∑j=0N∑k=1KC*xijk∑i=1Nmax⁡⁡{e*(ai-ti);0;f*(ti-bi)}+m*c,s
.
t.
(2.13)tijk=xijk(ti+dijv),(2.14)∑j=1Nxijk=∑j=1Nxjik≤1,(2.15)∑i=0N∑k=1Kxijk=1,(2.16)∑j=0N∑k=1Kxijk=1,(2.17)∑j=1Nyojk=∑i=0N∑j=0Nxijk*gi(k∈K,i≠j),(2.18)∑i=1Nzi0k=∑i=0N∑j=0Nxijk*pi(k∈K,i≠j),(2.19)∑j=1Ny0jk≤q,(2.20)∑j=0Nz0jk≤1,(2.21)∑i=0Nyijk+(pi-qi)∑i=0Nxijk=∑i=1Nyijk,(2.22)∑j=0N∑k=1Kx0jk=∑j=0N∑k=1Kxj0k=m≤K,(2.23)yijk≤q*xijk,(2.24)yijk≻0,(2.25)zijk≻0.In the model, formula (2.12) is the goal of the problem, that is minimum the cost. Formula (2.13) calculates the time when the vehicle k travels from customer i to j. Constrain (2.14) means that all the vehicles start from the depot and terminate at the depot. Constrain (2.15) and (2.16) represent that each customer can only be served by one vehicle. Constrain (2.17) and (2.19) stand for the capacity of each vehicle cannot exceed the specified capacity. Constrain (2.18) and (2.20) stand for the sum pickups of each vehicle cannot exceed the specified capacity. Constrain (2.21) represents the capacity of the vehicle after it has serviced the customer i, but before it services customer j. Constrain (2.22) means the number of the vehicles arranged should not exceed k. Constrain (2.23) represents that no matter where the vehicle is, the capacity cannot exceed the vehicle’s specified capacity. Constrain (2.24) and (2.25) mean that the capacity and pickups should not lower 0.


## 3. MCPSO Algorithm for VRPTW-SDP

### 3.1. Description of MCPSO Algorithm

The particle swarm optimization (PSO) was proposed by Kennedy and Eberhart in 1995. It was inspired by the social behavior of animals, such as bird flocking and fish schooling. The mechanism of PSO is as follows. Each member of PSO is seen as a particle, and each particle is a potential solution to the problem. The particle keeps memory of its previous locations, and evolves from one generation to another. The experiences of all the particles are transmitted to the global best particle. By using the communication mechanism, PSO optimum search is speeded up.However, PSO also has disadvantages. It sometimes converges to undesired local solution, thus the accuracy of the algorithm can achieve is limited. MultiSwarm Cooperative Particle Swarm Optimizer (MCPSO) proposed in [4] is a variant of PSO. It takes a multiswarm cooperative evolutionary strategy where the master swarms change their particles based on their own knowledge and the knowledge of the particles in the slave swarms, while the slave swarms carry out PSO independently to get the diversity of particles. The mechanism of information transfer for MCPSO is showed in Figure 1 (the arrows stand for information transfer).Figure 1
Master-slave model.Suppose that there areN swarms, among which there are Nm(Nm∈[1,N]) main swarms, the remains are slave swarms Ns (Nm+Ns=N). Each slave swarm evolves according to the PSO. When all the slave swarms are ready with the next new generation, each slave swarm then sends the best local individual to one of the master swarms. The master swarm chooses the best of all received individuals and other main swarms’ experiences and evolves according to the following equations:(3.1)vidM(t+1)=wvidM(t)+c1*r1*(pidM(t)-xidM(t))+c2*r2*(pgM(t)-xidM(t))+c3*r3*(pgQ(t)-xidM(t)),xidM(t+1)=xidM(t)+vidM(t+1),
where M represents one of the master swarm, Q represents the other master swarms, c1 and c2 are accelerating constant, w is inertial weight, c3 is migration coefficient, while r1, r2, and r3  are uniform random sequences in the range [0,1]. Note that the particle’s velocity update in the master swarm is associated with three factorspidM(t): previous best position of the master swarm in dth;pgM(t): best global of the master swarm in dth;pgQ(t): best global of the other master swarms in dth.In what follows we briefly sketch out the multiswarm cooperative particle swarm optimizer (MCPSO) step by step. See Algorithm1.Algorithm 1:Pseudocode for MCPSO algorithm.
Algorithm MCPSO
Begin
Initialize the population // master swarm and slave swarm
Evaluate the fitness value of each particle in the population
Repeat
Do in parallel
Swarmi, 1≤i≤Ns
End Do in parallel
Barrier synchronization  // wait for all processes to finish
Select the fittest global individual (pgS or pgM) from all the swarms
Do in parallel
Swarmj,  1≤j≤Nm
End Do in parallel
Barrier synchronization
Evolve the mast swarm // update the velocity and position using Eqs. (3.1)
Evaluate the fitness value of each particle
Until a terminate-condition is met
End


### 3.2. Particle Encoding Scheme

For our proposed vehicle route problem with time windows, the appropriate expression of particles in MCPSO algorithm is a key issue. After all, the essence of VRPTW is to determine the best route of vehicle from a series of potential solutions. The object is to minimize the sum of all costs. Based on the uncertain number of vehicles in this paper, for the VRPTW withN demand points, 2N dimensions encoding scheme is applied. For each demand point, two dimensions are involved. The first k is the serial number of the vehicle that serviced the demand point. The second r is the serial of the demand point that was serviced by k. Conveniently, each particle corresponds to a matrix of (2,N). The first dimension represents the serial number of vehicles, and the second dimension represents the serial number of customers.


### 3.3. Design of Fitness Function

In this paper, the capacity of vehicles is calculated by(3.2)∑i=0Nyijk+(pi-qi)∑i=0Nxijk=∑i=0Nyjik(k∈K,i≠j).
Both sides of the equation have unknown quantities. It is complicated to convert the constraint by penalty function method. Hence the constraint is dealt with when we program for the VRPTW-SDP. Thus the objective function corresponds to the fitness function:
(3.3)fitness=min⁡⁡z.


## 3.1. Description of MCPSO Algorithm

The particle swarm optimization (PSO) was proposed by Kennedy and Eberhart in 1995. It was inspired by the social behavior of animals, such as bird flocking and fish schooling. The mechanism of PSO is as follows. Each member of PSO is seen as a particle, and each particle is a potential solution to the problem. The particle keeps memory of its previous locations, and evolves from one generation to another. The experiences of all the particles are transmitted to the global best particle. By using the communication mechanism, PSO optimum search is speeded up.However, PSO also has disadvantages. It sometimes converges to undesired local solution, thus the accuracy of the algorithm can achieve is limited. MultiSwarm Cooperative Particle Swarm Optimizer (MCPSO) proposed in [4] is a variant of PSO. It takes a multiswarm cooperative evolutionary strategy where the master swarms change their particles based on their own knowledge and the knowledge of the particles in the slave swarms, while the slave swarms carry out PSO independently to get the diversity of particles. The mechanism of information transfer for MCPSO is showed in Figure 1 (the arrows stand for information transfer).Figure 1
Master-slave model.Suppose that there areN swarms, among which there are Nm(Nm∈[1,N]) main swarms, the remains are slave swarms Ns (Nm+Ns=N). Each slave swarm evolves according to the PSO. When all the slave swarms are ready with the next new generation, each slave swarm then sends the best local individual to one of the master swarms. The master swarm chooses the best of all received individuals and other main swarms’ experiences and evolves according to the following equations:(3.1)vidM(t+1)=wvidM(t)+c1*r1*(pidM(t)-xidM(t))+c2*r2*(pgM(t)-xidM(t))+c3*r3*(pgQ(t)-xidM(t)),xidM(t+1)=xidM(t)+vidM(t+1),
where M represents one of the master swarm, Q represents the other master swarms, c1 and c2 are accelerating constant, w is inertial weight, c3 is migration coefficient, while r1, r2, and r3  are uniform random sequences in the range [0,1]. Note that the particle’s velocity update in the master swarm is associated with three factorspidM(t): previous best position of the master swarm in dth;pgM(t): best global of the master swarm in dth;pgQ(t): best global of the other master swarms in dth.In what follows we briefly sketch out the multiswarm cooperative particle swarm optimizer (MCPSO) step by step. See Algorithm1.Algorithm 1:Pseudocode for MCPSO algorithm.
Algorithm MCPSO
Begin
Initialize the population // master swarm and slave swarm
Evaluate the fitness value of each particle in the population
Repeat
Do in parallel
Swarmi, 1≤i≤Ns
End Do in parallel
Barrier synchronization  // wait for all processes to finish
Select the fittest global individual (pgS or pgM) from all the swarms
Do in parallel
Swarmj,  1≤j≤Nm
End Do in parallel
Barrier synchronization
Evolve the mast swarm // update the velocity and position using Eqs. (3.1)
Evaluate the fitness value of each particle
Until a terminate-condition is met
End


## 3.2. Particle Encoding Scheme

For our proposed vehicle route problem with time windows, the appropriate expression of particles in MCPSO algorithm is a key issue. After all, the essence of VRPTW is to determine the best route of vehicle from a series of potential solutions. The object is to minimize the sum of all costs. Based on the uncertain number of vehicles in this paper, for the VRPTW withN demand points, 2N dimensions encoding scheme is applied. For each demand point, two dimensions are involved. The first k is the serial number of the vehicle that serviced the demand point. The second r is the serial of the demand point that was serviced by k. Conveniently, each particle corresponds to a matrix of (2,N). The first dimension represents the serial number of vehicles, and the second dimension represents the serial number of customers.


## 3.3. Design of Fitness Function

In this paper, the capacity of vehicles is calculated by(3.2)∑i=0Nyijk+(pi-qi)∑i=0Nxijk=∑i=0Nyjik(k∈K,i≠j).
Both sides of the equation have unknown quantities. It is complicated to convert the constraint by penalty function method. Hence the constraint is dealt with when we program for the VRPTW-SDP. Thus the objective function corresponds to the fitness function:
(3.3)fitness=min⁡⁡z.


## 4. Experimental Study

We consider the following situatoion. 8 customers need served by one depot. The depot has 5 vehicles. Each of the vehicles has a capacity of 3 units and velocity of 10 units. The distance between each customer and the demands supplies of each customer is shown in Table1. The distance between depot and each customer is shown in Table 2. The service time for each customer is 0.Table 1
The supplies and demands and time window of each customer.
Customeri
1
2
3
4
5
6
7
8
Demandsgi
0.8
1.5
0.8
0.7
1.4
1.5
0.6
0.8
Suppliespi
0.6
1
1
1
0.8
1.2
0.8
0.8
Time window
[1,3]
[2,4]
[1,2]
[3.5, 7]
[3.5, 5]
[2,5]
[2,5]
[1,4.5]Table 2
The distance between depot and each customer.
d
i
j
0
1
2
3
4
5
6
7
8
0
0
10
9
7
8
8
8
3
4
1
0
4
9
14
18
18
13
14
2
0
5
10
14
17
12
13
3
0
5
9
15
10
11
4
0
6
13
11
12
5
0
7
10
12
6
0
6
8
7
0
2
8
0Base on Matlab 7.0, we apply MCPSO algorithm to solve this VRPTW-SDP. In MCPSO algorithm, the number of swarms isn=4. The 4 swarms conclude 3 slave swarms and 1 master swarm. The migration coefficient of the master swarm is  c1=c2=c3=1.333. The migration coefficient of the slave swarm is c1=c2=2. To validate its performance, the MCPSO algorithm is also compared with GA and PSO.For fair comparison, in all cases the population sizenum of GA/PSO and MCPSO was set at 100 (all the swarms of MCPSO include the same particles) and a fixed number of maximum iterations 100 is set in GA/PSO and MCPSO. In the fitness function, we set four parameters as follows: a=1,b=5,s1=1,s2=2. The results are shown in Tables 3 and 4. To serve the 8 customers, 3 vehicles are arranged. The first vehicle serves customer 2, 8 and 4 in sequence. The second vehicle serves customer 3 and 6 in turn. The last one serves customer 1, 5, and 7 successively. The total fees are 549 units. From Table 4, we can see that the time points when vehicles arrive at customers are within the time windows. Figure 2 shows MCPSO algorithm, which converge on the minimum point quickly.Table 3
The route and fees.
Vehicle
Route of vehicle
Fee
1
0
→
2
→
8
→
4
→
0
549
2
0
→
3
→
6
→
0
3
0
→
1
→
5
→
7
→
0Table 4
The time the vehicle arrives at each customer.
Customer
1
2
3
4
5
6
7
8
Time
1
0.9
0.7
3.4
2.8
2.2
3.8
2.2Figure 2
The Convergence graph for MCPSO.The ten fitness values are obtained in ten runtimes in GA algorithm, PSO algorithm and MCPSO algorithm. The final result including maximum, minimum, mean, and variance is shown in Table5. It is apparently that the maximum, minimum and mean of MCPSO are smaller than PSO, whose are smaller than GA. Although the difference on maximum, minimum, and mean between PSO and MCPSO are small, the gap on variances is big. This means that the MCPSO is more stability. The convergence graph for MCPSO, PSO and GA can be seen in Figure 3. MCPSO achieves a minimal solution 549, which is lower than 570 (PSO) and 675 (GA). It is obvious that MCPSO is more suitable at least for our proposed vehicle route problem with time windows.Table 5
Summary result.
Algorithm
Maximum
Minimum
Mean
Variance
PSO
642
522
571.2
1296
GA
708
631
676.4
557.6
MCPSO
604
500
551.5
456.5Figure 3
The convergence for MCPSO, PSO, and GA.


## 5. Conclusions

The VRPTW-SDP studied in this paper is an extension of VRPTW, where the number of vehicles with simultaneous delivery and pickup is uncertain. The proposed problem is meaningful for transport companies in real world. To solve the proposed problem, we used MCPSO algorithm. A new encoding method is proposed, which may be somewhat benefit to deal with other problems. Furthermore, with comparison to PSO algorithm and GA algorithm, the results show that the simple and robust MCPSO algorithm is more effective for the proposed VRPTW-SDP.



---

*Source: 104279-2012-08-21.xml*