Received July 30, 2021, accepted August 5, 2021, date of publication August 9, 2021, date of current version August

17, 2021.
Digital Object Identifier 10.1109/ACCESS.2021.3103640

Multi-Depot Split-Delivery Vehicle

Routing Problem
HYUNPAE LIM1 , GYU M. LEE 2, AND IVAN KRISTIANTO SINGGIH 3
1 Schoolof Mechanical, Industrial and Manufacturing Engineering, Oregon State University, Corvallis, OR 97331, USA
2 Department of Industrial Engineering, Pusan National University, Busan 43241, Republic of Korea
3 Department of Industrial and Systems Engineering, Korea Advanced Institute of Science and Technology, Daejeon 34141, Republic of Korea

Corresponding author: Gyu M. Lee (glee@pnu.edu)

This work was supported by the National Research Foundation of Korea (NRF) grant funded by the Korea government (MSIT)
(No. 2018R1A2B3008890).

ABSTRACT The rapid advancements in information technologies and globalization change the way of
distributing goods to customers. Many enterprises have multiple factories, warehouses, and distribution
centers and strive for competitive efficiency in the distribution operations to minimize transportation costs.
This study proposed the mixed-integer programming (MIP) model for the multi-depot split-delivery vehicle
routing problems (MDSDVRPs) with hetero vehicles, allowing multiple visits to a customer. A genetic
algorithm (GA) with a novel two-dimensional chromosome representation has been proposed with dynamic
mutation policies. The process parameters of the proposed GA are optimized using the Taguchi method.
The proposed algorithms showed the benefits of split-delivery in MDSDVRPs and showed the competitive
performance even for the classical single-depot vehicle routing problems with no split-delivery.

INDEX TERMS Vehicle routing problem (VRP), multi-depot split-delivery VRP, genetic algorithm (GA),
taguchi method.

I. INTRODUCTION the advances in electronics and new technologies, e.g., VRPs

Logistics and distribution systems change our lives rapidly for electric vehicles ([1]–[3]) and drones ([4], [5]).
due to the changes in customer demands and enabling tech- The VRP was initially introduced by Dantzig and
nologies while saving significant distribution costs. How- Ramser [6], and it has been widely studied thereafter. They
ever, traditional logistics and distribution operations exhibit described a real-life application concerning the delivery of
an underlying limitation in modeling and implementing gasoline to service stations. Fisher [7] describes the problem
shared services from multiple manufacturers, warehouses, as finding the efficient use of a fleet of vehicles that must
and depots. make several stops to deliver passengers or goods. The term
The vehicle routing problem (VRP) is a problem in which ‘‘customer’’ is used to denote the stops to make. Every cus-
a set of routes for a fleet of vehicles based at one or several tomer has to be assigned to precisely one vehicle in a specific
depots must be determined for a certain number of geograph- order.
ically dispersed customers. The objective of the VRP is to A particular case of the VRP arising when only one
minimize the total distance traveled by all vehicles, which vehicle is available at a depot and no additional opera-
can be considered transportation or delivery costs. Recently, tional constraints are imposed, i.e., traveling salesman prob-
with the increase in fuel prices, the importance of minimizing lem (TSP), is extensively described by Lawler et al. [8],
delivery costs has been emphasized as a critical factor that can Knox [9], Barvinok et al. [10], Engebretsen and Karpin-
reduce the total costs of production and distribution. In recent ski [11], Ouaarab et al. [12], and Mahi et al. [13]. The TSP
decades, various engineering areas produce many variants of has one vehicle, one depot, and multiple customers. The
VRPs to utilize the theory and to optimize their systems with customer demands are assumed to be satisfied with one visit
to each customer by a vehicle. Another version of the VRP
is capacitated VRP (CVRP). The CVRP has n customers and
The associate editor coordinating the review of this manuscript and a single depot with several vehicles of an identical capacity.
approving it for publication was Kuo-Ching Ying . The vehicles must accomplish the delivery with the minimum
This work is licensed under a Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 License.
112206 For more information, see https://creativecommons.org/licenses/by-nc-nd/4.0/ VOLUME 9, 2021
H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

total travel cost, where the cost is the distance dij from nodes Golden et al. [50], Laporte [51], Ghorbani et al. [52], and
i to j (i, j ∈ {0, 1, . . . , n}), where 0 stands for a single depot Anuar et al. [53]. They could be divided into two main
and n is the number of customers. The applications of the classes: classical heuristics, mostly between 1960 and 1990,
CVRP can be found in Desrosiers et al. [14], Osman [15], and metaheuristics from 1990 ([44], [51]).
Laporte [16], Toth and Vigo [17], Lysgaard et al. [18], This paper addresses the VRP, which has heterogeneous
and Alssager et al. [19]. Some studies considered hetero- vehicles departing from multiple depots, allowing split deliv-
geneous vehicles to reduce delivery costs by dispatching eries to customers. This multi-depot split-delivery VRP
appropriate vehicles to the routes according to customer (MDSDVRP) was firstly studied by Lim [54] and then
demands ([20]–[22]). mentioned by Gulczynski et al. [55] independently. The
The vast majority of papers have been published on a clas- research by Gulczynski et al. [55] has been cited more than
sical single-depot capacitated VRP (SDCVRP) ([23], [24]) 100 times. They defined the MDSDVRP and developed an
and some others were found dealing with problems known integer programming-based heuristic. Their heuristic deter-
as multiple-depot capacitated VRP (MDCVRP) ([25]–[28]). mines the reduction in traveled distance, allowing split deliv-
The MDCVRP is an extension of SDCVRP with vehicles eries among vehicles based at the same depot and different
starting from different depots. The MDCVRP has similar depots.
constraints to those of the SDCVRP, except for the require- The study completed earlier by Lim [54] has not been
ment that each vehicle starts from and finishes the delivery published until now. It defined the MDSDVRP and proposed
at the same depot. In the MDCVRP, if customers are clus- the first mixed-integer programming model for MDSDVRP.
tered around depots, then the problem can be modeled as The model is solved optimally using a CPLEX solver. The
a set of independent SDCVRPs. However, if customers and addition of this study to MDSDVRP literature may enrich
depots are intermingled, the problem must be modeled as the future research.
MDCVRP. The applications of the MDCVRP can be found The classical heuristics can be divided into three groups:
in [23], [25], [26], and [29]–[32], all using adaptations of construction methods, two-phase methods, and improvement
classical SDCVRP procedures. methods [56]. Construction methods gradually build a feasi-
The VRP is one of the combinatorial optimization prob- ble solution by selecting arcs based on minimizing cost. The
lems belonging to the non-deterministic polynomial-time two-phase method divides the problem into two stages: clus-
hard (NP-hard) class [33] which cannot be solved to optimal- tering customers into feasible routes disregarding their order,
ity within polynomially-bounded computational time [34]. and constructing routes. One of the two-phase methods is the
The NP-hardness of VRPs can be proven via ‘‘proof by sweep algorithm in Laporte et al. [44]. Improvement methods
restriction’’. A problem is NP-complete when we can prove start with a feasible solution and improve it by exchanging
that it contains a known NP-complete problem as a spe- arcs or nodes within or between the routes. The local search
cial case [35]. Likewise, we can prove that VRP contains a algorithms developed by Aarts and Lenstra [57] belong to the
traveling salesman problem (TSP) when we limit the num- improvement heuristics. The advantage of classical heuristics
ber of salesmen to one. Since TSP is NP-hard ([36], [37]), is that they have a polynomial running time [44]. When using
then VRP (multiple TSP [38]) is NP-hard as well. Solving them, one can provide good solutions within a reasonable
NP-hard problems of large sizes using exact solution methods time [58]. On the other hand, they only perform a limited
is expensive in computational efforts. It is often impossible search in the solution space. Therefore, they have a risk of
due to the limited computational power. Hence, the extensive resulting in a local optimum.
efforts in studying the approximation algorithms or heuris- During the past few decades, there have been many
tics were justified by many researchers. The optimal solu- attempts to solve VRPs quickly and effectively by using meta-
tions for NP-hard problems can be obtained by using those heuristics such as tabu search (TS), simulated annealing (SA),
algorithms and matching their solutions with the theoretic genetic algorithm (GA), and ant colony optimization (ACO)
bounds ([39]–[42]). algorithm ([44], [59], [60]). Braekers et al. [61] reported
Many different approaches have been developed to solve that metaheuristic was the most applied method to solve
this NP-hard problem. The branch-and-bound method has VRP cases. Metaheuristics were used more often than exact
been used for small problems with only a few customers [43]. methods, classical heuristics, and simulation approaches. The
Most approaches for large problems are based on heuristics, TS and the SA are local search-based algorithms that move
i.e., approximation algorithms aiming to find good feasible from one solution to another in the neighborhood until a stop-
solutions quickly [44]. Many models and algorithms have ping criterion is met. Many different TS heuristics have been
been proposed to obtain the optimal or approximate solutions proposed with unequal success. Rochat and Taillard [62] used
for different variants of the VRP. A thorough classification the TS heuristic to solve some benchmark VRPs. Osman [63]
was given in Desrochers et al. [45]. Laporte and Nobert [46] obtained similar results using the SA. Unlike TS and SA,
presented an extensive survey devoted to exact methods GA maintains a population of good solutions that are recom-
for VRPs. Other surveys about VRP studies were reported bined to produce new solutions. A considerable amount of
by Laporte [16], Toth and Vigo [17], Bodin et al. [33], research on the GA has recently been done to solve many vari-
Christofides et al. [47], Magnanti [48], Christofides [49], ants of VRPs, including VRPs with time windows (VRPTW)

