https://doi.org/10.31449/inf.v45i1.2814 Informatica 45 (2021) 93–104 93
Penalty Variable Neighborhood Search for the Bounded Single-Depot Multiple
Traveling Repairmen Problem
Ha-Bang Ban
School of Information and Communication Technology
Hanoi University of Science and Technology, Hanoi, Vietnam
E-mail: BangBH@soict.hust.edu.vn
Keywords: nounded-mTRP, penalty variable neighborhood search, metaheuristic
Received: June 4, 2019
Multiple Traveling Repairmen Problem (mTRP) is a class of NP-hard combinatorial optimization problems
with many practical applications. In this paper, a general variant of mTRP, also known as the Bounded
Single-Depot Multiple Traveling Repairmen Problem (Bounded-mTRP), is introduced. In the Bounded-
mTRP problem, a fleet of identical vehicles is dispatched to serve a set of customers. Each vehicle that
starts from the depot is only allowed to visit the number of customers within a predetermined interval, and
each customer must be visited exactly once. Such restrictions appear in real-life applications where the
purpose is to have a good balance of workloads for the repairmen. The goal is to find the order of customer
visits that minimizes the sum of waiting times. In our work, the proposed algorithm is encouraged by the
efficiency of the algorithms in [15, 19, 20] that are mainly based on the principles of the VNS [14]. The
penalty VNS extends the well-known VNS [14] by including constraint penalization, to solve the Bounded-
mTRP effectively. Extensive numerical experiments on benchmark instances show that our algorithm
reaches the optimal solutions for the problem with 76 vertices at a reasonable amount of time. Moreover,
the new best-known solutions are found in comparison with the state-of-the-art metaheuristic algorithms.
Povzetek: Razvita je nova metoda za preiskovanje grafov - za reševanje naloge serviserja, tipiˇ cnega NP-
polnega problema.
1 Introduction
1.1 Motivation and definition
The Traveling Repairman Problem (TRP) has been stud-
ied in the number of previous work [1, 2, 5, 15, 19, 20]. It
is known as the Minimum Latency Problem (MLP), or the
Deliveryman Problem (DMP). These problems arise appli-
cations, e.g., whenever repairmen or servers have to ac-
commodate a set of requests to minimize their total (or
average) waiting times [1, 2, 5, 15, 19, 20]. A direct gen-
eralization of the TRP is the Multiple Traveling Repair-
men Problem (mTRP) that considers k vehicles simultane-
ously. Applications of the mTRP can be found in Routing
Pizza Deliverymen, or Scheduling Machines to minimize
mean flow time for jobs. Several prior studies that we can
find in the literature are [10, 12]. In this paper, we study
the Bounded Single-Depot Multiple Traveling Repairmen
Problem (Bounded-mTRP) by involving the restriction of
the number of vertices that a repairman must visit in his
tour. The restriction is defined by lower (denoted by K)
and upper (denoted by L) bounds regarding the traveled
vertices. Therefore, the number of vertices that a repair-
man can visit lies within a predetermined interval with the
aim of obtaining balanced solutions. The requirement of
the problem is to find a tour such that the above restric-
tion is satisfied, and the overall cost of visiting all vertices
is minimized. Such restriction appears in many real-life
applications whose purpose is to have a good balance of
workloads for the repairmen.
1.2 Approach and contributions
There are three approaches for solving the Bounded-
mTRP: 1) exact algorithms, 2) approximation algorithms,
and 3) heuristic algorithms. The exact algorithms find the
optimal solution with an exponential time in the worst case.
Therefore, the exact algorithm only solves the problem
with small sizes. To describe related works, we denote an
approximation algorithm asp-approximation when the al-
gorithm finds the solution at mostp times worse than the
optimal solution. Here p is an approximation ratio with
a constant value. In this approach, the best approxima-
tion ratio of 16.994 is for the mTRP [10, 12]; however, it
is still far from the optimal solution. Heuristic algorithms
perform well in practice, and their efficiency can be eval-
uated through experiments. Our algorithm falls into this
approach.
Previously, research on the Bounded-mTRP has not
studied much, and this work presents the first metaheuris-
tic approach for this problem. Our algorithm is encouraged
by the efficiency of the algorithms in [14, 19, 20] that are
mainly based on the principles of the VNS [14]. However,

