Available online at www.sciencedirect.com
Electronic Notes in Discrete Mathematics 69 (2018) 141–148 www.elsevier.com/locate/endm
An integer linear programming model for the constrained shortest path tour problem Rafael Castro de Andrade 1 Department of Statistics and Applied Mathematics Federal University of Cear´ a, Campus do Pici, BL 910 CEP 60.440-900, Fortaleza, Brazil
Rommel Dias Saraiva 2 Department of Computer Science Federal University of Cear´ a, Campus do Pici, BL 910 CEP 60.455-760, Fortaleza, Brazil
Abstract Let D = (V, A) be a directed graph with set of vertices V and set of arcs A, and let each arc (i, j) ∈ A, with i, j ∈ V , be associated with a non-negative cost. The constrained shortest path tour problem (CSPTP) is NP-Hard and consists in finding a shortest path between two distinct vertices s ∈ V and t ∈ V such that the path does not include repeated arcs and must visit a sequence of vertex disjoint subsets T1 , . . . , TN in this order. In this work, we formulate the CSPTP as an integer linear programming (ILP) model and present valid inequalities for the problem. Computational experiments performed on benchmark data sets from the literature show that our ILP model consistently outperforms existing exact algorithms for the CSPTP and finds optimal solutions for most instances. Keywords: Combinatorial optimization, Integer linear programming, Network flow problems, Constrained shortest path tour problem.
https://doi.org/10.1016/j.endm.2018.07.019 1571-0653/© 2018 Elsevier B.V. All rights reserved.
142
1
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
Introduction
Let D = (V, A) be a directed graph, where V is the set of vertices and A is the set of arcs. Let each arc (i, j) ∈ A, with i, j ∈ V , be associated with a nonnegative cost cij . The constrained shortest path tour problem (CSPTP) [3] consists in finding a shortest path between two distinct vertices s ∈ V and t ∈ V (or (s, t)-path, for short) such that at least one vertex of each disjoint subset T1 , . . . , TN is visited according to this order. At least one vertex belonging to Th must be visited after the path crosses at least one vertex belonging to Tl , with l = 1, . . . , N − 1, and h = l + 1, . . . , N . Any vertex can be visited more than once. However, each arc (i, j) ∈ A can be visited at most once. This (s, t)-path can be viewed as a sequence of arc disjoint subpaths, where s is the origin of the subpath with destination u1 ∈ T1 ; u1 is the origin of the subpath with destination u2 ∈ T2 ; and so on, until uN ∈ TN be the origin of the subpath with destination t. Practical applications of the problem arises in: (i) safe cargo delivery, where orders to N regions must be delivered by a vehicle that stores items following a priority queue; (ii) control of robot motions, where a robot performs ordered operations (from a set of N ) to manufacture workpieces; (iii) reliable network design [7], when one selects ways to travel arcs whose usage does not compromise the survivability of the flow that must be sent from a given origin to a given destination such that ordered sets of vital nodes must be crossed at least once [3], among others. The CSPTP can be viewed as a variant of some constrained shortest path problems. Indeed, the shortest path tour problem (SPTP) [4,5], which follows a similar definition of the CSPTP except for allowing repeated arcs; and the shortest path visiting specified vertices (SPVSV) [1], where the path must visit only once all vertices of a given set in any order and without repeated arcs. The CSPTP was introduced in [3], where the proposed exact branch-andbound (B&B) method solves instances of small size when taking advantage of the specific structure of the problem. A greedy randomized adaptive search procedure (GRASP) is also developed to address medium to large instances. The reader is also referred to [3] in order to assess theoretical properties of the CSPTP, including its membership to the NP-Hard class. Hamiltonian path problem (HPP) [2] is polynomially Karp-reducible to the CSPTP. Thus, it is rather unlikely to find a polynomial-time algorithm that solves the problem. Figure 1 shows an example of a digraph D and its arc costs. Consider s = 1 and t = 4 as source and destination vertices, respectively. Clearly, arcs 1 2
Email:
[email protected] Email:
[email protected]
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
143
1 1
1
2
1
3
1
4
5
Fig. 1. Digraph reproduced from [3].
(1, 2), (2, 3), and (3, 4) compose the classical shortest path from s to t with cost equal to 3. Now, let us define T1 = {3} and T2 = {2} as disjoint sets to be visited in this order, i.e., T1 first, then T2 . An optimal SPTP solution from s to t is composed by arcs (1, 2), (2, 3), (3, 2), (2, 3), and (3, 4), with cost equal to 5, where arc (2, 3) is used twice. Contrary to the SPTP, a CSPTP solution does not allow repeated arcs. Then, a corresponding CSPTP optimal solution includes arcs (1, 3), (3, 2), (2, 3), and (3, 4), of cost equal to 8. In this work, we formulate the CSPTP as an integer linear programming (ILP) model. We also present some valid inequalities that showed to be very helpful in solving the problem. To do so, we first add N new vertices and N i=1 2Ti new arcs to the original digraph. We show how this modification allows us to tackle the CSPTP by means of mathematical programming techniques. Computational experiments performed on benchmark instances from the literature show the quality of the proposed ILP model when compared to results obtained by existing state-of-the-art exact algorithms for the problem. The remainder of the text is organized as follows. In Section 2, we conceive an ILP model for the CSPTP and present valid inequalities for the problem. In Section 3, we perform computational experiments on benchmark instances from the literature and present comparisons with existing exact algorithms. In Section 4, we give a brief conclusion and future perspectives.
2
Integer linear programming model
To the best of our knowledge, there is no mathematical programming model for the CSPTP in the literature so far. In this section, we formulate the CSPTP as an ILP model. Initially, we adapt each input instance as follows. For every Tk , k = 1, . . . , N , we create a dummy vertex vk and add it to V . For all u ∈ Tk , we add arcs (u, vk ) and (vk , u) to A, each of them with associated cost equal to zero, as shown in Fig. 2. We refer to V (resp. A ) as the augmented set of vertices (resp. arcs) obtained after the instance transformation. Now, define a set P consisting of N + 1 ordered pairs of vertices as P = {(s, v1 ), (v1 , v2 ), (v2 , v3 ), . . . , (vN , t)}. Let any k-th element of P , 1 ≤ k ≤
144
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
v1
v2
0 0 0 0 s
t1a
T1
vN
0
0 (. . . )
T2
t1b
0
0
TN
t
Fig. 2. Adapting a general input instance.
N + 1, be represented by (ak , bk ). We refer to ak (resp. bk ) as its first (resp. second) vertex in (ak , bk ). For instance, for the second element of P above, a2 = v1 and b2 = v2 . The idea is to find N + 1 subpaths, each one having origin ak and destination bk related to an element (ak , bk ) of P , and use them to compute the constrained shortest (s, t)-path tour. More precisely, we have to construct a path between s0 and v1 , which naturally visits one vertex of T1 in the transformed instance; in a second phase, we determine a path between v1 and v2 , which visits one vertex of T2 ; and so on, until determining the last path from vN to t, which has no vertex to visit. To formulate the CSPTP as an ILP model, we define the following binary decision variables: 1, if arc (i, j) belongs to the (ak , bk )-subpath, xkij = 0, otherwise. for all (i, j) ∈ A , and for k = 1, . . . , N + 1. We first define the objective function (1) of the problem, which minimizes the cost of the (s, t)-path tour. Min
N +1
cij xkij
(1)
k=1 (i,j)∈A
Constraints (2) are arc flow conservation constraints needed for each subpath k = 1, . . . , N + 1. ⎧ ⎪ ⎨1, if i = ak , k k (2) xij − xji = −1, if i = bk , ∀k = 1, . . . , N + 1, ∀i ∈ V ⎪ ⎩ (i,j)∈A (j,i)∈A 0, otherwise. Constraints (3) are forcing constraints guaranteeing that the last arc belonging to subpath k defines the first arc belonging to subpath k + 1, for each subpath k = 1, . . . , N . xkivk ≤ xk+1 vk i ,
k = 1, . . . , N, ∀(i, vk ) ∈ A
(3)
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
145
Constraints (4) impose that each arc can be used at most once. N +1
xkij ≤ 1,
∀(i, j) ∈ A
(4)
k=1
Constraints (5) state that x variables referred to arcs (i, vk ), with i ∈ Tk , are null in any subpath p = k. xkivp = 0,
∀(i, vp ) ∈ A , k = 1, . . . , N + 1 : k = p
(5)
Finally, constraints (6) define the domain of the decision variables. xkij ∈ {0, 1},
∀(i, j) ∈ A , k = 1, . . . , N + 1
(6)
The ILP model (1)–(6), which we refer to as (Q), contains O(N × A ) binary decision variables, and O(N ×V +A ) constraints. Model (Q) can be enhanced with the valid inequalities (7)–(9). x1is = 0 (7) (i,s)∈A
+1 xN =0 ti
(8)
(t,i)∈A
xkij ≤ 1,
k = 1, . . . , N + 1, ∀j ∈ V
(9)
(i,j)∈A
Constraint (7) (resp. (8)) states that no arc enters (resp. leaves) the source vertex s (resp. the destination vertex t) in the first (resp. last) subpath, while constraints (9) impose that, for each subpath k = 1, . . . , N + 1, the number of arcs entering vertex j ∈ V is at most one. In what follows, we refer to (Q)+ as model (Q) with the addition of inequalities (7)–(9). Proposition 2.1 If constraints (5) are relaxed, disconnected subpaths in the original digraph can appear in a feasible solution. Proof. We resort to the example in Fig. 3 to prove it.
3
2
Computational experiments
We implemented model (Q)+ using the C++ API of CPLEX 12.7. Experiments were performed on a PC Intel Pentium i7, 8 × 3.60 GHz, 16 GB
146
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148 v 0
1
1
2
0 0
0 0
0
M
3
M
4
1
5
(a) Digraph with s = 1, t = 5, T1 = {2, 3, 4}, v as dummy vertex associated with T1 , and M as a big non-negative value. Values crossing arcs represent their costs. Let P = {(1, v), (v, 5)}, and let variable x1 (resp. x2 ) refer to the subpath from vertex 1 to v (resp. from v to 5). v 0 0 1
1
M
2
3
M
4
1
5
(b) An optimal solution represented by variables x112 , x12v , x2v2 , x223 , x234 , and x245 , all equal to 1, with cost equal to 2M + 2. v 0
1
1
2
0
0
M
3
0 4
1
5
(c) When relaxing constraints (5), the optimal solution is represented by x112 , x12v , x2v2 , x223 , x23v , x2v4 , and x245 , all equal to 1, with cost equal to M + 2. It is infeasible, since the dummy vertex v is visited more than once in the second subpath. Consequently, the path tour from 1 to 5 is disconnected in the original digraph. Fig. 3. Possible consequence when relaxing constraints (5).
DDR3 RAM under Linux Ubuntu 14.05 LTS. The CPU time limit for each instance was set to 3600 seconds. The instances adopted in our experiments are from [3], which includes complete and grid digraphs. They are pseudorandomly generated by using the algorithm in [6]. See [4] for more details. We conducted comparative analysis between model (Q)+ and the B&B method of [3] implemented using two different strategies to build the branching tree: depth first (BBdf) and best bound first (BBbf). Results are detailed in Table 1. For each set of instances, we report the kind of digraph under consideration (type), the number of vertices (|V |), and the number of arcs (|A|) For the number of sets to visit (|T |), it varies according to each set of considered instances. Thus, we inform it as an interval containing the minimum and the maximum number of sets to be visited for instances of same dimensions. For each exact solution approach, we present the number (opt) of instances
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
Instance type
|V |
complete complete
BBbf
(Q)+
BBdf
|A|
|T |
100
9900
[20, 80]
16/16
3.55
16/16
3.98
16/16
3.52
0.00
0.18
150
22350
[20, 75]
20/20
27.69
20/20
33.37
20/20
17.01
0.00
0.75
complete
200
39800
[100, 180]
8/20
n/a
8/20
n/a
20/20
73.08
0.00
3.52
grid
9×9
288
15
4/4
0.03
4/4
0.37
4/4
0.15
0.00
0.01
grid
5 × 20
350
19
2/4
n/a∗
2/4
n/a
4/4
0.31
0.03
0.04
grid
10 × 10
360
19
3/4
n/a
3/4
n/a
4/4
0.28
0.04
0.03
grid
7 × 15
376
19
2/4
n/a
2/4
n/a
4/4
0.30
0.02
0.04
grid
12 × 12
528
27
0/4
n/a
0/4
n/a
4/4
0.92
0.07
0.13
grid
14 × 14
728
37
0/4
n/a
0/4
n/a
4/4
6.26
0.38
1.07
grid
20 × 20
1520
76
0/4
n/a
0/4
n/a
1/4
n/a
n/a
102.35
∗
opt
cpu(s)
opt
147
cpu(s)
opt
cpu(s)
gap(%)
t(s)
n/a means not available due to either time or memory limitation.
Table 1 Summary of results for complete and grid digraphs of all sets of instances.
of that type solved to optimality, and the average execution time (cpu) in seconds to obtain their optimal solutions. For model (Q)+ , in particular, we also give the average integrality (gap) in percentage between the linear relaxed solutions and the optimal ones, and the average time (t) in seconds needed to solve their linear relaxations. Concerning complete digraphs, the proposed model found optimal solutions for all instances, while both BBbf and BBdf failed to solve 12 (out of 56) instances. To solve instances with 100 vertices, (Q)+ , BBbf, and BBdf needed, respectively, 3.52, 3.55, and 3.98 seconds (in average). For instances with 150 vertices, they took, respectively, 17.01, 27.69, and 33.37 seconds (in average). The average time needed by (Q)+ to solve instances with 200 vertices was 73.08 seconds, and the average time spent to solve linear relaxations is between 0.18 and 3.52 seconds for all sets of instances. The integrality gap produced by (Q)+ for these instances was equal to zero. Regarding grid digraphs, model (Q)+ found optimal solutions for almost the whole set of grid digraphs, except for three instances. For this data set, both BBbf and BBdf were not able to solve 17 (out of 28) instances, mainly those with number of vertices greater than 100. Actually, the obtained results confirm that the B&B method can only be used to solve small instances [3]. Our enhanced ILP model took, in average, between 0.15 and 6.26 seconds to solve instances with up to 196 vertices. The average time to solve the linear relaxation of instances with 400 vertices was 102.35 seconds.
148
4
R.C. de Andrade, R.D. Saraiva / Electronic Notes in Discrete Mathematics 69 (2018) 141–148
Final remarks
We addressed in this work the CSPTP. To the best of our knowledge, there was no mathematical programming model for the CSPTP in the literature, what encouraged us to formulate the problem as an ILP model and propose valid inequalities. Experiments performed on benchmark instances from the literature showed that our enhanced ILP model outperforms existing exact approaches for the CSPTP while providing optimal solutions for most cases. Future directions include the development of a Lagrangian relaxation for the problem as well as decomposition methods and additional valid inequalities.
Acknowledgement The authors sincerely thank Daniele Ferone for making available both the executable files of BBdf and BBbf solution approaches, and the benchmark instances reported in [3]. The authors are also grateful to the anonymous referees for their valuable comments and suggestions.
References [1] Andrade, R. C., New formulations for the elementary shortest-path problem visiting a given set of nodes, European Journal of Operational Research 254 (2016), pp. 755–768. [2] Bertossi, A. A., The edge hamiltonian path problem is NP-complete, Information Processing Letters 13 (1981), pp. 157–159. [3] Ferone, D., P. Festa, F. Guerriero and D. Lagan`a, The constrained shortest path tour problem, Computers & Operations Research 74 (2016), pp. 64–77. [4] Festa, P., Complexity analysis and optimization of the shortest path tour problem, Optimization Letters 6 (2012), pp. 163–175. [5] Festa, P., F. Guerriero, D. Lagan`a and R. Musmanno, Solving the shortest path tour problem, European Journal of Operational Research 230 (2013), pp. 464– 474. [6] Festa, P. and S. Pallottino, A pseudo-random networks generator, Technical report, Department of Mathematics and Applications “Renato Caccioppoli”, University of Napoli Federico II, Italy (2003). [7] Heegaard, P. E. and K. S. Trivedi, Network survivability modeling, Computer Networks 53 (2009), pp. 1215–1234.