VOLUME 9, 2021 112207

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

([64], [65]), where each customer has a time window for Kang [72] introduced the VRP, allowing multiple visits to a
which the vehicle has to arrive. Berger and Barkaoui [65] customer using a heuristic method. Some other split delivery
presented a hybrid GA to solve the CVRP. Their HGA uses VRP studies are Archetti et al. [73], Berbotto et al. [74],
two populations of solutions that periodically exchange some Bianchessi et al. [75], and Chen et al. [76].
chromosomes, which are the feasible solutions to the CVRP. This study differs significantly from Gulczynski et al. [55].
The algorithm is competitive in comparison to the best TS They used pre-generated routes using a combined algorithm
heuristics. However, Renaud et al. [26] reported that such as input data. Then, the integer programming model was
heuristics require substantial computing times and several used to move some delivery tasks from one route to another.
parameter settings. When dealing with VRP, the GA was The routes were then updated only considering pre-defined
applied more often than most local search-based metaheuris- input sets. In contrast to Gulczynski et al. [55], this study
tics [66], and GA was the most applied method among formulates a complete MDSDVRP mixed-integer program-
population-based metaheuristics ([66], [67]). It could be con- ming (MIP) model and optimally solves it using the commer-
cluded that the GA has been proven to perform effectively cially available solver. This study also proposes a GA to solve
compared with other methods used for solving VRPs. the MDSDVRPs of large sizes effectively and efficiently. The
The GA is a randomized, global search algorithm that proposed GA achieves optimality for small-sized problems
solves problems by imitating genetic processes observed and demonstrates its effectiveness by solving single-depot
during natural evolution and has been extensively used benchmark VRPs to the best-known solutions for comparison
to tackle many combinatorial problems, including various purposes.
VRPs ([60], [68]–[70]). In the GA, a population of chromo- The objectives of this paper are as follows: 1) to generalize
somes (individuals) or solutions is maintained during the evo- the MDSDVRPs by removing the constraints of the number
lution, in which solution evaluation, selection, crossover, and of depots and the number of visits allowed to each customer;
mutation occur. The quality of the solution is evaluated by its 2) to develop and validate a MIP model to achieve the opti-
fitness function, which represents an individual’s survivabil- mality; 3) to propose and validate a GA to solve effectively
ity in the wild. This fitness determines the individuals for the and efficiently the medium or large VRPs with heterogeneous
crossover or mating, which produces offsprings for the next vehicles from multiple depots, allowing splits deliveries; and
generation. The mutation is also used to prevent local conver- 4) to optimize the parameters of the proposed GA using the
gence by diversifying the search space. The average quality of Taguchi method and understand the effects of parameters on
the population gradually improves as new and better solutions the performance.
are generated and worse solutions are removed. Analogous To our best knowledge, this study presents the first MIP
to biological processes, offspring with relatively good fitness formulation for MDSDVRP and a novel GA. An interesting
levels are more likely to survive and reproduce, expecting that study similar to the problem in this study is Ray et al. [77].
fitness levels throughout the population may improve as they Their research is closer to the location problem as a variant
evolve. More details can be found in Reeves [71]. of our proposed MDSDVRP. They focused on determining
Most prior research on the VRP has considered capacitated the depot locations in which multi-depot shared commodity
vehicles from a depot while only allowing a visit to each cus- delivery by vehicles was completed. They proposed a heuris-
tomer. Vehicles dispatched from a single depot must deliver tic to determine where to locate the depots among customer
the required amount to customers, satisfy all demands, and node candidates, while this study focuses on the efficient
finally return to the depot. The vehicle routes are designed routing among given depot locations. Our proposed model
so that each customer is visited only once by precisely one for the MDSDVRPs can limit the maximum number of visits
vehicle, and the total demands of all customers on a partic- to each customer, and it provides the decision-makers more
ular route must not exceed the vehicle’s capacity. However, control in planning the distribution.
the constraints of homogeneously capacitated vehicles, a sin- This study also proposes a GA to solve the proposed MDS-
gle depot, and one allowed visit to customers are unrealistic DVRP with the limited number of vehicle visits to a customer
in the real world. Therefore, this research presents a genetic by developing a novel two-dimensional chromosome design.
algorithm (GA) to find all routes that minimize the total A chromosome representation is designed to tackle the com-
distance traveled by heterogeneous vehicles from multiple plexity of MDSDVRPs and maintain the efficiency of the GA
depots with split deliveries. Split deliveries enable multiple simultaneously. The chromosome design must have a way to
visits to a customer to satisfy his/her demand. If multiple represent the sequence of vehicles visiting customers, and the
visits are allowed, it is conjectured that they may reduce the sequence of customers vehicles visit simultaneously while its
number of vehicles, travel distances, and green gas emissions crossover and mutation maintain their efficiency. Hence, this
to satisfy the customer demands. study can play an important role in facilitating the research in
Assume that there are five customers, each of whom MDSDVRPs.
demands a little more than half of the homogeneous vehicle’s An essential characteristic of the GA is its non-
capacity. To satisfy all demands by only one visit to each cus- deterministic evolution, i.e., it is stochastic in natural deci-
tomer, 5 vehicles are necessary. Using split deliveries, only sions, making the GA more robust than other heuristics.
three vehicles might be needed to solve the problem. Shin and This evolution by GA can be determined by a set of

112208 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