94 Informatica 45 (2021) 93–104 H.-B. Ban
the difference between the Penalty VNS (P-VNS) and their
VNS is that our algorithm builds up penalty value during a
search. The proposed algorithm includes two phases. The
algorithm is developed based on the GRASP [9] to build an
initial solution in the construction phase. In the improve-
ment phase, the P-VNS combined with shaking techniques
not only exploits good local solution space but also pre-
vents the search from escaping from local optimal. More-
over, several novel neighborhoods’ structure as well as a
constant time operation for calculating the cost of each
neighboring solution is also introduced. The main problem
is that there exists no other metaheuristic reported in the
literature for this problem; this is, we found no previous at-
tempts to solve this problem, neither exact nor heuristically,
to compare with. Therefore, we adapt the metaheuris-
tic algorithms in [17] to solve the Bounded-mTRP, and
choose several state-of-the-art metaheuristic algorithms for
the mTRP [8, 18], and Bounded-mTSP [17] as a baseline in
our research. Extensive numerical experiments on bench-
mark instances show that our algorithm reaches the optimal
solutions for the problems with up to 76 vertices at a rea-
sonable amount of time. Moreover, the new best-known
solutions are found in comparison with the state-of-the-art
metaheuristic algorithms.
The rest of this paper is organized as follows. Section
2, and 3 present literature review, and neighborhood struc-
ture, respectively. The proposed algorithm is described in
Section 4. Computational evaluations and discussions are
reported in Section 5. Finally, Section 7 concludes the pa-
per.
2 Literature review
The Bounded-mTRP has, as we know, not been studied
much, although it is a natural extension of the mTRP prob-
lem. In the literature, several variants of the problem are
introduced as follows:
– The mTRP is a popular case since no constraint is
considered. Numerous works for the mTRP can be
found in [8, 10, 12, 18]. Some metaheuristic algo-
rithms [10, 12] can give good solutions fast for large
instances.
– The mTRP with Profits (mTRPP) finds a travel plan
for server that maximizes the total revenue. Meta-
heuristic algorithm [3] produces solutions well.
– Another variant of the mMLP is mMLP with distance
constraints [4, 16]. Lou et al. [16] proposed an exact
algorithm that reaches the optimal solutions for the in-
stances with up to 50 vertices. Ban et al. [4] then
presented a metaheuristic algorithm based on VNS.
The experimental results concluded that the algorithm
found good-quality solutions for small and medium-
size instances.
– mTRPD [6] finds a tour with minimum latency sum in
post-disaster road clearance. Unlike mMLP, in disas-
ter situations, travel costs need to be added to debris
removal times. Their metaheuristic obtained the opti-
mal or near-optimal solutions on Istanbul data within
seconds.
– The TRP is a particular case where there is only
a repairman. Numerous metaheuristic algorithms
[1, 2, 5, 15, 19, 20] for the problem have proposed in
the literature. The experimental results showed that
their algorithms obtain good solutions fast for the in-
stances with up to 500 vertices.
These algorithms are the best algorithms for some variants
of the Bounded-mTRP problem. However, they do not in-
volve the bounded constraint. Therefore, they cannot be
used directly to solve the Bounded-mTRP.
3 Mathematical formulation
The formulation is obtained from the formulation proposed
by Christofides et al. [7] for the Capacitated Vehicle Rout-
ing Problem (CVRP). Letu
i
be a non-negative real variable
representing the length of the route from the depot 0 to the
vertex i. Letx
r
ij
be the following binary variables:
x
r
ij
=
(
1 if if edge(i;j) is used in router
0 otherwise
minz =
P
n
i=1
u
i
Subject to:
n
X
i=0
i6=j
m
X
r=1
x
r
ij
= 1; (j = 1; 2;:::;n)
(1)
n
X
i=0
i6=l
x
r
il
 n
X
j=0
j6=l
x
r
lj
= 0; (l = 0; 1;:::;n;
r = 1; 2;:::;m) (2)
n
X
j=1
x
r
0j
= 1; (r = 1; 2;:::;m) (3)
K n
X
i=0
n
X
j=0
x
r
ij
+ 1 L; (r = 1; 2;:::;m) (4)
u
i
 u
j
+ (T +c
ij
) m
X
r=1
x
r
ij
+(T c
ji
) m
X
r=1
x
r
ji
 T ;
(i;j = 1;:::;n;i6=j) (5)

Penalty Variable Neighborhood Search for. . . Informatica 45 (2021) 93–104 95
u
i
 c
0i
 m
X
r=1
x
r
0i
; (i = 1; 2;:::;n) (6)
x
r
ij
2f0; 1g; (i;j = 0; 1;:::;n;j6=i;
r = 1;:::;m) (7)
u
i
 0; (i = 1;:::;n) (8)
Constraints (1) show that each vertex is contained in only
one route. Constraints (2) indicate that when a vertex is in
a route, then it has a predecessor and a successor in that
route. Constraints (3) ensure that one vertex is sequenced
as first in each route and constraints (4) guarantee that and
the number of vertices visited of each repairman must less
than L and more than K. Constraints (5) demonstrate that
u
j
=u
i
+c
ij
when
P
m
r=1
x
r
ji
= 1 and they are redundant
when
P
m
r=1
x
r
ji
= 0. Constraints (6) initialize the latency
of the first vertex in each route. Finally, constraints (7) and
(8) establish the nature of the variables.
4 Neighborhood structure
Seven neighborhoods investigated are divided into
two categories: intro-route and intra-route. Now, let
T = (R
1
;R
2
;:::;R
l
;:::;R
k
) (l = 1;:::;k) be a tour, we
introduce a novel neighborhoods’ structure and complexity
of their exploration. Note that: in [15, 20], the complexity
of some neighborhoods is already mentioned. Therefore,
we only introduce the complexity of new ones.
For Intro-route: Intro-route is used to optimize on a
single route. Assume that, R and m (m < n) are a
route and its length, respectively. We then introduce five
neighborhoods’ structure in turn.
Remove-insert neighborhood considers each vertex v
i
in the route at the end of it. This neighborhood of R is
defined as a set N
1
(R) =fR
i
= (v
1
;v
2
;:::;v
i 1
;v
i+1
;
:::;v
m
;v
i
) : i = 2; 3;:::;m 1g. Obviously, the size of
N
1
(R) isO(m).
Property 1. The time complexity of exploring N
1
(R) is
O(m
2
).
Swap adjacent neighborhood attempts to swap each pair
of adjacent vertices in the route. This neighborhood of
R is defined as a set N
2
(R) = fR
i
= (v
1
;v
2
;:::;v
i 2
;
v
i
;v
i 1
;v
i+1
;:::;v
m
) : i = 3; 4;:::;m 1g. The size of
the neighborhood isO(m).
Property 2. The time complexity of exploring N
2
(R) is
O(m).
Swap neighborhood attempts to swap the po-
sitions of each pair of vertices in the route.
This neighborhood of R is defined as a set
N
3
(R) = fR
ij
= (v
1
;v
2
;:::;v
i 1
;v
j
;v
i+1
;:::;v
j 1
;
v
i
;v
j+1
;:::;v
m
) : i = 2; 3;:::;m 3;j = i + 3;:::;mg.
The size of the neighborhood isO(m
2
).
Property 3. The complexity of exploring N
3
(R) is
O(m
2
).
3-opt neighborhood attempts to reallocate three ad-
jacent vertices to another position of the route.
This neighborhood of R is defined as a set
N
4
(R) = fR
i
= (v
1
;v
2
;:::;v
i 1
;v
i
;v
j+1
;:::;v
k
;
v
i+1
;:::;v
j
;v
k+1
::::;v
m
) : i = 2; 3;:::;m 5;j =
4;:::;m 3;k = 6;:::;m 1g. The size of the neighbor-
hood isO(m
3
).
Property 4. The complexity of exploring N
4
(R) is
O(m
3
).
2-opt neighborhood removes each pair of edges from
the solution and reconnects the vertices. This neigh-
borhood of T is defined as a set N
5
(T ) = fT
ij
=
(v
1
;v
2
;:::;v
i
;v
j
;v
j 1
;:::; v
i+2
;v
i+1
;v
j+1
;:::;v
m
) : i =
1;:::;n 4;j =i+4;:::;mg. The size of the neighborhood
isO(m
2
).
Property 5. The complexity of exploring N
5
(T ) is
O(m
2
).
It is realized that the calculation of a neighboring solution’s
cost by using the known cost of the current solution can be
done in constant time [15, 20]. As a result, the algorithm
spendsO(m
3
) operations for a full neighborhood search.
For intra-route: LetR
l
,R
h
,ml, andmh be two different
routes and their sizes in T , respectively. Intra-route is
used to exchange vertices between two different routes
or remove vertices from a route and then insert them to
another as followings:
The swap-intra-routes neighborhood tries to exchange
the positions of each pair of vertices in R
l
and R
h
in
turn. The neighborhood of R
l
and R
h
is defined as a
set N
8
(T ) = fT
i
= (R
1
;:::;R
2
;:::;R
l
= (v
1l
;v
2l
;:::;
v
ih
;v
il+1
;:::;v
ml
);:::;R
h
= (v
1h
;v
2h
;:::;v
il
;v
ih+1
;
:::;v
mh
);:::;R
k
) : il = 2; 3;:::;ml 1;ih = 2; 3;
:::;mh 1g. The size of the neighborhood isO(ml mh).
Property 6. The complexity of exploring N
6
(T ) is
O(ml mh).
Proof. For a tourT2N
6
(R), we have
L(Ti) = L(T) (ml i+1)c(v
il 1
;v
il
)
 (ml i)c(v
il
;v
il+1
)
 (mh i 1)c(v
ih 1
;v
ih
)
 (mh i)c(v
ih
;v
ih+1
)
+(ml i+1)c(v
il 1
;v
ih
)
+(ml i)c(v
ih
;v
il+1
)
+(mh i+1)c(v
ih 1
;v
il
)
+(mh i)c(v
il
;v
ih+1
): (9)
Hence, we can calculateL(T
i
) by the formulation (10) in
O(1) time. Therefore, the complexity of exploringN
6
(T )
isO(ml mh).
The insert-intra-routes neighborhood considers each
vertex v
i
in R
l
and insert it into each position
in R
h
. The neighborhood of R
l
and R
h
is de-
fined as a set N
7
(T ) = fT
i
= (R
1
;:::;R
l
=
(v
1l
;v
2l
;:::;v
ih 1
;v
ih
;v
il+1
;:::;v
ml
);:::;R
h
= (v
1h
;
v
2h
;:::;v
ih 1
;v
ih+1
;:::;v
mh
);:::;R
k
) :il = 2; 3;:::;ml

96 Informatica 45 (2021) 93–104 H.-B. Ban
1;ih = 2; 3;:::;mh 1g. The size of the neighborhood is
O(ml mh).
Property 7. The complexity of exploringN
7
(T ) isO(ml mh).
Proof. For a tourT2N
7
(R), we have
L(Ti) = L(T) (ml i)c(v
il
;v
il+1
)
 (mh i+1)c(v
ih 1
;v
ih
)
 (mh i)c(v
ih
;v
ih+1
) ih 1
X
k=1h
c(v
k
;v
k+1
)
+(ml i)c(v
il
;v
ih
)
+(mh i 1)c(v
ih
;v
il+1
)
+(mh i+1)c(v
ih 1
;v
ih+1
)
+
il 1
X
k=1l
c(v
k
;v
k+1
): (10)
Hence, we can calculateL(T
i
) by the formulation (11)
in O(max(mh;ml)) time. Therefore, the complexity of
exploringN
7
(T ) isO(max(mh;ml) ml mh).
5 Algorithmic design
The proposed algorithm includes two phases as follows:
In the construction phase, the algorithm [9] allows a con-
trolled amount of randomness to overcome the behavior of
a purely greedy heuristic. It is used to build an initial so-
lution for our algorithm. In the improvement phase, the
penalty VNS [14] is combined with shaking operators to
escape from local optima. The proposed algorithm is re-
peated a number of times, and the best solution found is
reported. An outline of the algorithm is shown in Algo-
rithm 1. In Step 1, the algorithm starts with an initial so-
lution. In Step 2, it is explored switches between different
neighborhoods. To explore new promising solution spaces,
a diversification step is added in Step 3. In the remaining
of this section, more details about the three steps of our
algorithm are given.
5.1 Feasible solution space
Penalty method is a technique to solve optimization prob-
lem with constraints. It adds a penalty value to the original
objective function. The advantage of penalty technique is
simple to implement. It is used to solve successfully many
problem [21]. All infeasible solutions are penalized by a
value. With a tourT , letV (T ) be the violation. The viola-
tion valueV (T ) is computed as follows:
V (T ) =
k
X
l=1
maxfjR
l
j L; 0g +
k
X
l=1
maxfKj R
l
j; 0g:
Solutions are then evaluated according to the weighted fit-
ness functionL
0
(T ) = L(T ) +  V (T ), where is the
penalty parameter
Algorithm 1 The Proposed Algorithm
Input: v1;V;Ni(T)(i=1;:::;7);level are a starting vertex, the
set of vertices inKn, the set of neighborhoods and the number
of swap, respectively.
Output: The best solutionT
 .
Step 1 (the construction phase):
{Initially,T is an empty tour}
repeat
T = ;
for(l =1;l<k;l++) do
R
l
=R
l
[v1; {main depot isv1}
while (all vertices are not visited) do
{Pick a random route that still satisfies the constraint if
the new insertion is occurred}
R =fR
l
jR
l
2T andL 1j R
l
j&&jR
l
j K 1g;
if9!R
l
then
R = Choose a random route(R
l
2T) with minimum
cost; {accept an invalid route}
Create a RCL ofve;{ve is the last vertex ofR
l
} Select a
randomly vertexv =fvijvi2RCL andvi is not visited
g to add toR
l
;
for(l =1;l k;l++) do
T =T[R
l
;{update the tourT }
LT =LT[T ;{LT is stored a list of solutions}
until iter
T = The best feasible solution if any. Otherwise the best infea-
sible one in LT;
while stop criteria not met do
Step 2 (the improvement phase):
fori:1!7 do
T
0
 argmin
T
00
2N
i
(T)
L(T
00
)
if ((L(T
0
)<L(T)) or(L(T
0
)<L(T
 ))) then
T T
0
if (L(T
0
)<L(T
 )) and(T
0
is feasible) then
T
 T
0
else
i++
Step 3 (Diversification):
type = rand(2);{Select randomly a number from 1 to 2}.
if type==1 then
R
l
= Select randomly a route2T ;
R
l
= shaking-single-route(R
l
;level);
else
(R
l
;R
h
) = Select randomly two routes ofT ;
(R
l
;R
h
)= shaking-multi-routes(R
l
;R
h
;level) ;
UpdateT ;
return T
 ;
Algorithm 2 shaking-single-route(R
l
;level)
Input: R
l
;level are thel th route, and the number of swap,
respectively.
Output: a new solutionR
l
.
while(level>0) do
selecti;j positions fromR
l
at random
if(i6=j) then
InsertR
l
[i] betweenR
l
[j] andR
l
[j+1];
level  level 1;
return R
l
;

Penalty Variable Neighborhood Search for. . . Informatica 45 (2021) 93–104 97
Algorithm 3 shaking-multi-routes(R
l
;R
h
;level)
Input: R
l
;R
h
;level are thel th,h th route, and the number
of swap, respectively.
Output: a new solutionR
l
andR
h
.
while(level>0) do
selecti th andj th positions fromR
l
andR
h
at random,
respectively;
swapR
l
[i] betweenR
h
[j];
level  level 1;
return R
l
andR
h
;
5.2 The construction phase
Our construction phase is developed on the GRASP
scheme in [9]. Only one iteration is performed, and one
solution is found which is either feasible or not. Its steps
is described in Algorithm 1. All routes are initialized with
v
1
because it is a starting vertex. Each vertex of K
n
is
then added to the tour by using a Restricted Candidate List
(RCL). The RCL of each vertex includes a number of ver-
tices that are the closest to it. At an iteration, we find a
routeR
l
that does not violate the constraint if a new inser-
tion occurs. Otherwise, if we cannot find any route, then
we accept the infeasibility. That means a route R
l
with
minimum cost in the tour is picked. Letv
e
be the current
last vertex of the route R
l
. An unvisited vertex v is then
picked randomly from the RCL ofv
e
to add toR
l
. A solu-
tion is generated when all vertices ofK
n
are routed. The
above steps are executed iter times to create iter solutions.
They are stored in a LT list. The procedure then returns
the feasible solution with minimum cost in the list if any.
If it cannot produce any feasible solution, the solution with
minimum cost is penalized by adding a value to the objec-
tive function.
5.3 The improvement phase
For a given current solutionT , the neighborhood explores
the neighboring solution space set N(T ) of T iteratively
and tries to replaceT by the best solutionT
0
2N(T ). The
main operation in exploring the neighborhood is the calcu-
lation of a neighboring solution’s cost. In straightforward
implementation, this operation requires Tsol = O(n).
However, by using the known cost of the current solution,
we show that this operation can be done in constant time
for considered neighborhoods. Thus, we speed up the run-
ning time of exploring these neighborhoods.
In a preliminary study, we realize that the efficiency of
VNS algorithm relatively depends on the order in which
the neighborhoods are used. Therefore, the neighborhoods
are explored in a specific order based on the size of their
structure, namely, from the small to large, such as the swap-
adjacent, remove-insert, swap, 2-opt, or, swap-intra-route,
and insert-intra-route. The time complexity of exploring
the neighborhoods is reduced by choosing a random ver-
tex, and then we are only interested in neighborhoods gen-
erated from this vertex’s moves. The strategy is called “re-
stricted". As a result, the size of the neighborhood is re-
duced by a factor O(n). Another is “without restricted"
strategy when the entire neighborhood is explored without
first fixing a random vertex. The reduction of neighborhood
size is also used in [19]. The aim of using two strategies is
to introduce several options to run the proposed algorithm
effectively.
5.4 Diversification
Shaking procedure allows to guide the search towards an
unexplored part of the solution space. In this work, two
types of shaking are used to give a new solution: shaking
in a single route (shaking-single-route) and shaking in two
(shaking-multi-routes). In shaking procedure in a single
route, it selects thel-th routeR
l
ofT and then swaps ran-
domly several vertices for each other. In the rest, it picks
two routesR
l
andR
h
in a random manner, and after that,
exchanges randomly some several vertices in them. We fi-
nally return to Step 2 with the new solution. The shaking
procedure is described in Algorithm 2 and 3.
The last aspect to discuss is the stop criterium of our al-
gorithm. A balance must be made between computation
time and efficiency. Here, the algorithm stops if no im-
provement is found after the number of loop (NL).
5.5 The time complexity
The running time of our algorithm mainly spends on ex-
ploring in VNS. In the VNS step, insert-intra-routes neigh-
borhood consumes time, at least as well as the others. As-
sume that if these neighborhoods are invokedk
1
times, then
the complexity of neighborhoods’ exploration is O(k
1
 max(mh;ml) ml mh) O(k
1
 n
3
) (in the worst
case the size of mh or ml is n). It is also the theoretical
complexity of our algorithm.
6 Computational results
The experiments are conducted on a personal computer,
which is equipped with an Intel Pentium core i7 duo 2.10
Ghz CPU and 4 GB bytes RAM memory.
6.1 Datasets
The numerical analysis is performed on a set of benchmark
problems for the mTRP and Bounded-mTSP [8, 17, 18].
As testing our algorithm on all instances would have been
computationally too expensive, we implement our numeri-
cal analysis of some selected instances. In [17], R. Necula
et al. propose the Bounded-mTSP instances based on the
TSPLIB benchmark. Specifically, they transform TSPLIB
four instances (eil51, berlin52, eil76, and rat99) by setting
the number of salesmen to be, by turn, 2, 3, 5, and 7. With
the aim of obtaining balanced solutions, they choose to set
the bounds (L and K) on the number of vertices in each

98 Informatica 45 (2021) 93–104 H.-B. Ban
route, by running the k-means clustering algorithm. Be-
sides that, we add more real instances by randomly choos-
ing some instances from TSPLIB. We divide the instances
into two groups: 1) in group one (G1), the vertices are con-
centrated; 2) in the other (G2) the vertices are scattered.
6.2 Metrics
To evaluate our algorithm’s solution quality, we need to
compare it with the other metaheuristics. The main prob-
lem is that there exists no other metaheuristic reported
in the literature for this problem. That means we found
no previous attempts to solve the problem, neither ex-
act nor heuristic (or metaheuristic), to compare. We try
to solve exactly several small instances using some state-
of-art solvers and use those results to evaluate the per-
formance of the proposed algorithm. However, the ap-
proach only solve the problem with small instances, while
metaheuristics is a suitable approach for the problem with
large sizes. Therefore, we adapt the existing algorithms
in [17] to compare with the proposed algorithm for the
Bounded-mTRP. We define the improvement of our algo-
rithm with respect to Best.Sol (Best.Sol is the best solu-
tion found by our algorithm) in comparison with the ini-
tial solution (Init.Sol), upper bound (UB) obtained by the
GRASP and the adapted algorithm in [17], respectively.
Improv[%] =
Best:Sol Init:Sol
Init:Sol
 100%, andGap[%] =
Best:Sol UB
UB
 100%. In addition, we choose several state-
of-the-art metaheuristic algorithms for the Bounded-mTSP
[17] (Bounded Multiple Traveling Salesman Problem) and
mTRP (Multiple Traveling Repairmen Problem) in [8, 18]
as a baseline in our research.
6.3 Results and discussions
Through preliminary experiments, we observe that the val-
uesiter = 10; = 10; = 5, level=5, PF=100, and NL
= 100 resulted in a good trade-off between solution quality
and run time. In this paper, the neighborhoods’ order is as
follows: swap adjacent, remove-insert, swap, 2-opt, or-opt,
swap-intra-routes, and insert-intro-routes. These settings
have thus been used in the following experiments.
In the tables, Init.Sol, Best.Sol, Aver.Sol, and T corre-
spond to the initial, best, and average solution, and the av-
erage time in seconds of ten executions obtained by our
algorithm, respectively. The column ACS in Tables 1 and
2 describe the best results obtained from the adapted algo-
rithms in [17]. Figure 1 and 2 shows the evolution of the av-
erage improvement in two strategies. The values in Figures
are extracted from Table 4. In Table 5, kM-ACS, g-ACS,
s-ACS, gb-ACS, and sb-ACS [17] are developed on Ant
Colony System (ACS) with different strategies. The pro-
posed algorithm is tested by selecting a fixed random ver-
tex that is labeled “restricted". Runs in which the search ex-
plores all possible moves are labeled “without restricted".
6.3.1 Experimental results for the Bounded-mTRP
The experimental results in Table 3 are the average values
calculated from Table 1 and 2. In Table 3, for all instances,
it can be observed that our algorithm is capable of improv-
ing the solutions in comparison with Init.Sol. The average
improvement of our algorithm with the two strategies is
about 16:94% and 13:69%, respectively. Obviously, our al-
gorithm can obtain a significant improvement for almost
instances and required small-scaled running time. Both
strategies seem to work well. The neighborhood imple-
mentation with fixed random vertex uses significantly less
computing time, combined with a slight loss of solution
quality (about 3:25%). However, the strategy proves use-
ful for the larger instances, for which full neighborhood
search is too time-consuming. Moreover, in comparison
with the algorithms in [17], the proposed algorithm also
outperforms for most instances.
From Tables 5 to 6 we can draw some conclusions about
the working of our algorithm. Unsurprisingly, the multi-
start version of our algorithm (algorithm settings 5 to 8) re-
quires a much larger computation time than the single-start
version (settings 1 to 4). However, the quality improvement
obtained by this method is relatively small. The perturba-
tions in the GRASP+VNS algorithm (Table 6) seem to help
marginally, as the solutions obtained by this algorithm are
usually slightly better. This may indicate that the GRASP
multi-start is not able to provide enough diversification, and
that the perturbation move is useful.
For two strategies, Figure 1 and 2 shows the evolution of
the average deviation to the initial solutions with respect to
improv andT during the iterations in some instances. The
deviations in two strategies are 14.61% (10.38%), 16.25%
(11.63%), 16.48% (11.79%), 16.66% (11.94%), 16.77%
(12.03%), 16.94% (13.69%), and 16.94% (13.69%) for
the first local optimum, obtained by one, ten, twenty,
thirty, fifty, one-hundred, and two-hundred iterations, re-
spectively. A major part of the descent obtained by from
fifty to one-hundred iterations. As can be observed, addi-
tional iterations give a minor improvement with the large
running time. Hence, the first way to reduce the large run-
ning time is to use no more than one-hundred iterations,
and the improvement of the proposed algorithm is about
16.94% (13.69%) for two strategies, respectively. A much
faster option is to run the initial construction phase then
improve it by using a single iteration, which obtains an av-
erage deviation of 14.61% (10.38%) and an average time
of 0.28 (0.21) seconds.
6.3.2 Experimental results for some variants
To the best of our knowledge, most algorithms are devel-
oped for a specific variant that is not applicable to other
variants. Our algorithm can be applicable to the Bounded-
mTSP, although it was not designed for solving them. In
comparison with the state of the art algorithms for the
Bounded-mTSP, and mTRP in [8, 17, 18], our algorithm’s
solutions are better than the other algorithms. Specifically,

Penalty Variable Neighborhood Search for. . . Informatica 45 (2021) 93–104 99
 Table 1: The experimental results for Bounded-mTRP without restricted neighborhood.

100 Informatica 45 (2021) 93–104 H.-B. Ban
 Table 2: The experimental results for Bounded-mTRP with restricted neighborhood.

Penalty Variable Neighborhood Search for. . . Informatica 45 (2021) 93–104 101
 Gap Time Gap Time
 3.84 2.45
 5.46 4.46 4.64 3.45 Table 3: The average results for two schemes.
 Init.Sol Improv T Improv T Improv T Improv T Improv T Improv T Improv T
 3 0.64 0.93 1.55 3.22 3.84 14.15 0.33 0.93 1.35 2.24 4.47 5.46 19.81
 0.28 0.78 1.14 1.89 3.85 7.65 16.99
 0.15 0.41 0.60 0.99 2.05 2.45 9.03
 0.27 0.76 1.11 1.83 3.65 4.46 16.17
 0.21 0.58 0.85 4 2.85 3.45 12.60 Table 4: Evolution of average deviation to Init.Sol without restricted neighborhood.
Figure1.Evolutionof averageImprove[%]
Figure2.Evolutionof averagerunningtime
1 10 20 30 50 100 200
withoutrestricted 14.61 16.25 16.48 16.66 16.77 16.94 16.94
restricted 10.38 11.63 11.79 11.94 12.03 13.69 13.69
0
4
8
12
16
20
Improve[%]
Iterations
1 10 20 30 50 100 200
withoutrestricted 0.28 0.78 1.14 1.89 3.85 7.65 16.99
restricted 0.21 0.58 0.85 1.41 2.85 3.45 12.6
0
4
8
12
16
20
Time(seconds)
Iterations
Figure 1: Evolution of average Improve [%]. Figure1.Evolutionof averageImprove[%]
Figure2.Evolutionof averagerunningtime
1 10 20 30 50 100 200
withoutrestricted 14.61 16.25 16.48 16.66 16.77 16.94 16.94
restricted 10.38 11.63 11.79 11.94 12.03 13.69 13.69
0
4
8
12
16
20
Improve[%]
Iterations
1 10 20 30 50 100 200
withoutrestricted 0.28 0.78 1.14 1.89 3.85 7.65 16.99
restricted 0.21 0.58 0.85 1.41 2.85 3.45 12.6
0
4
8
12
16
20
Time(seconds)
Iterations
Figure 2: Evolution of average running time.
for the Bounded-mTSP in [17], our algorithm reaches bet-
ter solutions for 12 out of 16 tested instances at a rea-
sonable computational time. In addition, our algorithm
can find the optimal solutions (eil51-2-23-27, eil51-3-15-
20) or near-optimal solutions (berlin52-2-10-41, berlin52-
3-10-27, berlin52-5-6-17) for the problems with 50 vertices

102 Informatica 45 (2021) 93–104 H.-B. Ban
 OPT
 LB, UB Best.Sol Best.Sol Best.Sol Best.Sol Best.Sol Best.Sol Aver.Sol Time 
 2.36 1.26 5 0.33 0.93  1.32 0.34
 0.51 5.52 4.04 2.00 2.05 14.55 9.30 4.38 2.52
Table 5: Comparisons with the state of the art metaheuristics for Bounded-mTSP without restricted neighborhood.
 Best.Sol T Best.Sol T Best.Sol T
 Table 6: Comparisons with the state of the art metaheuristics for mTRP.
in several seconds in Table 5. For the mTRP in [8, 18] in
Table 6, the quality of our solutions is much better than I.
O. Ezzine et al.’s algorithm (IOE) in [8] and every compa-
rable with S. Nucamendi-Guillen et al.’s algorithm (SNG)
in [18]. Moreover, our algorithm can find the optimal solu-
tions for the problems for the mTRP instance with up to 76
vertices in several seconds.
6.4 Discussions
Metaheuristic approach is a suitable approach to solve
the large sized-problem. The VNS [14] is a popular
schemes used widely to solve NP-hard problems. They
are very effective for some variants of the mTRP problem
[15, 19, 20]. The proposed algorithm is encouraged by the
efficiency of the algorithms in [15, 19, 20]. However, the
difference between the proposed VNS and their VNS is that
our algorithm builds up penalty value during a search. Our
main contribution is to adapt the VNS scheme, that extends
the well known VNS by including constraint penalization,
to solve the Bounded-mTRP effectively.
A good metaheuristic must balance between exploration
and exploitation. Exploration is to create diverse solu-
tions on a global space, while exploitation is to focus on
the search good current local regions. In the proposed al-
gorithm, the VNS implements exploitation while shaking
maintains exploration. In terms of experiments, the pro-
posed algorithm obtains better solutions than the adapted
algorithms in [17] in many cases. We also implement two
strategies that provide several choices: 1) The first choice
is to run the proposed algorithm with one iteration in “re-
stricted" strategy that obtains 10:38% solution quality on
average. The running time is very fast; 2) The second is to
run the proposed algorithm with one iteration in "without
restricted" strategy that obtains an average improvement of
14:61%. The second option trades off solution quality and
running time; 3) The last is to run the proposed algorithm
no more than 100 iterations. The average improvement is
16:94%. The option is the best in terms of solution quality.
Moreover, the proposed algorithm obtains comparable
or better solutions than the algorithms for the mTRP and
Bounded-mTSP in [8, 17, 18]. It shows that our algorithm
is still effective for various problems.