parameters, including population size, crossover rate, muta- Cmt Capacity of vehicle t from depot m for
tion rate, elitism rate, terminal condition, etc., which signif- dij Distance between nodes i and j (1 ≤ i, j ≤ N + M)
icantly affects the performance of GA. This paper optimizes V Maximum number of visits to a customer
the parameters of the proposed GA using the Taguchi method B A large number
to make the proposed GA robust to different problems in
MDSDVRP. Decision Variables:
Our paper is presented as follows. In Section 2, we describe Ujmt Unloaded amount by vehicle t from depot m at
the MDSDVRP completely and propose a MIP model. customer j, where 1 ≤ j ≤ N .
In Section 3, we introduce the details of the proposed GA. (Ujmt = 0 for N + 1 ≤ j ≤ N + M , Ujmt = 0 for
In Section 4, the parameters of the proposed GA have been t∈/ STm )
optimized using the Taguchi method. Section 5 provides the xijmt 1, if vehicle t from depot m travels from node i
computational results to demonstrate the effectiveness of our to j, where t ∈ STm ; 0, otherwise.
GA through some numerical experiments and the compar- yimt auxiliary variable for sub-tour elimination
isons with optimal solutions obtained by the proposed MIP
model. Finally, we conclude our study and list future research The MDSDVRP with heterogeneous vehicles from multi-
topics in Section 6. ple depots, allowing split deliveries to the customers, is for-
mulated as a MIP model in the following:
II. PROBLEM DESCRIPTION AND MATHEMATICAL XL XL XM XTm
Min dij xijmt
MODEL i=1 j=1 m=1 t=1
This section presents a MIP model for the MDSDVRP with Subject to
heterogeneous vehicles from multiple depots, allowing split xijmt = 0, ∀i 6= j ∈ SN , ∀m ∈ SM , ∀m ∈ STm
deliveries. The VRP under consideration can be represented
(1)
as a network, where nodes are customers or depots, and the XL XL
links between pairs of nodes are the roads. In a network, there xijmt = xjimt
i=1 i=1
are N customers with known demands Di (i = 1, . . . , N ), and ∀j ∈ SN , ∀m ∈ SM , ∀m ∈ STm (2)
M depots, each of which has Tm (m = N + 1, . . . , N + M ) XTm XM XL
vehicles. All vehicles may have homogenous or heteroge- xijmt ≤ V ∀j ∈ SN (3)
t=1 m=1 i=1
neous capacities. XL
xiimt = 0 ∀m ∈ SM , ∀m ∈ STm (4)
The assumptions for the MDSDVRP in this study are i=1
detailed in the following. XL
B xiimt ≥ Ujmt ∀j ∈ SN , ∀m ∈ SM ,
i=1
• Each vehicle must start and finish its route at a single
∀m ∈ STm (5)
depot.
• Customer demands must be satisfied by vehicles within xijmt ≤ Ujmt ∀i, j ∈ SN , ∀m ∈ SM , ∀m ∈ STm
the given maximum number of visits. (6)
• The vehicles’ capacities are known (homogenous or
XTm XM
Ujmt = Dj ∀j ∈ SN , ∀m ∈ SM ,
heterogeneous). t=1 m=1
• The sum of unloaded amounts at customers from a vehi- ∀m ∈ STm (7)
XN
cle must not exceed the vehicle’s capacity. Ujmt ≤ Cmt ∀m ∈ SM , ∀m ∈ STm (8)
• The locations of all customers and depots are given. j=1
• The distances between all pairs of locations are known. yimt − yjmt + Lxijmt ≤ L − 1 ∀i 6= j ∈ SN ,
The definitions of constants, variables, and sets used in the ∀m ∈ SM , ∀m ∈ STm (9)
MIP formulation are given as follows. xijmt = {0, 1} ∀i, j ∈ SN , ∀m ∈ SM , ∀m ∈ STm
Sets: (10)
SN Set of customer indices yjmt > 0, integer ∀j ∈ SN , ∀m ∈ SM , ∀m ∈ STm
SM Set of depot indices
(11)
S Set of indices for all customers and depots; S =
SN ∪ SM The objective function of the MIP model is to minimize the
STm Set of all vehicle indices at depot m total traveled distance by all vehicles to satisfy all customers’
demands. Constraint (1) ensures that each vehicle starts from
Parameters: its origin depot and terminates its route at the same depot.
N Number of customers In other words, each vehicle cannot visit the depots other
M Number of depots than its origin depot. Constraint (2) ensures that all vehicles
L Number of customers and depots (L = N + M) visiting a node must leave that node. The number of visits
Tm Number of vehicles at depot m must be the same as the number of departures for each vehicle
Di Demand of customer i (1 ≤ i ≤ N) at each node. It ensures the continuous flow of vehicles in the

VOLUME 9, 2021 112209

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

the route of V21 is [D2 − C2 − C3 − D2 ], and the route of V31

is [D3 − C5 − D3 ], where D2 and D3 represent depots 2 and
3, respectively.

B. INITIAL POPULATION AND FITNESS FUNCTION

The initial population of the pre-determined size is randomly
generated. However, a way to obtain a good initial population
FIGURE 1. The proposed representation of a chromosome, where
5 customers and 3 vehicles in depot 1 (V11 , V12 and V13 ), 1 vehicle in
is desired since it significantly impacts the performance of
depot 2 (V21 ) and 1 vehicle in depot 3 (V31 ), allowing up to 3 visits to a the proposed GA but it is beyond the scope of this study.
customer. The objective of the proposed GA is to minimize the overall
distance traveled by all vehicles. The fitness, which is the
network. Constraint (3) ensures that each customer node can survival chance of a feasible solution satisfying all customer
have up to V visits by vehicles to satisfy customer demand. demands, is an inverse of the sum of the total distance traveled
Constraint (4) prevents the loop of any route at a node. Con- by all vehicles.
straints (5) and (6) ensure that if vehicle t from depot m travels
from node i to node j, the vehicle should unload Ujmt at node C. PARENTS SELECTION
j. Constraint (7) ensures that the sum of the unloaded amounts In GA, an appropriate method to select chromosomes for
at a customer node j should be the same as its demand. crossover must be employed to give more chances to those
Constraint (8) ensures that the total unloaded amounts of each fittest chromosomes in a population. The genetic search
vehicle over its route cannot exceed the vehicle’s capacity. terminates prematurely; with too little chance, evolutionary
Constraint (9) presents the sub-tour elimination constraint. progress is slower than necessary. Typically, lower selection
Constraints (10) and (11) are the binary and integer variable pressure is desirable at the start of the genetic search in favor
constraints. of a broad exploration of the search space, while a higher
selection pressure is recommended at the end to converge
III. PROPOSED SOLUTION METHOD efficiently. The roulette wheel selection method and linear
As forementioned, the VRP is known as an NP-hard com- scaling method have been used during the selection process
binatorial problem. It is difficult to solve even small prob- in the proposed GA [78]. Based on several test runs of the
lems optimally in a reasonable amount of time. The GA proposed GA, the linear scaling function has been chosen as
has been applied successfully in many combinatorial opti- fi0 = 0.1fi + 1, where fi0 and fi are the scaled and raw finesses
mization problems. The GA does not guarantee optimality for chromosome i.
because of its stochastic nature, but it finds a good near-
optimal solution in significantly less time. In this section, D. CROSSOVER
the proposed GA implemented for MDSDVRPs in this paper According to the crossover rate, two chromosomes from the
is described in detail. current population are selected for mating through the selec-
tion process, i.e., a probability of crossover. If a randomly
A. CHROMOSOME REPRESENTATION generated number between 0 and 1 is smaller than the given
A way to encode a solution of the problem into a chromosome crossover rate, these chromosomes reproduce to form new
has a high impact on the GA’s performance. The proposed offsprings to be included in the next generation. Otherwise,
representation is a 2-dimensional (V + 1) × N matrix. Its the crossover does not take place. An appropriate method
columns represent N customers. The first row contains ran- to crossover a pair of selected chromosomes to improve the
domly generated sequences of visiting orders, and the other finesses of offsprings has been proposed for the MDSDVRP
V rows contain the vehicles visiting each customer. The under consideration in this paper.
maximum number of visits to a customer is limited to V . Two different crossover methods are used to produce off-
A chromosome representation is illustrated as a 4 × 5 matrix springs in this paper. One is the position-based crossover
in Figure 1, where V = 3 and N = 5. method, which is applied to the first row. Because the first
The way to interpret the chromosome representation in Fig- row represents the visiting order of vehicles, the row must
ure 1 is explained in the following. There are 3 vehicles in be ensured to have no same gene. The other is the uniform
depot 1 (V11 , V12 and V13 ), 1 vehicle in depot 2 (V21 ) and crossover method, which is superior to other crossover strate-
1 vehicle in depot 3 (V31 ). The vehicles serve five customers, gies for combinatorial problems [79], and it is applied to
C1 , C2 , C3 , C4 and C5 . The V11 will visit C1 , C3 and C4 , the other rows. These crossovers are described below and
respectively. The visiting order of the V11 depends on the illustrated in Figures 2, 3, and 4.
values in the first row. Since the corresponding values in Step 1: Select a set of cells from Parent 1 randomly, which
the first row for C1 , C3 and C4 are 3, 4 and 5, respectively, are shaded in Figure 2.
the route of V11 is [D1 − C1 − C3 − C4 − D1 ], where D1 Step 2: For the first row, copy the selected cells into
represents depot 1. In the same manner, the route of V12 is offspring at the corresponding locations and delete the cor-
[D1 − C2 − C5 − C1 − D1 ], the route of V13 is [D1 − C3 − D1 ], responding values from Parent 2 (position-based crossover).

112210 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

FIGURE 2. Step 1 of the crossover.

FIGURE 3. Step 2 of the crossover.

random changes prevent pre-mature local convergence. The

proposed GA has a relatively simple mutation procedure with
elitism. All chromosomes in the population except the elites
are subject to mutation at the mutation rate. The elitism rate
is set to 10%, which is the percentage of best solutions in
each population and these elites are immune to the mutation.
FIGURE 4. Step 3 of the crossover.
While most GA implementations use the static mutation rate,
the proposed GA introduces the dynamic self-adapting muta-
tion rate. The mutation rate has been dynamically adapted
during the evolutions. It starts with an initial value and
For the other rows, copy the selected cells into offspring then increases by the fixed or logarithmic amount when-
at the corresponding locations and delete the cells at the ever no improvement is observed over a certain number of
corresponding locations from Parent 2 (uniform crossover). generations. This number of generations is called mutation
See Figure 3. crank-up interval in this paper. If the best solution improves,
Step 3: For the first row, copy the remaining cells from the mutation rate drops to the initial mutation rate. It increases
Parent 2 into empty cells of offspring from left to right in the global search capability to escape from the pre-mature
the order of cells showing up in Parent 2 (position-based local optima and to search for better solutions from diverse
crossover). For the other rows, copy the remaining cells from directions. The implementation of adaptive mutation rate
parent 2 into the empty cells of offspring at the corresponding has been proved effective in solving the various VRPs
locations (uniform crossover). See Figure 4. ([80]–[82]). The proposed GA uses the inversion mutation
The crossover procedure described in this study may con- method, which is explained in the following steps and illus-
tain the same vehicle more than once in a column, i.e., a trated in Figure 5.
vehicle unnecessarily visits a customer more than once. Then, Step 1: Select a pair of columns randomly in a selected
a repair rule has been performed on the generated offspring chromosome, which are shared in Figure 5.
by removing redundant visits for the same customer. Step 2: Swap the corresponding cells between these two
columns except the cells on the first row. Note that the cells
E. MUTATION of the first rows are not swapped.
The mutation is another important operator in the GA and is
applied to a chromosome at a mutation rate. Syswerda [79] F. TERMINATION CONDITIONS
proved that the mutation operator could sometimes play a The proposed GA uses two termination criteria. One is that
more crucial role than the crossover. Therefore, the crossover the proposed GA completes a specified maximum num-
and mutation operators need to well-designed per the problem ber of generations in this implementation. The other is
on hand. The mutation operator brings random changes into a that the proposed GA can also be terminated due to no
single chromosome. If a randomly generated number between improvement over a specified number of generations. After
0 and 1 is smaller than the mutation rate, a chromosome the proposed GA terminates, the chromosome with the
reproduces a new member to be included in the next gen- highest fitness is interpreted as the best solution in that
eration. Otherwise, the mutation does not take place. These run.

VOLUME 9, 2021 112211

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

FIGURE 5. Mutation of the proposed GA.

IV. TAGUCHI METHOD FOR PARAMETER TUNING TABLE 1. Four factors with three levels per factor.
Many GA process parameters have been defined to solve the
combinatorial optimization problems effectively. The values
of those process parameters need to be carefully selected.
Good process parameter setting is important for the GA
to obtain a good final solution. Usually, it is difficult to
determine a good set of process parameters because their
relationships can be rather complicated and unclear.
To fine-tune the performance of algorithms or processes,
many parameters must be set carefully. The Taguchi method
for parameter tuning is an important tool for robust design.
Robust design is an engineering methodology for optimizing
the product and process conditions that are minimally sensi-
tive to the causes of variation. The orthogonal array and the
signal-to-noise ratio (SNR) are two primary tools used in the
Taguchi method. Additional details can be found in the books
(100, 150, 200) for the population size, (0.6, 0.7, 0.8) for
presented by Taguchi et al. [83] and Wu [84].
the crossover rate, (50, 75, 100) for the mutation crank-up
An orthogonal array is a fractional factorial matrix, which
interval, and (Mutation Policy 1, Mutation Policy 2, Muta-
assures a balanced comparison of levels of any factor. It is
tion Policy 3) for the mutation policy. In Mutation Policy 1,
a matrix of numbers arranged in rows and columns where
the static mutation rate of 0.05 has been used. Mutation Policy
each row represents the level of the factors in each run, and
2 increases the mutation rate linearly from 0.05 by 0.10
each column represents a specific factor that can be changed
whenever no improvement is made for the mutation crank-
from each run. The symbol of three-level orthogonal arrays
up interval. Mutation Policy 3 uses the logarithmic increase
is Ln(3k), where n is the number of experimental runs, 3 is
of the mutation rate from 0.05. The logarithmic increase of
the number of levels for each factor, and k is the number of
Mutation Policy 3 instead of linear increase in MP 2 (possible
factors. The letter L comes from Latin since the orthogonal
mutation rates are 0.05, 0.15, 0.25) can be calculated as
arrays were associated with Latin square designs from the
outset. An = 0.05 + ln (1 + 0.15n),
SNR is the signal ratio over the noise, which measures the
where n = {1, 2, 3} (possible mutation rate can be 0.05,
strength of the signal with the existence of noises. A higher
0.19, 0.31). Elitist rate is used to retain 10% of the best
SNR means that a process or a design is the more robust.
chromosomes at each generation.
Suppose that we have a set of experiment runs. Since the
Table 1 presents four process parameters with three levels
objective of MDSDVRPs is to minimize the total traveled
in the proposed GA. To conduct the full factorial experiment
distance, the smaller-the-better is an appropriate measure in
with all factors, 34 (or 81) experiments are necessary to deter-
this study. The following formulation for smaller-the-better
mine the optimal process parameters. However, the Taguchi
characteristics is used;
method only requires 9 runs to optimize the process parame-
n
1X 2 ters when L9 orthogonal array is used.
SNR = −10 log( xi ). Each row in Table 2 shows 9 experiments with process
n
i=1
parameters A, B, C, and D at their corresponding levels.
Through extensive computational experiments for various To account for the characteristics of stochastic disturbance in
problems, four process parameters are identified as important the proposed GA, each experiment has been tested 40 times,
design factors for the performance of the proposed GA; pop- i.e., 10 runs for each of 4 test problems. These 4 test problems
ulation size, crossover rate, mutation crank-up interval, and are introduced in the following sections. Various sizes and
mutation policy. For each design factor, three levels of the structures of VRPs can be considered noise factors. The
process parameters are chosen from the previous research; relative gap between the i-th experimental solution and the

112212 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

TABLE 2. SNR values of the L9 experiments.

best-known solution, Gi, is calculated as TABLE 3. SNR values of the process parameters.

ith solution − best known solution

Gi =
best known solutoin
where i = {1, 2, . . . , 40}.
Since the smaller Gi is more desirable, the smaller-the-
better SNR calculation has been used as
n
1X 2
SNR = −10 log( xi ),
n
i=1
where n = 40. The SNR values of process parameter A are
calculated as
SNRA1 = SNR1 + SNR2 + SNR3 ,
SNRA2 = SNR4 + SNR5 + SNR6 ,
SNRA3 = SNR7 + SNR8 + SNR9 , by noise. Using the optimal process parameters suggested
where SNRi represents the SNR value of the ith run and by the Taguchi method, the proposed GA’s searchability has
SNRAi denotes the aggregated SNR values of level i of process been improved and the proposed GA has generated better
parameter A. solutions.
The optimal levels of process parameters A, B, C, and D
are the level with the largest SNR value, and the calculated V. COMPUTATIONAL RESULTS
SNRAi SNRBi SNRCi and SNRDi are shown in Table 3. SNR In this section, computational experiences with the MIP
values at the optimal level for each process parameter are in model in Section 2 and the proposed GA in Section 3 are
bold. According to Table 3, the optimal parameter settings of presented. All computational experiments are carried out on a
robust design are 200 for population size, 0.7 for the crossover PC with a Pentium IV CPU at 3.4 GHz and 2.0 GB RAM. The
rate, 75 for mutation crank-up interval, and Mutation MIP model is solved using CPLEX in an OPL-Studio envi-
Policy 2. ronment. The program for the proposed GA is implemented
From the contribution of each parameter on the per- in C++ programming language.
formance, process parameter A impacts 77%, and process
parameter D influences 15% on the performance of the pro- A. EXPERIMENT OF AN EXEMPLIFIED MDSDVRP
posed GA, as shown in Table 3. In other words, the large size A hypothetical test problem defined in Tables 6 and 7 has
of the population and the mutation policy of linear increase been generated to demonstrate the effectiveness of the
in mutation rate dominate other process parameters on the proposed GA. The problem has 6 customers with known
performance of the proposed GA. demands and two depots with 3 and 2 vehicles, respectively,
The process parameters optimized by the Taguchi method of heterogeneous capacities. The locations and demands of
are robust, so the signal or performance measure always customers are given in Table 4. Table 5 shows the locations
centralizes to the optimal expected values and is less affected of two depots and their vehicles with the given capacities.