Penalty Variable Neighborhood Search for. . . Informatica 45 (2021) 93–104 103
7 Conclusions
In this paper, we propose the first metaheuristic algorithm,
which is mainly based on the VNS and GRASP’s prin-
ciples to solve the problem. Extensive numerical experi-
ments on benchmark instances show that, on average, our
algorithm leads to significant improvement. For small in-
stacnes, our algorithm obtains the optimal solutions for the
problem with 76 vertices at a reasonable a mount of time.
For larger instances, the proposed algorithm reaches better
solutions than the state-of-the-art algorithms in many cases.
References
[1] A. Archer, A. Levin, and D. Williamson, “A Faster,
Better Approximation Algorithm For The Minimum
Latency Problem", J. SIAM, V ol. 37, No. 1, 2007, pp.
1472-1498. https://doi.org/10.1137/07068151x.
[2] F. Afrati, S. Cosmadakis, C. Papadimitriou, G. Papa-
georgiou, and N. Papakostantinou, “The Complexity
Of The Travelling Repairmen Problem", J. Informa-
tique Theorique Et Applications, V ol. 20, pp.79–87.
https://doi.org/10.1051/ita/1986200100791
[3] M. Avci, M.G. Avci, 2017, “A GRASP with iter-
ated local search for the Traveling Repairman Problem
with Profits", J. Comput. Ind. Eng. 113, PP. 323–332.
https://doi.org/10.1016/j.cie.2017.09.032.
[4] Ha-Bang Ban, Duc-Nghia Nguyen, and Kien-
Nguyen, “An Effective Metaheuristic for Multiple
Traveling Repairman Problem with Distance Con-
straints", J. CAI, V ol. 38, 2019, pp. 1001-1034.
https://doi.org/10.31577/cai_2019_4_883.
[5] A. Blum, P. Chalasani, D. Coppersmith, W. Pulley-
blank, P. Raghavan, and M. Sudan, “The Minimum
Latency Problem", Proc. STOC, 1994, pp.163-171.
https://doi.org/10.1145/195058.195125
[6] M.E. Bruni, P. Beraldi, S. Khodaparasti, “A Fast
Heuristic for Routing in Post-Disaster Humani-
tarian Relief Logistics", J. Computers and Op-
erations Research, V ol. 30, pp. 304-313, 2018.
https://doi.org/10.1016/j.trpro.2018.09.033.
[7] N. Christofides, A. Mingozzi, p. Toth, “Exact Al-
gorithms for the Vehicle Routing Problem based
on Spanning Tree and Shortest Path Relaxations".
j. Math Program, V ol. 20, 1981, pp. 255–282.
https://doi.org/10.1007/bf01589353.
[8] I. O. Ezzine, and Sonda Elloumi, “Polynomial For-
mulation and Heuristic Based Approach for the k-
Travelling Repairmen Problem", Int. J. Mathematics In
Operational Research, V ol. 4, No. 5, 2012, pp. 503-
514. https://doi.org/10.1504/ijmor.2012.048928.
[9] T.A. Feo, and M.G.C. Resende, “Greedy Randomized
Adaptive Search Procedures", J. Global Opt., 1995, pp.
109–133.
[10] F. Jittat, C. Harrelson, and S. Rao, “The k-
Traveling Repairmen Problem", Proc. ACM-SIAM,
2003, pp.655-664.
[11] D. S. Johnson, and L. A. Mcgeoch, “The Travel-
ing Salesman Problem: A Case Study In Local Opti-
mization In Local Search In Combinatorial Optimiza-
tion", E. Aarts and J. K. Lenstra, Eds., pp. 215-310.
https://doi.org/10.2307/j.ctv346t9c.13
[12] R. Jothi, and B. Raghavachari, “Minimum La-
tency Tours and The k-Traveling Repairmen
Problem", Proc. LATIN, 2004, pp. 423–433.
https://doi.org/10.1007/978-3-540-24698-5_46:
[13] O. Martin, S. W. Otto, and E. W. Felten, “Large-Step
Markov Chains For The Traveling Salesman Problem",
J. Complex Systems, V ol. 5, No. 3, 1991, pp. 299-326.
[14] N. Mladenovic, and P. Hansen, “Variable
Neighborhood Search", J. Operations Re-
search, V ol.24, No. 11 24, 1997, pp.1097-1100.
https://doi.org/10.1016/s0305-0548(97)00031-2.
[15] N. Mladenovic, D. Urosevi, and S. Hanafi, “Vari-
able Neighborhood Search for the Travelling De-
liveryman Problem", J. 4OR, 2012; 11: 1-17.
https://doi.org/10.1007/s10288-012-0212-1.
[16] Z. Luo, H. Qin, and A. Lim, “Branch-and-
Price-and-Cut for the Multiple Traveling Repair-
man Problem with Distance Constraints", J. Opera-
tions Research, V ol. 234, No. 1, 2013, pp. 49-60.
https://doi.org/10.1016/j.ejor.2013.09.014.
[17] R. Necula, M. Breaban, M. Raschip, “Perfor-
mance Evaluation of Ant Colony Systems for the
Single-Depot Multiple Traveling Salesman Prob-
lem", Proc. HAIS, vol. 9121, pp. 257-268, 2015.
https://doi.org/10.1007/978-3-319-19644-2_22.
[18] S. Nucamendi-Guillén, I. Martínez-Salazar, F. Angel-
Bello, and J. M. Moreno-Vega, “A Mixed Inte-
ger Formulation and An Efficient Metaheuristic Pro-
cedure for the k-Travelling Repairmen Problem",
J. JORS, V ol. 67, No. 8, 2016, pp. 1121-1134.
https://doi.org/10.1057/jors.2015.113.
[19] A. Salehipour, K. Sorensen, P. Goos, and O.Braysy,
“Efficient GRASP+VND and GRASP+VNS meta-
heuristics for the Traveling Repairman Problem", J.
Operations Research, 2011, V ol. 9, No. 2, 189-209.
https://doi.org/10.1007/s10288-011-0153-0.
[20] M. Silva, A. Subramanian, T. Vidal, and L.
Ochi, “A simple and effective metaheuristic for

104 Informatica 45 (2021) 93–104 H.-B. Ban
the Minimum Latency Problem", J. Operations Re-
search, V ol. 221, No. 3, 2012, pp.513-520. DOI:
10.1016/j.ejor.2012.03.044
[21] A. Piwonska, F. Seredynski, “A Genetic Al-
gorithm with a Penalty Function in the Selec-
tive Travelling Salesman Problem on a Road
Network. Proc. IPDPS, 2011. pp. 381-387.
https://doi.org/10.1109/ipdps.2011.177.