VOLUME 9, 2021 112213

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

TABLE 4. Locations and demands of 6 customers.

TABLE 5. Vehicles available at two depots and their capacities.

TABLE 6. Results from the MIP model.

TABLE 7. Results from the proposed GA.

Customers and depots are on the 2-dimensional Euclidean distances as the maximum number of visits increases (See
space. Tables 6 and 7).
The optimal solutions from the MIP model presented in The best solutions from the proposed GA are shown for
Section 2 are shown for V = 1, 2, and 3, in Table 6. V = 1, 2, and 3 in Table 7. Tables 6 and 7 show that
The numbers in the parenthesis in Tables 6 and 7 are the the total traveled distances from the MIP model and the
unloaded amounts by the corresponding vehicles. As conjec- proposed GA are identical, and it indicates that the pro-
tured earlier, split deliveries to a customer reduce the total posed GA achieves the optimality for this small hypothetical

112214 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

FIGURE 6. Vehicle routes for one allowed visit (a), two allowed visits (b), and three allowed visits (c).

TABLE 8. Performance comparison for SDCVRP benchmark problems.

problem. Note that the routes of V13 and V22 in bold from B. EXPERIMENT OF BENCHMARK SDCVRPS
Tables 6 and 7 are identical, but the unloaded amounts at The CPLEX solver can only solve the MIP model of small
C6 are different. The proposed GA generated an alternative sizes because of the memory limitation. The proposed GA
solution. can generate good solutions for the MDSDVRPs of larger
From Tables 6 and 7, the proposed GA effectively solves sizes. In the next computational experiment, the proposed
the MDSDVRP with heterogeneous vehicles from multiple GA has been applied to benchmark VRPs available at
depots, allowing split deliveries. Note that the increase in the VRPLIB repository on the website (https://www.coin-
the maximum number of visits leads to a shorter travel or.org/SYMPHONY/branchandcut/VRP/data/index.htm.old).
distance. These problems have been widely used as benchmarks
The vehicle routes with different maximum numbers of ([18], [85]–[87]) and they are derived from Eilon et al. [88].
visits from the MIP model are illustrated in Figure 6. The The benchmark VRPs under comparison are SDCVRP,
circled Ci stands for customer i. Three vehicles V11 , V12 , and not allowing multiple visits. Since there is no benchmark
V13 are housed in depot 1 (D1 ) and two vehicles, V21 and V22 , problem for MDSDVRPs yet, we assessed the compet-
in depot 2 (D2 ). In Figure 6, the routes of different vehicles itiveness of our proposed GA by comparison with the
are represented by different types of arrows. Customer C6 has Eilon instances. Much research has used these benchmark
two vehicle visits for V = 2 in Figure 6(b) and three-vehicle instances to demonstrate the effectiveness of their algo-
visits for V = 3 in Figure 6(c). rithms ([89]–[91]), considering the NP-hardness of VRPs.

VOLUME 9, 2021 112215

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

TABLE 10. Demands of the retailer cities.

FIGURE 7. The MDSDVRP with 35 US cities and 3 depots.

TABLE 9. Vehicle capacities in each warehouse city.

TABLE 11. Results of the proposed GA for the VRP with one or two
allowed visit(s).

As explained in Section 1, many approximation algorithms

and heuristics keep updating the best solutions for these
benchmark instances. For some instances, those solutions
were matched with the upper bounds, proving that they are
optimal. depots, allowing split deliveries to customers and solving the
The effectiveness of the proposed GA for MDSDVRPs is classical SDCVRP well. Route of all vehicles in the bench-
already verified in Table 7. By limiting the proposed GA’s mark problems are reported in Appendix for the archival
ability with the maximum number of visits to a customer purpose.
set to 1 and a single depot, the effectiveness of the pro-
posed GA was compared with the best-known solution for C. EXPERIMENT OF A REAL-LIFE SCALE MDSDVRP
SDCVRP. A real-life scale MDSDVRP with heterogeneous vehicles
Note that the proposed GA is proposed to solve the het- from multiple depots, allowing split deliveries, is presented
erogeneous VRPs with multi-depot, allowing split-deliveries, and solved by the proposed GA. The problem has 35 US
but it also shows good effectiveness for the benchmark cities (customers), 3 depots, and 9 heterogeneous vehicles
SDCVRPs in Table 8. Table 8 shows the number of (3 vehicles at each depot). The distances between all cities
customers, the number of vehicles at a depot, the vehi- are the approximated driving distances on the road, obtained
cle capacity, the best-known solution by previous works, by database on Google Maps, instead of the Euclidean
and the solution obtained by the proposed GA. For this distances.
benchmark library, customers and a single depot are with The 38 nodes of the proposed problem are shown in Fig-
coordinates in a network, and Euclidean distances are ure 7, where the cities in red circles are 35 retailer cities and
used. the ones in blue rectangles are 3 warehouse cities located
Table 8 shows that the proposed GA is also effec- in Denver, Chicago, and Atlanta. There are three vehicles
tive in solving SDCVRPs and shows comparable perfor- in each warehouse city. The capacities of vehicles in each
mances for the benchmark problems. The proposed GA does warehouse city are given in Table 9. The demands of the
solve the VRP with heterogeneous vehicles from multiple retailer cities are given in Table 10.

112216 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

TABLE 12. Vehicle routes with two allowed visits. TABLE 13. Best solutions obtained by the proposed GA for the
benchmark problems.

FIGURE 8. The convergence of the proposed GA.

The results of the proposed GA for this problem with V =

1 and 2 are shown in Table 11. Note that routes with V = 2 are
better than ones with V = 2. The vehicle routes with V =
2 are given in Table 12. Two visits to customer 11 by vehicles
V22 and V33, shown in bold, contribute to the reduction of
the total traveled distance. All other customer’s demands are
satisfied only with one visit in Table 12. into a MIP formulation and tested to solve small problems.
The evolution of the best solution during 1776 generations As the motivation of this paper conjectures, it has been
by the proposed GA is shown in Figure 8. The graph shows identified that the introduction of split deliveries or multi-
the convergence of the best solution over the generations. ple visits to each customer leads to a reduction of delivery
After the 976th generation, the best solution of the problem cost.
is finally obtained. A GA to effectively and efficiently solve the medium
or large VRPs with heterogeneous vehicles from multiple
VI. CONCLUSION AND FUTURE RESEARCH depots, allowing split deliveries, has been proposed and vali-
A generalized VRP with heterogeneous vehicles from mul- dated successfully. The proposed algorithm has produced the
tiple depots, allowing split deliveries, has been identified solutions that are equal or close to the best-known solutions
by relaxing from the classical VRPs the constraints of for the benchmark SDCVRPs for which the proposed algo-
the number of depots and the number of visits allowed rithm has been executed with the restriction of one depot, one
to each customer. The identified VRP has been modeled allowed visit, and capacitated vehicles.

VOLUME 9, 2021 112217

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

Taguchi’s robust design method has been introduced and [10] A. Barvinok, S. P. Fekete, D. S. Johnson, A. Tamir, G. J. Woeginger, and
applied in optimizing the process parameters of the proposed R. Woodroofe, ‘‘The geometric maximum traveling salesman problem,’’
J. ACM, vol. 50, no. 5, pp. 641–664, Sep. 2003.
GA. Using the optimized parameters, the proposed GA shows [11] L. Engebretsen and M. Karpinski, ‘‘TSP with bounded metrics,’’ J. Com-
robust performance regardless of the size or structure of the put. Syst. Sci., vol. 72, no. 4, pp. 509–546, Jun. 2006.
problems. The proposed algorithm has effectively solved the [12] A. Ouaarab, B. Ahiod, and X.-S. Yang, ‘‘Discrete cuckoo search algorithm
for the travelling salesman problem,’’ Neural Comput. Appl., vol. 24,
test problem with 35 US cities and 3 depots, allowing multiple nos. 7–8, pp. 1659–1669, 2014.
visits and heterogeneous vehicles. [13] M. Mahi, Ö. K. Baykan, and H. Kodaz, ‘‘A new hybrid method based on
A new mutation policy has been developed for the pro- particle swarm optimization, ant colony optimization and 3-opt algorithms
posed GA. The existing GA implementations use the fixed for traveling salesman problem,’’ Appl. Soft Comput., vol. 30, pp. 484–490,
May 2015.
mutation rate. This thesis proposed the idea of self-adapting [14] J. Desrosiers, Y. Dumas, M. Solomon, and F. Soumis, ‘‘Time constrained
mutation rate, which enables dynamic speeds of evolutions routing and scheduling in: Network routing,’’ in Handbooks in Operations
in nature. The proposed mutation policy has proven effective Research and Management Science, vol. 8. Amsterdam, The Netherlands:
North Holland, 1995.
from the Taguchi method analysis and has generated consis- [15] I. H. Osman, ‘‘Vehicle routing and scheduling: Applications, algorithms
tently good solutions. and developments,’’ in Proc. Int. Conf. Ind. Logistics, Rennes, France,
As for future research, it may be helpful to investigate 1993, p. 277.
the issue where there is a restriction on the driving distance [16] G. Laporte, ‘‘The vehicle routing problem: An overview of exact and
approximate algorithms,’’ Eur. J. Oper. Res., vol. 23, pp. 631–640,
of vehicles available at each depot. The problem where the Jun. 1992.
customers have different time windows as their requirements [17] P. Toth and D. Vigo, ‘‘Exact algorithms for vehicle routing,’’ in
for delivered goods may also be worth considering. In addi- Fleet Management and Logistics. Boston, MA, USA: Kluwer, 1998,
pp. 1–31.
tion, future research can be conducted to improve the pro- [18] J. Lysgaard, A. N. Letchford, and R. W. Eglese, ‘‘A new branch-and-cut
posed algorithm and the performance comparison with other algorithm for the capacitated vehicle routing problem,’’ Math. Program.,
solution methods Additional improvements might lie in the vol. 100, no. 2, pp. 423–445, Jun. 2004.
combination of various selection and population replacement [19] M. Alssager, Z. A. Othman, M. Ayob, R. Mohemad, and H. Yuliansyah,
‘‘Hybrid cuckoo search for the capacitated vehicle routing problem,’’
schemes and new fitness models. Further investigation in Symmetry, vol. 12, no. 12, p. 2088, Dec. 2020.
optimizing the performance of the proposed algorithm and [20] M. Gendreau, G. Laporte, C. Musaraganyi, and É. D. Taillard, ‘‘A tabu
developing other solution methods to solve MDSDVRP can search heuristic for the heterogeneous fleet vehicle routing problem,’’
Comput. Oper. Res., vol. 26, no. 12, pp. 1153–1173, Oct. 1999.
be conducted to compare the comparison in performance and [21] H.-Y. Kang and A. Lee, ‘‘An enhanced approach for the multiple vehicle
the solution quality. Applications of the approach to related routing problem with heterogeneous vehicles and a soft time window,’’
problems can be explored as well. Symmetry, vol. 10, no. 11, p. 650, Nov. 2018.
[22] D. Li, Q. Cao, M. Zuo, and F. Xu, ‘‘Optimization of green fresh food
logistics with heterogeneous fleet vehicle route problem by improved
APPENDIX genetic algorithm,’’ Sustainability, vol. 12, no. 5, p. 1946, Mar. 2020.
See Table 13. [23] G. Laporte, Y. Nobert, and D. Arpin, ‘‘Optimal solutions to capacitated
multidepot vehicle routing problems,’’ Congressus Numeratium, vol. 44,
pp. 283–292, Dec. 1984.
REFERENCES [24] L. Zhishuo and C. Yueting, ‘‘Sweep based multiple ant colonies algorithm
for capacitated vehicle routing problem,’’ in Proc. IEEE Int. Conf. e-Bus.
[1] S. Hulagu and H. B. Celikoglu, ‘‘An electric vehicle routing problem with Eng. (ICEBE), Beijing, China, Oct. 2005, pp. 387–394.
intermediate nodes for shuttle fleets,’’ IEEE Trans. Intell. Transp. Syst.,
[25] I.-M. Chao, B. L. Golden, and E. Wasil, ‘‘A new heuristic for the multi-
early access, Sep. 23, 2020, doi: 10.1109/TITS.2020.3023673.
depot vehicle routing problem that improves upon best-known solu-
[2] S. Hulagu and H. B. Celikoglu, ‘‘A multiple objective formulation of tions,’’ Amer. J. Math. Manage. Sci., vol. 13, nos. 3–4, pp. 371–406,
an electric vehicle routing problem for shuttle bus fleet at a university Feb. 1993.
campus,’’ in Proc. 6th Int. Conf. Models Technol. Intell. Transp. Syst. (MT- [26] J. Renaud, F. F. Boctor, and G. Laporte, ‘‘An improved petal heuristic for
ITS), Cracow, Poland, Jun. 2019, pp. 1–5. the vehicle routeing problem,’’ J. Oper. Res. Soc., vol. 47, pp. 329–336,
[3] A. Ham and M.-J. Park, ‘‘Electric vehicle route optimization under Feb. 1996.
time-of-use electricity pricing,’’ IEEE Access, vol. 9, pp. 37220–37228, [27] Y. Shi, L. Lv, F. Hu, and Q. Han, ‘‘A heuristic solution method for multi-
2021. depot vehicle routing-based waste collection problems,’’ Appl. Sci., vol. 10,
[4] K. Dorling, J. Heinrichs, G. G. Messier, and S. Magierowski, ‘‘Vehicle no. 7, p. 2403, Apr. 2020.
routing problems for drone delivery,’’ IEEE Trans. Syst., Man, Cybern., [28] D. Hou, H. Fan, X. Ren, P. Tian, and Y. Lv, ‘‘Time-dependent multi-
Syst., vol. 47, no. 1, pp. 70–85, Jan. 2017. depot heterogeneous vehicle routing problem considering temporal–spatial
[5] R. G. Ribeiro, L. P. Cota, T. A. M. Euzebio, J. A. Ramirez, and distance,’’ Sustainability, vol. 13, no. 9, p. 4674, Apr. 2021.
F. G. Guimaraes, ‘‘Unmanned-aerial-vehicle routing problem with mobile [29] F. A. Tillman and R. W. Hering, ‘‘A study of a look-ahead procedure for
charging stations for assisting search and rescue missions in postdis- solving the multiterminal delivery problem,’’ Transp. Res., vol. 5, no. 3,
aster scenarios,’’ IEEE Trans. Syst., Man, Cybern., Syst., early access, pp. 225–229, Aug. 1971.
Jun. 22, 2021, doi: 10.1109/TSMC.2021.3088776. [30] R. Liu, Z. Jiang, R. Y. K. Fung, F. Chen, and X. Liu, ‘‘Two-phase heuris-
[6] G. B. Dantzig and J. H. Ramser, ‘‘The truck dispatching problem,’’ Man- tic algorithms for full truckloads multi-depot capacitated vehicle routing
age. Sci., vol. 6, no. 1, pp. 80–91, 1959. problem in carrier collaboration,’’ Comput. Oper. Res., vol. 37, no. 5,
[7] M. L. Fisher, ‘‘Optimal solution of vehicle routing problems using mini- pp. 950–959, May 2010.
mum K-trees,’’ Oper. Res., vol. 42, no. 4, pp. 626–642, Aug. 1994. [31] H. Bae and I. Moon, ‘‘Multi-depot vehicle routing problem with time win-
[8] E. L. Lawler, J. K. Lenstra, A. H. G. R. Kan, and D. B. Shmoys, The Trav- dows considering delivery and installation vehicles,’’ Appl. Math. Model.,
eling Salesman Problem. Chichester, U.K.: Wiley, 1985. vol. 40, nos. 13–14, pp. 6536–6549, Jul. 2016.
[9] J. Knox, ‘‘Tabu search performance on the symmetric traveling sales- [32] P. Stodola and J. Mazal, ‘‘Applying the ant colony optimisation algo-
man problem,’’ Comput. Oper. Res., vol. 21, no. 8, pp. 867–876, rithm to the capacitated multi-depot vehicle routing problem,’’ Int. J. Bio-
Oct. 1994. Inspired Comput., vol. 8, no. 4, pp. 228–233, 2016.

112218 VOLUME 9, 2021

H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

[33] L. D. Bodin, B. L. Golden, A. A. Assad, and M. Ball, ‘‘Routing and [59] G. Barbarosoglu and D. Ozgur, ‘‘A tabu search algorithm for the vehi-
scheduling of vehicles and crews, the state of the art,’’ Comput. Oper. Res., cle routing problem,’’ Comput. Oper. Res., vol. 26, no. 3, pp. 255–270,
vol. 10, no. 2, pp. 212–263, 1983. Mar. 1999.
[34] E. Falkenauer, ‘‘A hybrid grouping genetic algorithm for bin packing,’’ [60] B. M. Baker and M. A. Ayechew, ‘‘A genetic algorithm for the vehi-
J. Heuristics, vol. 2, pp. 5–30, Jun. 1996. cle routing problem,’’ Comput. Oper. Res., vol. 30, no. 5, pp. 787–800,
[35] M. R. Garey and D. S. Johnson, Computers and Intractability: A Guide to Apr. 2003.
the Theory of NP-Completeness. New York, NY, USA: Freeman, 1979. [61] K. Braekers, K. Ramaekers, and I. Van Nieuwenhuyse, ‘‘The vehicle
[36] L. Wang, R. Cai, M. Lin, and Y. Zhong, ‘‘Enhanced list-based simulated routing problem: State of the art classification and review,’’ Comput. Ind.
annealing algorithm for large-scale traveling salesman problem,’’ IEEE Eng., vol. 99, pp. 300–313, Sep. 2016.
Access, vol. 7, pp. 144366–144380, 2019. [62] Y. Rochat and É. D. Taillard, ‘‘Probabilistic diversification and intensi-
[37] Y. Cui, J. Zhong, F. Yang, S. Li, and P. Li, ‘‘Multi-subdomain grouping- fication in local search for vehicle routing,’’ J. Heuristics, vol. 1, no. 1,
based particle swarm optimization for the traveling salesman problem,’’ pp. 147–167, Sep. 1995.
IEEE Access, vol. 8, pp. 227497–227510, 2020. [63] I. H. Osman, ‘‘Metastrategy simulated annealing and tabu search algo-
[38] M. Wang, T. Ma, G. Li, X. Zhai, and S. Qiao, ‘‘Ant colony optimization rithms for the vehicle routing problem,’’ Ann. Oper. Res., vol. 41, no. 4,
with an improved pheromone model for solving MTSP with capacity and pp. 421–451, 1993.
time window constraint,’’ IEEE Access, vol. 8, pp. 106872–106879, 2020. [64] R. Cordone and R. W. Calvo, ‘‘A heuristic for the vehicle routing problem
[39] I. K. Singgih, J. Lee, and B.-I. Kim, ‘‘Node and edge drone surveillance with time windows,’’ J. Heuristics, vol. 7, pp. 107–129, Mar. 2001.
problem with consideration of required observation quality and battery [65] J. Berger and M. Barkaoui, ‘‘A new hybrid genetic algorithm for the
replacement,’’ IEEE Access, vol. 8, pp. 44125–44139, 2020. capacitated vehicle routing problem,’’ J. Oper. Res. Soc., vol. 41, no. 2,
[40] I. K. Singgih, O. Yu, B.-I. Kim, J. Koo, and S. Lee, ‘‘Production scheduling pp. 179–194, 2003.
problem in a factory of automobile component primer painting,’’ J. Intell. [66] G. D. Konstantakopoulos, S. P. Gayialis, and E. P. Kechagias, ‘‘Vehi-
Manuf., vol. 31, no. 6, pp. 1483–1496, Aug. 2020. cle routing problem and related algorithms for logistics distribution: A
[41] S. Lee and J. A. Ventura, ‘‘An algorithm for constructing minimal c- literature review and classification,’’ Oper. Res., pp. 1–30, Sep. 2020.
broadcast networks,’’ Networks, vol. 38, no. 1, pp. 6–21, 2001. [Online]. Available: https://link.springer.com/article/10.1007/s12351-020-
[42] I. K. Singgih and G. M. Lee, ‘‘Lower and upper bounds for military 00600-7#citeas, doi: 10.1007/s12351-020-00600-7.
deployment planning considering combat,’’ Int. J. Ind. Eng., Theory, Appl., [67] R. Elshaer and H. Awad, ‘‘A taxonomic review of metaheuristic algorithms
Pract., vol. 21, no. 6, pp. 408–420, 2014. for solving the vehicle routing problem and its variants,’’ Comput. Ind.
[43] F. B. Pereira, J. Tavares, P. Machado, and E. Costa, ‘‘GVR: A new genetic Eng., vol. 140, Feb. 2020, Art. no. 106242.
representation for the vehicle routing problem,’’ in Proc. 13th Irish Int. [68] C. J. Malmborg, ‘‘A genetic algorithm for service level based vehi-
Conf. Artif. Intell. Cogn. Sci., 2002, pp. 95–102. cle scheduling,’’ Eur. J. Oper. Res., vol. 93, no. 1, pp. 121–134,
[44] G. Laporte, M. Gendreau, J.-Y. Potvin, and F. Semet, ‘‘Classical and Aug. 1996.
modern heuristics for the vehicle routing problem,’’ Int. Trans. Oper. Res., [69] S. Jung and A. Haghani, ‘‘Genetic algorithm for a pickup and delivery
vol. 7, nos. 4–5, pp. 285–300, 2000. problem with time Windows,’’ Transp. Res. Rec., vol. 1733, no. 1, pp. 1–7,
[45] M. Desrochers, J. K. Lenstra, and M. W. P. Savelsbergh, ‘‘A classification Jan. 2000.
scheme for vehicle routing and scheduling problems,’’ J. Oper. Res. Soc., [70] D. K. Gupta, ‘‘A new algorithm to solve vehicle routing problems (VRPs),’’
vol. 46, pp. 322–332, Jun. 1990. Int. J. Comput. Math., vol. 80, no. 3, pp. 267–274, Mar. 2003.
[46] G. Laporte and Y. Nobert, ‘‘Exact algorithms for the vehicle routing [71] C. R. Reeves, Modern Heuristic Techniques for Combinatorial Problems.
problem,’’ Ann. Discrete Math., vol. 31, pp. 147–184, Jan. 1987. Oxford, U.K.: Blackwell Scientific Press, 1993.
[47] N. Christofides, A. Mingozzi, and P. Toth, ‘‘The vehicle routing problem,’’ [72] H. W. Shin and M. K. Kang, ‘‘A heuristic for vehicle routing problem
Combinatorial Optimization. Chichester, U.K.: Wiley, 1979, pp. 315–338. allowing multiple visits,’’ Soc. Korea Ind. Syst. Eng., vol. 14, no. 24,
[48] T. L. Magnanti, ‘‘Combinatorial optimization and vehicle fleet planning: pp. 141–147, 1991.
Perspectives and prospects,’’ Networks, vol. 11, pp. 179–214, Jun. 1981. [73] C. Archetti, N. Bianchessi, and M. G. Speranza, ‘‘A column generation
[49] N. Christofides, Vehicle Routing: The Traveling Salesman Problem: approach for the split delivery vehicle routing problem,’’ Networks, vol. 58,
A Guided Tour of Combinatorial Optimization. Chichester, U.K.: Wiley, no. 4, pp. 241–254, Dec. 2011.
1985, pp. 431–448. [74] L. Berbotto, S. García, and F. J. Nogales, ‘‘A randomized granular tabu
[50] B. L. Golden, E. A. Wasil, J. P. Kelly, and I. M. Chao, ‘‘Metaheuristics in search heuristic for the split delivery vehicle routing problem,’’ Ann. Oper.
vehicle routing,’’ in Fleet Management and Logistics. Boston, MA, USA: Res., vol. 222, no. 1, pp. 153–173, Nov. 2014.
Kluwer, 1998, pp. 33–56. [75] N. Bianchessi, M. Drexl, and S. Irnich, ‘‘The split delivery vehicle routing
[51] G. Laporte, ‘‘What you should know about the vehicle routing problem,’’ problem with time windows and customer inconvenience constraints,’’
Nav. Res. Logistics, vol. 54, no. 8, pp. 811–819, 2007. Transp. Sci., vol. 53, no. 4, pp. 1067–1084, Jul. 2019.
[52] E. Ghorbani, M. Alinaghian, G. B. Gharehpetian, S. Mohammadi, and [76] Z. Chen, M. Yang, Y. Guo, Y. Liang, Y. Ding, and L. Wang, ‘‘The
G. Perboli, ‘‘A survey on environmentally friendly vehicle routing problem split delivery vehicle routing problem with three-dimensional loading
and a proposal of its classification,’’ Sustainability, vol. 12, no. 21, p. 9079, and time Windows constraints,’’ Sustainability, vol. 12, no. 17, p. 6987,
Oct. 2020. Aug. 2020.
[53] W. K. Anuar, L. S. Lee, S. Pickl, and H.-V. Seow, ‘‘Vehicle routing optimi- [77] S. Ray, A. Soeanu, J. Berger, and M. Debbabi, ‘‘The multi-depot split-
sation in humanitarian operations: A survey on modelling and optimisation delivery vehicle routing problem: Model and solution algorithm,’’ Knowl.-
approaches,’’ Appl. Sci., vol. 11, no. 2, p. 667, Jan. 2021. Based Syst., vol. 71, pp. 238–265, Nov. 2014.
[54] H. Lim, ‘‘A genetic algorithm for the vehicle routing problem with [78] M. Gen and R. Cheng, Genetic Algorithms and Engineering Optimization.
heterogeneous vehicles from multiple depots, allowing multiple visits,’’ Hoboken, NJ, USA: Wiley, 2000.
M.S. thesis, Dept. Ind. Eng., Oregon State Univ., Corvallis, OR, USA, [79] G. Syswerda, ‘‘Uniform crossover in genetic algorithms,’’ in Proc. 3rd
Aug. 2007. Int. Conf. Genetic Algorithms, J. Schaffer, Ed. San Mateo, CA, USA:
[55] D. Gulczynski, B. Golden, and E. Wasil, ‘‘The multi-depot split delivery Morgan Kaufmann, 1989, pp. 2–9.
vehicle routing problem: An integer programming-based heuristic, new test [80] K. Q. Zhu, ‘‘A diversity-controlling adaptive genetic algorithm for the
problems, and computational results,’’ Comput. Ind. Eng., vol. 61, no. 3, vehicle routing problem with time windows,’’ in Proc. 15th IEEE Int. Conf.
pp. 794–804, Oct. 2011. Tools Artif. Intell., Sacramento, CA, USA, Nov. 2003, pp. 176–183.
[56] G. Laporte and F. Semet, ‘‘Classical heuristics for the vehicle routing [81] S. Wang, Z. Lu, L. Wei, G. Ji, and J. Yang, ‘‘Fitness-scaling adap-
problem,’’ GERAD, Montréal, QC, Canada, Tech. Rep. G-98-54, 1999. tive genetic algorithm with local search for solving the multiple depot
[57] E. Aarts and J. K. Lenstra, Local Search in Combinatorial Optimiza- vehicle routing problem,’’ Simulation, vol. 92, no. 7, pp. 601–616,
tion. Princeton, NJ, USA: Princeton Univ. Press, 2003. [Online]. Avail- Jul. 2016.
able: https://press.princeton.edu/books/paperback/9780691115221/local- [82] J. Li, F. Wang, and Y. He, ‘‘Electric vehicle routing problem with battery
search-in-combinatorial-optimization swapping considering energy consumption and carbon emissions,’’ Sus-
[58] J. F. Cordeau, M. Gendreau, G. Laporte, P. J-Y, and F. Semet, ‘‘A guide tainability, vol. 12, no. 24, p. 10537, Dec. 2020.
to vehicle routing heuristics,’’ J. Oper. Res. Soc., vol. 53, pp. 512–522, [83] G. Taguchi, S. Chowdhury, and S. Taguchi, Robust Engineering. New York,
May 2002. NY, USA: McGraw-Hill, 2000.

VOLUME 9, 2021 112219

 H. P. Lim et al.: Multi-Depot Split-Delivery Vehicle Routing Problem

[84] Y. Wu, Taguchi Methods for Robust Design. New York, NY, USA: ASME, GYU M. LEE received the B.S. and M.S. degrees
2000. in industrial engineering from Seoul National Uni-
[85] L. C. T. Gomes and F. J. Von Zuben, ‘‘Multiple criteria optimization based versity and the Ph.D. degree in industrial engi-
on unsupervised learning and fuzzy inference applied to the vehicle routing neering from The Pennsylvania State University.
problem,’’ J. Intell. Fuzzy Syst., vol. 13, pp. 143–154, Jan. 2002. He held a Professorship with the University of Ari-
[86] R. Fukasawa, H. Longo, J. Lysgaard, M. P. D. Aragão, M. Reis, E. Uchoa, zona, Northeastern University, the University of
and R. F. Werneck, ‘‘Robust branch-and-cut-and-price for the capacitated Illinois, Oregon State University, and POSTECH,
vehicle routing problem,’’ Math. Program., vol. 106, no. 3, pp. 491–511,
before he joined Pusan National University, where
May 2006.
he is currently a Professor with the Department
[87] K. Ganesh and T. T. Narendran, ‘‘CLOVES: A cluster-and-search heuristic
to solve the vehicle routing problem with delivery and pick-up,’’ Eur. J. of Industrial Engineering. His research interests
Oper. Res., vol. 178, no. 3, pp. 699–717, May 2007. include the operations research, combinatorial optimization, data analytics,
[88] S. Eilon, C. D. T. Watson-Gandy, and N. Christofides, Distribution Man- logistics and manufacturing, computational engineering, and information
agement: Mathematical Modelling and Practical Analysis. London, U.K.: technology. He is currently the Editor-in-Chief of International Journal of
Griffin, 1971. Industrial Engineering: Theory, Applications and Practice (SCI) and an
[89] A. M. Altabeeb, A. M. Mohsen, and A. Ghallab, ‘‘An improved hybrid Associate Editor of Decisions Science (SSCI).
firefly algorithm for capacitated vehicle routing problem,’’ Appl. Soft
Comput., vol. 84, Nov. 2019, Art. no. 105728.
[90] A. M. Altabeeb, A. M. Mohsen, L. Abualigah, and A. Ghallab, ‘‘Solving
capacitated vehicle routing problem using cooperative firefly algorithm,’’
Appl. Soft Comput., vol. 108, Sep. 2021, Art. no. 107403.
[91] I. Sbai, S. Krichen, and O. Limam, ‘‘Two meta-heuristics for solving IVAN KRISTIANTO SINGGIH received the
the capacitated vehicle routing problem: The case of the Tunisian B.S. and M.S. degrees in industrial engineering
Post Office,’’ Oper. Res., pp. 1–43, Jan. 2020. [Online]. Available: from Bandung Institute of Technology, Indone-
https://link.springer.com/article/10.1007/s12351-019-00543-8#citeas, sia, in 2009 and 2010, respectively, and the
doi: 10.1007/s12351-019-00543-8. Ph.D. degree in industrial engineering from Pusan
National University, Busan, South Korea, in 2017.
He was a Postdoctoral Researcher and a Research
Assistant Professor with the Department of Indus-
trial and Management Engineering, Pohang Uni-
HYUNPAE LIM received the B.S. degree from versity of Science and Technology (POSTECH).
Korea Military Academy and the M.S. degree in He is currently a Postdoctoral Researcher with the Department of Indus-
industrial engineering from Oregon State Univer- trial and Systems Engineering, Korea Advanced Institute of Science and
sity under Prof. Gyu M. Lee’s guidance. He is cur- Technology (KAIST), South Korea. His research interests include container
rently a Lieutenant Commander of Korean Army. terminal logistics, vehicle routing problems, scheduling, and serious game
He served for the logistics operations in Korean development. He is a member of KIIE, T-LOG, and Korean AI Association.
Army. His research interests include logistics, mil- He received the Best Student Paper Award at the IJIE 2013 Conference
itary defense operations, metaheuristics, and oper- and the Korean Government Scholarship Excellent Academic Achievement
ations research. Award in 2014.

112220 VOLUME 9, 2021