EUROPEAN JOURNAL OF OPERATIONAL RESEARCH
ELSEVIER
European Joumal of Operational Research 104 (1998) 139-153
Theory and Methodology
Minimum directed 1-subtree relaxation for score orienteering problem Seiji Kataoka a, *, Takeo Yamada a, Susumu Morito b a Department of Computer Science, The National Defense Academy, Yokosuka, Kanagawa 239, Japan b Department of Industrial Engineering and Management, Waseda Unit,ersity, Ohkubo 3-4-1, Shinjuku-ku. Tokyo 169, Japan
Received 21 March 1995; accepted 15 August 1996
Abstract
Score orienteering is a variant of the orienteering game where competitors start at a specified control point, visit as many control points as possible and thereby collect prizes, and return to the starting point within a prescribed amount of time. The competitor with the highest amount of prizes wins the game. In this paper, we propose the minimum directed l-subtree problem as a new relaxation to this problem and develop two methods to improve its lower bound. We describe a cut and dual simplex method and a Lagrangian relaxation method, and construct an algorithm that combines these two methods in an appropriate way. Computational experiments are carried out to analyse the behavior of the proposed relaxation with respect to the characteristics of the test problems. Especially for large-scale instances, we find that the proposed relaxation is superior to the ordinary assignment relaxation. © 1998 Elsevier Science B.V. Keywords: Orienteering problem; Minimum directed l-subtree relaxation; Traveling salesman problem
1. I n t r o d u c t i o n
Score orienteering is a variant of the orienteering game usually played in a mountainous terrain, where we are given a set of control points with each having an associated score. Competitors are required to start from a specified control point, which is referred to as the center, visit as many control points as possible, and return to the center within a prescribed amount o f time. The score gained by a competitor is the sum of scores of the control points visited. Competitors returning to the center after the time limit are disqualified, and the competitor with the highest score wins the game. Due to the time limit, competitors usually cannot visit all the control points and hence need to select the points to visit. At the same time, they must determine the sequence of their visits. The optimization problem underlying this sport game is called the score orienteering problem (SOP) [12-14,20,25,27,28,31 ] or maximum collection problem [17,19]. The prize collecting traveling salesman
* E-mail:
[email protected]
0377-2217/98/$19.00 © 1998 Elsevier Science B.V. All rights reserved. PH S0377-2217(96)00309-8
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
140
problem [2,3,6,26], the traveling salesman's subtour problem [10], the generalized traveling salesman problem [11], and the selective travelling salesman problem [22] can be thought as variants of SOP. They are proved to be JF.ga-hard [9]. In this paper, we propose the minimum directed l-subtree problem (MD1SP) as a new relaxation to SOP and develop two methods to improve its lower bound. In the first method we start from a trivial solution to a linearly relaxed MD1SP, and try to improve the lower bound by introducing, one by one, a cut and by performing the dual simplex procedure. Three types of cuts are proposed for this purpose. Specifically, we show that one of the cuts gives a facet of the directed 1-subtree polytope. The second method is an application of the Lagrangian relaxation to our problem. We propose a strategy that combines these two methods in an appropriate way. The organization of this paper is as follows. Section 2 presents an integer programming formulation of SOP and introduce MD1SP as its relaxation. In Section 3, based on this relaxation, two methods are presented to improve the lower bound of this problem, and we combine these two to propose an overall strategy. Computational experiments are carried out in Section 4.
2. Problem formulation and its relaxation Let G = (N, A) be a directed graph, where N = {1, 2 . . . . . n} is a set of nodes (control points) and A g N X N is a set of arcs (possible links between control points). Associated with each arc (i, j ) e A is tij(>1 0) representing the units of time to traverse it. Also, each node i ~ N has the associated prize Pi(t> 0). Node 1 represents the center, and we assume p~ = 0. Our problem is to find a subtour (a circuit in G) starting from and returning to the center such that the total traverse time along the subtour is within the prescribed amount of time T, and the total prize collected en route is maximized. To formulate SOP, for an arbitrary subtour in G, we define 0-1 variables with the following meaning:
xiJ= Yi=
1 0,
if arc (i, j ) is in the subtour, otherwise. if node i is not in the subtour, otherwise.
1 0,
It is convenient to regard y~ as corresponding to the self-loop at node i. In the sequel, we assume that G has a self-loop at each node. Denoting P = Ei~ NPi, SOP is formulated as follows:
(soP) Max
E Pi( 1 - Y i ) = P i~N
s.t.
~-, PlY,
(2.1)
iEN
Y~.xij+yj= 1 V j E N ,
(2.2)
i~N
~-, x , i + Y i = l
Vi~N,
(2.3)
j~N
Y'. xij~ I V I - 1
VVgN\{1},
Vv~l,
(2.4)
i,j~ V
Y~. ti]x q <~T,
(2.5)
(i, j ) E A
xij, yiE{0, 1} V ( i , j ) ~ A , i~N.
(2.6)
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
141
Here the objective may be switched to minimizing Y'.,~ NPiYi. Constraints (2.2) and (2.3) are the so-called assignment constraints. Constraint (2.4) eliminates subtours which do not involve the center, and (2.5) is for time restriction. In the traveling salesman problem (TSP), the assignment problem and the minimum 1-tree problem are usually used for the relaxations [15,16,18,21,24]. Both of these can be tailored for SOP. However, the assignment relaxation tends to produce subtours that are far from the center, that are sometimes hard to eliminate in the subsequent branch-and-bound procedure. This is especially the case for SOP defined on a plane where t,j represents the Euclidean distance between nodes i and j. On the other hand, if we apply the l-tree relaxation and obtain a tree consisting of only one connected component rooted at the center, the above-mentioned difficulty is unlikely to occur and we can expect a good lower bound from this relaxation. We therefore try to extend the 1-tree relaxation to SOP and thereby define the minimum directed 1-subtree relaxation. Directed 1-subtree is a directed subtree rooted at the center plus an additional arc back to the center, and thus has only one cycle in it [29]. Before formulating the minimum directed 1-subtree problem, we note that adding the constraints
Xijq-yi<~ 1 V(i, j) ~ a
(2.7)
to SOP does not change the problem at all, since (2.7) are valid for any feasible solutions to SOP. Therefore, in what follows by SOP we mean minimizing Ei ~ u Pi Yi subject to (2.2)-(2.7). With the Lagrangian multipliers w~ for (2.3) and ~b >/0 for (2.5), we obtain the following Lagrangian relaxation problem R(~r, + ) for SOP: R ( ~ , ~b) Min
Y'.
(d?t~j+lri)xij+ Y'~ (pi+~r~)y,-[dpT+ ~ / v ~ ]
(i, j ) ~ A
s.t.
i~N
(2.8)
i
Y'. x,i + yj = 1 Vj ~ N,
(2.9)
i~N
xO+yi<~ 1 V(i, j) E A ,
E %~< IVl-1
(2.10)
VV_.N\{1},
V4:~,
(2.11)
i,jEV
x,j, yi~{O, 1} V ( i , j ) ~ a ,
i~N.
(2.12)
For each node i constraint (2.10) forbids i to have a self-loop and an outgoing arc (i, j ) simultaneously. Because of these constraints, it is clear that a feasible solution to R ( ~ , qb) corresponds to a directed l-subtree. Thus for a fixed (~r, ~b), R ( ~ , ~b) is also referred to as the Minimum Directed l-Subtree Problem (MD1SP). MD1SP belongs to one of AP~-hard problems [9], even if minimum directed spanning tree problem, that is also called arborescence, is solved easily [4,5,8,23,30].
3. Two methods to improve lower bounds
For simplicity, we change the notation for the objective function of R ( ~ , ~b) as
E (i, j ) E A
c,)xij + E s i Y i ,
(3.1)
iEN
where c o = qbtij + ~, and s~ = p~ + "tr,. In this section, we propose two methods to obtain stronger lower bound to MD1SP as well as to SOP. Specifically in Section 3.2, we explicitly discuss strategies to improve lower bound to SOP.
142
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
© 1
a cycle
SA pairs
Fig. 1. A cycle and SA pairs.
3.1. Cut and dual simplex method Eliminating (2.10) and (2.11) and replacing (2.12) with 0 ~< xi:, Yi ~< 1 defines a relaxation of MD1SP. The resulting problem can be solved easily as follows: for each j E N, we find the smaller between ci,j = min i ~ N{Co} and s t, and set the corresponding variable, either x~,j or yj, to 1. All the other variables are set to 0. If this solution is also feasible to M D I S P , this clearly solves MD1SP. Otherwise, we look for a valid inequality, often referred to as a cut, add this to the relaxed problem, and solve the enlarged problem using the dual simplex method. Repeating this procedure we either obtain an optimal solution to MD1SP, or even if the process is truncated before an optimal solution is obtained, the final objective value gives a stronger lower bound to MD1SP. We present three kinds of cuts for this purpose.
3.1.1. Cut 1 The set of non-zero variables in the solution to the relaxed problem defines a subgraph in G. Such a subgraph may involve a cycle since (2.11) has been dropped. Also, since (2.10) has been relaxed, it may have a self-loop and an outgoing arc simultaneously at a same node. The latter is said to be an SA pair. This is shown in Fig. 1. If a cycle that does not include the center is found, we add
x,/< I V I - 1
(3.2)
i,j~ V
to the relaxed problem, where V is the set of nodes that constitute the cycle. Similarly, if an SA pair is found we adopt the following cut: Xij "Jr-Yi ~ l .
(3.3)
Inequalities (3.2) and (3.3) are referred to as Cut 1. We can also apply Cut 1 when the linearly relaxed solution has fractional values. In such a case, we prepare a directed graph consisting of the arcs with xi~ > 0, and take a strongly connected component in this graph as V in (3.2).
3.1.2. Cut 2 If the linearly relaxed solution happens to have all 0 - 1 values after applying Cut 1, a basic variable with value 1 can be written as 1 minus a linear combination of nonbasic variables x N. Let B denote the set of basic variables that forms a cycle or an SA pair. Then we obtain the following:
~_~ xij+ ~., yj= I n l - ~ ' x s . xij~B
(3.4)
)'j~B
Note that some components of tj may be zero. Let R be the set of nonbasic variables that correspond to the positive elements of ~ and ~k denote the slack variable associated with Cut 1. Then we obtain the inequality
E xij~R
xij+ E Yj+ E yj~R
~ER
~k >~ 1.
(3.5)
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
143
This is Cut 2. Note that when R is the subset of all nonbasic variables, Cut 2 is stronger than the Dantzig cut [21]. Cut 1 can be regarded as a special case of Cut 2 since the latter coincides with the former when every element of ~ is either 0 or 1. 3.1.3. Cut 3 Let V be a set of nodes other than the center, i.e., V c N \ {1}, V 4: 0. From (2.9) we obtain Ei ~ N]~i ~ v xij + Ej~ v Yj = I V 1. Restricting the first summation to V, we have ~ , ~ v E j ~ vXii + E j ~ v Yj <~ ] V I. Note that for a feasible solution to MD1SP, we have Ej ~ v Yj ~< I V I and Ei j ~ v x q ~< ] V [ - 1. Then for a set S c V, we obtain
Y'. E x , j + iEV j ~ V
Ey~
(3.6)
jES
which is termed as Cut 3. To detect a solution violating Cut 3, we consider a graph consisting arcs and self-loops with x,j > 0 and y~ > 0. Then we first diminish the graph by recursively deleting nodes each of which has only one arc. In the remaining graph, we can take the set V in (3.6) with some minute rules. The subset S is set to V \ { i } , where y i ( ~ V) is the smallest value. Cut 3 has some nice properties, but before describing these we need some notations. For a combinatorial problem Q described in terms of an integer programming problem, the convex hull of the feasible region of Q is denoted as conv(Q), while its dimension is written as dim(Q). The proofs to the following propositions are given in Appendix A. Proposition 3,1. dim(MD1SP) = n(n - 1). Proposition 3.2. When I SI = I V I - 1. Cut 3 gives a f a c e t o f conv(MDlSP). 3.2. Lagrangian relaxation method
The Lagrangian relaxation method [7,32] of this section is based on the bounding procedures initially developed by Balas and Christofides [1] for TSP. Here we take the same approach to solve M D I S P as well as SOP. 3.2.1. Description o f the method In the following integer programming problem:
Min s.t.
cx Aix = b I, A2x = b 2, x : integer vector,
we assume A l X = b 1 consists of a relatively small number of constraints. Also A 1 is assumed to have some special structure, such as totally unimodular, which is likely to yield integer solutions. As an example of such equations, we mention the assignment constraints. We call these constraints easy. On the other hand, A 2 x = b 2 are assumed to be very large in size, and hard to solve. In the case of TSP, the subtour elimination constraints may be regarded of this type. These are referred to as hard constraints. In the Lagrangian relaxation method, we first ignore the hard constraints and solve the remaining as a linear programming problem. The optimal objective value obtained herein is clearly a lower bound to the original problem, but usually the solution violates some of the hard constraints. Let a t x ~- b be a violated constraint. At
144
S. Kataoka et a l . / European Journal of Operational Research 104 (1998) 139-153
Fig. 2. Initial admissible graph G'.
this stage, the reduced cost vector is ~" = c - A t l ul(t> 0), where u I is the Lagrangian multiplier vector. Then we add the violated constraint to the relaxed problem with a Lagrangian multiplier u = min {~Jaj}.
(3.7)
aj> 0
This changes the objective value to btlul + bu, while u I remains unchanged and therefore without destroying dual feasibility. The reduced cost vector g'new is now ~+.ew --- g" - ua.
(3.8)
Repeating this procedure, we obtain improved lower bounds.
3.2.2. Application to SOP In applying the Lagrangian relaxation method to M D 1 S P or SOP, we take (2.9) as the easy constraints, and regard all others as hard ones. An admissible graph consists of the arcs and self-loops that correspond to the variables with zero reduced costs. This is denoted as G' = (N, A'). By the procedure in Section 3.1 we have an initial admissible graph as shown in Fig. 2. To this graph we add some arcs or self-loops, one by one, to construct a feasible solution to M D 1 S P as well as to SOP, observing the following four conditions. C1. There exist at least one arc originating from the center. C2. F o r each i ~ V\{1}, there exists either a self-loop at i or a path from the center to i. C3. For each i ~ V\{1}, there exists either a self-loop at i or a path from i to the center. C4. The graph either does not contain an articulation point, or otherwise each node beyond such a point from the center has a self-loop. Note that conditions C I and C2 are necessary for a solution to be feasible to M D I S P . In addition C3 and C4 are required for SOP. Here we prepare some notations. In an admissible graph when there is a path from node i to j , we say that i is connected to j, and denote this as i --*j. When i ~ j and i ~ j , these nodes are strongly connected, which is denoted as i ~ j. Note that each node is strongly connected to itself, i.e., i ~ i. A set of nodes S is said to be strongly connected, if any two nodes in S are strongly connected. A maximal strongly connected set is referred to as a strongly connected component or component for short. Generally, G' can be divided into several components. To indicate relations between the components, we use the same notations as for nodes; e.g., S~ ~ S 2, etc. A component S is called inactive if each node in the component has a self-loop; i.e., Y"i ~ S Yi = [ S [. Otherwise the component S is called active; i.e., ~.i e s Yi < [ S [. Root component R is an active component that is not connected to any other active components; if there exists such a component S as S ~ R, then S must be inactive. Head component H is an active component which does not connect to any other active component; if there exists a component S that H ~ S, then S must be inactive.
3.2.3. Four procedures Corresponding to conditions C 1 - C 4 , this section presents four procedures to augment the admissible graph G' in terms of a series of examples. W e start after solving relaxation problem with easy constraints (2.9). The
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
U
145
" Fig. 3. G' after Procedure 1.
solution has an incoming arc or a self-loop at each node. Therefore, in an initial reduced cost matrix, there are at least one 0 in each column. The arcs in the admissible graph coincide with 0s in the following reduced cost matrix. Example 0. In the reduced cost matrix shown below, the brackets indicate the candidates of arcs or self-loops which might be added to the admissible graph. The element in bold typeface represents the arc actually selected in the procedures: (18)
( 12)
(3)
(11)
(18)
(17)'
0
0
0
10
10
5
9
7
8
5
0
13
9
8
1
5
2
6
0
5
12
4 6
14 11
6 12
8 5
2 12
12 0
15 0m
11
12
6
9
15
14
3
The corresponding G' is shown in Fig. 2. Procedure 1. In the initial G' if no arc originates from node 1, then due to C1 we revive the following constraint that has been removed in the relaxation ]~ xlj >/1.
(3.9)
j~Kt
Here K t = N \ { I } , and the arc (I, Jl) is added to G', where ?U~ = minje rf?lj}. The reduced costs are updated according to (3.8).
Example 1. The heavy arc (1, 4) has been added to G' in Fig. 3. Correspondingly, the reduced cost matrix is updated to
(15) 0 7
1 4 6 11
(8) (5) (14) (11) (12)
(9)
Q
8
(15)
(14)
Q
10
10
(5)
(9)
(5) (2) (6) (12) (6)
0 6 8 5 9
13 Q 2 12 15
(9) (5) (12) Q (14)
(8) (12) (15) Q (3)
146
s. Kataoka et al./ European Journal of Operational Research 104 (1998) 139-153
P r o c e d u r e 2. Next, we look for a root component R which does not include the center and define the set of nodes K 2 = R t3 {il i ~ R}. Note that, by definition of the root component, each node in {i [ i ~ R} belongs to an inactive component. Note also that 1 ~ K 2. Let K 2 = N \ K 2 and r be a node in R without a self-loop. Then C2 is represented as xij + y,>~ l ,
(3.10)
r~R.
i~K2,jEK 2
Either an arc (i 2, J2) or a self-loop at r 2, whichever of the smaller reduced cost between 3i,~ = mini~F,~.j~x:{?ij} and ~ = min~R{~r}, is added to G'. The reduced costs are updated using (3.8). After Procedure 2, either one of the following occurs: (i) the set of nodes in R are reachable from the center, (ii) R becomes inactive, or (iii) new root components are created. This procedure is repeated until there is no root components remaining in G'. E x a m p l e 2. In Fig. 3, {3} and {7} are root components. For each of these, K 2 is {2, 3} and {6, 7} respectively. First, we take {2, 3} as K 2 and pick up arc (4, 3). Then remaining root component is only {7}, and adding self-loop at node 7 we obtain G' of Fig. 4. Correspondingly the matrix is updated to 13
7
0
8
12
11
0
0
0
10
10
2
6
(7)
(6)
(3)
Q
(13)
(6)
(5)
(1)
(3)
0
(6)
Q
(2)
(9)
(4)
(12)
(4)
(8)
(2)
(9)
(12)
6
9
lO
5
12
0
0
ll
lO
4
9
15
11
0
P r o c e d u r e 3. In a feasible solution to SOP, if a node h does not have a self-loop, we are able to find a path from h back to the center. To satisfy C3, we take a head component H in G' and define K 3 = H t3 {i [ i *-- H}. Note that nodes in {i J i *-- H} belong to inactive components. Also note that g 3 -,'->|. Let K3 = N \ K 3 and h ~ H be a node that does not have a self-loop. Then we add the following inequality: xij+yh>~l,
(3.11)
h~H,
i~K3,j~K, 3
to the set of constraints. In a similar way as in Procedure 2, either an arc (i 3, J3) or a self-loop at h 3, whichever of the smaller reduced cost between ?i3j3 =mini~g,.jEg~{?ij} and ~h, =minh~/4{Sh}, is added to G'. The reduced costs are updated by (3.8). Resultant admissible graph G' has either a path from the node h to the
Fig. 4. G' after Procedure 2.
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
147
Fig. 5. G' after Procedure 3.
center, a self-loop at h, or new head components. This procedure is repeated until there is no head components remaining in G'. E x a m p l e 3. In Fig. 4 the head component is only {5}, with K 3 = {5}, too. We pick up the self-loop at node 5. After applying Procedure 3, a new head component {3, 4} is created and we have K 3 = {3, 4, 5}. We repeat Procedure 3, taking this time arc (4, 1) to add to G'. Then no head component is found. The resulting G' is shown in Fig. 5, with the reduced cost matrix (13)
(7)
Q
8
12
ll'
0
0
0
10
10
2
6
(6)
(5)
(2)
0
12
5
4
0
2
0
5
0
1
8
1
9
1
5
0
6
9
6
9
10
5
12
0
0
II
lO
4
9
15
11
0
Procedure 4. After applying the above three procedures, we can decompose the nodes of G' into the following four sets; K s = { i l l ~ i}, K R = {ill ~ i, 1 ~ i}, K u = {ill ~ i, 1 ~ i}, and K 0 = {ill -~ i, 1 o- i}. Note that each node in K R, K n, K o has a self-loop, and clearly K s includes the center. Hence the nodes in K s are candidate for an optimal subtour of SOP. For each node d ~ K s that does not have self-loop, there exist a pair of paths from the center to d and from d back to the center. If these paths meet in between, we do not obtain a simple cycle in K s . This will happen when K s has an articulation point. Suppose there exists an articulation point ap ~ K s. Then we can divide K s into 2-connected components. Define K o as the component which does not include the center and K 4 = K s \ { K D U {ap}}. Here d is a node in K o without a self-loop. We also distinguish nodes in K R depending on whether the path from node 1 to that node includes the ap or not. The same discussion is applied to K n. That is: KR1 = { i l l ~ i, 1 ~- i, ap is not included in the path 1 ~ i}, KRd = {i I d ~ i, d ~ i, ap is not included in the path d ~ i}, K n , = { i l l ~ i, 1 *-- i, a n is not included in the path 1 ~ i}, K nd = { i I d -~ i, d ~ i, a e is not included in the path d ~ i}.
Define also the set of arcs Acut as follows: a c , , = {( i, j ) l i ~ { g 4 U KRI } , j ~ { g o 0 grid}} t.3 {( i, j ) li ~ { g o U gRd } , j E { g 4 U KH1}}.
148
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
Fig. 6. G' after Procedure 4.
Then we obtain X i j ' q - y d ~ 1,
E
d ~ K o.
(3.12)
(i, j)~Acu t
Similarly, either an arc (i 4, j4 ) o r a self-loop a t d 4 , whichever of the smaller reduced cost between ?i,j4 = minfi.j)~ A,~{F.ij} and sa4 = minae Ko{sd}, is added to G'. The reduced costs are updated by (3.8). This procedure is repeated until there is no articulation point remaining in G'. Example 4. In Fig. 5, node 4 is an articulation point. Hence, K 4 = {1}, K o = {3}, and K 0 = {6, 7}. Then, KR1 = KRd = ¢, and K n j = KHd = {2}. We take the self-loop at node 3 and add it to G'. The resultant G' is shown in Fig. 6, with the corresponding reduced cost matrix 11
5
0
8
12
11
0
0
10
10
2
6
4
3
0
0
12
5
4
0
2
0
5
0
1
8
1
9
1
5
0
6
9
6
9
10
5
12
0
0
11
10
4
9
15
11
0
0
-
-
B
After these four procedures, we find in G' a feasible solution to SOP. 3.3. A combination o f two methods In the previous subsections, we have discussed two methods to improve lower bounds. We now combine these two methods. In applying the Lagrangian relaxation method, we first solve a linearly relaxed problem with easy constraints. However, there is no clear-cut distinction between easy and hard constraints. In Section 3.2, we regarded (2.9) as easy constraints, while Fischetti and Toth [6] took the assignment constraints (2.2) and (2.3) for this purpose. The strategy we take here is to adopt (2.9) together with some cuts discussed in Section 3.1 as easy constraints. The strategy consists of two phases. In Phase 1, we add cuts and perform dual simplex method as in Section 3.1. Upon completion of this procedure, we switch to Phase 2 and apply the Lagrangian relaxation method of Section 3.2. In Phase 1, we apply Cut 1, Cut 2, and Cut 3, and keep lastly obtained integer solution. In Phase 2, we use the integer solution inherited from Phase 1, since in the Lagrangian relaxation integrality of solution is preferred. Thus the timing of switching from Phase 1 to Phase 2 is of crucial importance. We propose the following three criteria for this purpose. That is, we switch from Phase 1 to Phase 2 when: (a) an optimal solution to MD1SP is obtained, (b) a valid solution to the proposed three cuts is obtained, or
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
149
(c) the prescribed number of cuts have been added. In the following computational experiments, we examine the performance and analyse the behavior of these methods in some detail.
4. Computational experiments Although SOP is a variant of TSP or the vehicle routing problem, there exist few bench mark such as TSPLIB for TSP. Therefore, in our experiments we use the following two types of instances. The first is the random instance, where each coefficient ci/ is generated by uniformly random integer in [1, 1000]. In TSPs with random coefficients, assignment relaxation is frequently used and succeeds to solve large scale instances. But in actual applications, random instances appears to be rather rare. The other type is the 'slope' instance of Kataoka and Morito [18], which is a slight modification of the instance with two dimensional Euclidean distance. Here with (xi, Yi) denoting the location of the ith point in the plane (a square of 1000 × 1000), c# is given by:
~( xi _ x))2 + ( Y, _ yj)Z
if x~ >1xdj,
Cp2(xi_xj)2+(yi_yj)2
if x i < x v
(4.1)
ciJ=
In the experiment we set p = 2.0, which appears to produce the most intractable instances for TSP [18]. The objective of our experiments is to examine the efficiency of proposed lower bounds. The assignment optimal Table 1 The random instance Size
Easy constraints
n
Assignment
MD1SP(0)
MDIST(0.5)
MD1SP(1)
MDISP(2)
MDISP(3)
MDISP(4)
I0
1.0000 1.0181 0.0020
0-0-50 0.6868 1.0433 0.0014
19-0-31 0.8633 0.9835 0.0304
27-0-23 0.9095 0.9915 0.0511
44-0-6 0.9450 1.0072 0.1277
49-0-1 0.9553 1.0084 0.1767
50-0-0 0.9553 1.0109 0.1970
20
1.0000 1.0305 0.0060
0-0-50 0.6784 1.0609 0.0076
14-0-36 0.8082 0.9397 0.4883
32-0-18 0.8367 0.9245 1.0195
46-0-4 0.8477 0.9215 1.6310
50-0-0 0.8508 0.9217 2.9108
50-0-0 0.8508 0.9217 2.9054
30
1.0000 1.0176 0.0258
0-0-50 0.6560 1.0444 0.2500
23-0-27 0.7562 0.9128 1.8873
34-0-16 0.7673 0.9015 3.4930
46-1-3 0.7756 0.8962 11.2431
49-1-0 0.7770 0.8966 19.1100
49-1-0 0.7770 0.8966 19.1202
40
1.0000 1.0157 0.0840
0-0-50 0.6446 1.0313 0.0080
20-1-29 0.7213 0.8756 4.7523
32-1-17 0.7274 0.8666 8.5159
44-4-2 0.7331 0.8572 57.5056
44-6-0 0.7335 0.8572 82.6192
44-6-0 0.7335 0.8572 82.5546
50
1.0000 1.0081 0.1164
0-0-50 0.6444 1.0158 0.1300
19-0-21 0.7086 0.8393 12.3800
31-1-18 0.7136 0.8248 47.4659
46-2-2 0.7167 0.8171 189.0671
46-2-2 0.7186 0.8171 189.1118
46-4-0 0.7185 0.8171 148.3779
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
150 Table
2
The Slope instance
Size
Easy
n
Assignment
MDISP(0)
MDIST(0.5)
MDISP(1)
MDISP(2)
MDISP(3)
MDISP(4)
10
1.0000 1.0092 0.0008
0-0-50 0.7323 1.0099 0.00 i 2
23-0-27 0.9182 1.0364 0.0262
30-0-20 0.9542 1.0447 0.0418
40-0-10 0.9830 1.0469 0.0851
45-0-5 0.9983 1.0535 0.1554
48-0-2 1.0064 1.0536 0.2310
20
1.0000 1.0371 0.0050
0-0-50 0.8488 1.0426 0.0040
24-0-26 0.9843 1.0310 0.3205
43-0-7 1.0078 1.0384 0.6718
46-0-4 1.0165 1.0378 0.8205
47-0-3 1.0208 1.0380 1.0501
49-0-1 1.0232 1.0384 2.2097
30
-
0-0-50 0.8940 1.0702 0.0134
26-0-24 1.0130 ! .0594 1.9972
40-0-10 1.0333 1.0600 3.0063
47-0-3 1.0457 1.0621 5.9036
48-0-2 1.0485 1.0621 6.8541
40-0-1 1.0500 1.0621 9.8505
1.0000 1.0773 0.0320
0-0-50 0.9104 1.0889 0.0344
18-0- 32 1.0251 1.0715 6.3937
35-0-15 1.0526 1.0723 10.3405
47-0-3 1.0657 1.0740 24.0233
49-0-1 1.0684 1.0745 32.7855
49-0-1 1.0695 1.0745 32.7333
1.0000 1.1038 0.0770
0-0-50 0.9063 1.1224 0.0854
5 -0-45 1.0269 i .0931 20.1812
28-0-22 1.0787 1.0985 41.2416
47-0-3 1.0990 1.1019 91.1177
50-0-0 1.1008 1.1021 114.1868
50-0-0 1.1008 1.1021 114.1178
constraints
1.0000 1.0599 0.0124 40
50
-
value is taken as a standard. In Section 4.1, we inquire into the behavior of instance types and the number of added cuts. Section 4.2 shows the results for large scale instances. 4.1. Instance types Table 1 and Table give a summary of the experiments carried out for the random and slope instances, respectively. The coefficient s i in (3.1) associated at each node is randomly distributed in [I, 100]. In these tables each column represents the set of inequalities taken as 'easy' in Phase 1. Therefore, 'assignment' stands for the assignment constraints which coincides with those used by Fischetti and Toth [6], while in M D 1 S P ( f ) , in addition to these at most f n cuts are considered as 'easy'. Each row is average of 50 independent runs. Four entries in each cell show: (i) The number of switching from Phase 1 to Phase 2, shown as a-b-c, where a, b, c show the number of cases (out of 50) in which a: MD1SP is solved to optimality; b: a valid solution to the proposed cuts is obtained; c: the number of cuts to be added reaches to f n . This is not shown for assignment cases. (ii) The ratio of the lower bound after Phase 1 over the optimal assignment value. (iii) The same ratio as (ii) after Phase 2. (iv) The CPU time in seconds on a microSPARC II (85 MHz) using SPARComplier C 3.0.1.
S. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
Table 3 The larger case Size n
Random instance Assignment
MD1SP(0)
Slope instance Assignment
151
MDISP(0)
100
1.0000 1.0031 0.2594
0.6345 1.0064 0.3452
1.0000 1.2093 1.8934
0.9018 1.2472 1.9680
200
1.0000 1.0017 0.8886
0.6739 1.0014 5.2608
1.0000 1.3898 130.5480
0.8679 1.4714 129.6056
500
1.0000 1.0011 4.7852
0.6738 1.0014 5.2744
1000
1.0000 1.0004 11.4852
0.7328 1.0013 30.1532
In random instances, it is believed that assignment relaxation produces good lower bounds quickly. In our experiments, the optimal value thus obtained was generally quite satisfactory. These findings coincide with the experimental results obtained for TSP. However, in MD1SP relaxation adding cuts increases CPU time excessively, and still the resultant lower bound do not compete with the assignment value. In addition, Phase I often ends up with obtaining valid solutions for the cuts, which is obviously undesirable. For slope instances, too much cuts are again time consuming, However, in our experiments no runs ended up with obtaining valid solution for the cuts in Phase 1. Thus, the proposed cuts are useful if our purpose is to obtain a minimum directed l-subtree. In this case, the assignment solution is usually unsatisfactory as a lower bound. In M D 1 S P ( / ) , applying Lagrangian relaxation in Phase 2 is therefore valuable. However, adding cuts and performing dual simplex method often decrease each value of reduced cost matrix. Thus, increasing lower bound too much in Phase 1 does not necessarily lead to the higher lower bound in Phase 2. We can conclude by saying that adding cuts and performing the dual simplex method is usually a good strategy to solve the directed 1-subtree problem, but not necessarily for SOP. Tables 1 and 2 show that in MD1SP(0) we usually obtain a better lower bound than in 'assignment' both for random and slope instances. In the next subsection, we check if MD1ST(0) is superior to 'assignment' even for larger sized problems. 4.2. L a r g e s c a l e i n s t a n c e s
Table 3 gives the average of 50 runs for each n from 100 to 1000. Each row shows (ii), (iii) and (iv) of 4.1. In random instance, the assignment relaxation is satisfactory even if we do not apply Phase 2. On the other hand, for slope instance, Phase 2 works well for increasing lower bound dramatically, but it consumes too much CPU time both in assignment and MS1SP(0) strategies. Also we observe that more than 5% higher lower bounds can be obtained by MD1SP(0) than assignment. We see that for large n the slope instances are much more intractable than the random, even if we are to calculate lower bound only.
5. Conclusions
We have presented the minimum directed 1-subtree problem as a new relaxation to SOP. To obtain a stronger lower bound to this problem, we proposed the method of adding cuts and performing dual simplex procedure in
152
s. Kataoka et al. / European Journal of Operational Research 104 (1998) 139-153
Phase 1, and applying the Lagrangian relaxation in Phase 2. This strategy has been examined in series of computational experiments, and for large scale instances of both random and slope type, we ascertained that the lower bound from our strategy usually outperforms those obtained from the ordinary assignment relaxation. To solve SOP in the Lagrangian framework, we need to maximize R ( ~ , ~b) with respect to xt and d~. This may be accomplished by improving xt and ~b using the subgradient method. W e need to specify the parameters, such as the direction and the step-size of each iteration. This is left for future investigations.
Appendix A Proof of Proposition 3.1. In the formulation of MD1SP, we have n independent equality constraints with n 2 variables. Thus, dim(MD1SP) .<< n 2 - n -- n(n - 1). Therefore the proposition is proved if we are able to find n ( n - 1) linearly independent points in MD1SP. F o r the case of Yi = 0 and xil = 0 for all i ~ N \ { 1 } , the problem can be regarded as the minimum arborescence problem, whose dimension is known to be n 2 - 3n + 2 [8]. A feasible solution to 1-arborescence is defined as an arborescence plus an incoming arc to the center. W e can see the points given by x n = 1 for all i E N \ { 1 } are linearly independent, and show that dim(1-arborescence) = n 2 - 3n + 2 + (n - 1) = (n - 1) 2. Besides, we have n - 1 linearly independent points given by ys = 1 and Yi = O, i ~ N \ { 1 , s} ( s = 2, 3 . . . . . n). Thus in total M D 1 S P has at least ( n - 1) 2 + ( n - 1) = n ( n - 1) linearly independent points. [] Proof of P r o p o s i t i o n 3.2. As stated in the above proof, the minimum l-arborescence problem is equivalent to M D 1 S P if Yi = 0 for all i ~ N \ { 1 } . Also the dimension of the 1-arborescence problem is (n - 1) 2 and that of M D 1 S P is n(n - 1). Then for a given V and S = V \ { k } , if we are able to find n - 1 more linearly independent points in such a way that they make up a directed 1-subtree and satisfy (3.6) with equality, the proposition is proved. Case 1. W h e n each node in V has a self-loop: Suppose ( x v, y P ) ~ conv(MD1SP). Clearly E s ~ sY p = I SI = [ V I - 1. Then ( x p, y V ) satisfies (3.6) with equality. Case 2. When node k does not have a self-loop, and one of the nodes in S has a self-loop: Suppose ( x q, y q ) ~ conv(MD1SP), yq = 0, and yq = 1 for an s ~ S. Also suppose there exist I V I - 2 arcs in V. Then the solution satisfies (3.6) with equality. Case 3. W h e n no node in V has a self-loop, and one of the nodes in N \ { V U {1}} has a self-loop: Suppose ( x r, y r ) ~ conv(MD1SP), and these exist I V I - 1 arcs in V. Then ( x ~, y r ) satisfies (3.6) with equality. Clearly the points in these three cases are linearly independent and the number of such points in each case is 1, I S I, and (n - 1) - ( I S I + 1) respectively. Therefore, in addition to those in the case of 1-arborescence, directed 1-subtree has n - 1 more linearly independent points. []
References [1] Balas, E., and Christofides, N., "A restricted Lagrangean approach to the traveling salesman problem", Mathematical Programming 21 (1981) 19-46. [2] Balas, E., and Matin, G., "ROLL-A-ROUND software package for scheduling the rounds of a rolling mill", Research Report No. MSRR-539, GSlA, Carnegie Mellon University, Pittsburgh, 1985. [3] Balas, E., "The prize collecting traveling salesman problem", Networks 19 (1989) 621-636. [4] Chu, Y., and Lin, T., "On the shortest arborescence of a directed graph", Scientia Sinica 4 (1965) 1396-1400. [5] Edmonds, J., "Optimum branchings", Journal of Research of the National Bureau of Standards - B. Mathematics and Mathematical Physics 71/4 (1967) 233-240. [6] Fischetti, M., and Toth, P., "An additive approach for the optimal solution to the prize collecting TSP", in: B.L. Golden and A.A. Assad (eds.), Vehicle Routing: Methods and Studies, North-Holland, Amsterdam, 1988, 319-343.
S. Kataoka et a l . / European Journal of Operational Research 104 (1998) 139-153
153
[7] Fischetti, M., and Toth, P., "An bounding procedure for combinatorial optimization problem", Operations Research 37/2 (1989) 319-328. [8] Fulkerson, D.R., "Packing rooted directed cuts in a weighted directed graph", Mathematical Programming 6 (1974) I-13. [9] Garey, M.R., and Johnson, D.S., Computers and Intractability: A Guide to the Theory of NP-Completeness, Freeman, San Francisco, CA, 1978. [10] Gensch, D.H., "An industrial application of the traveling salesman's subtour problem", AIIE Transaction 10/4 (1978) 362-370. [11] Golden, B.L., Levy, L., and Dahl, R., "Two generalizations of the traveling salesman problem", OMEGA 9 / 4 (1981) 439-445. [12] Golden, B.L., Levy, L., and Vohra, R., "The orienteering problem", Naval Research Logistics 34 (1987) 307-318. [13] Golden, B.L., Wang, Q., and Liu, L., "A multifaceted heuristic for the orienteering problem", Naval Research Logistics 35 (1988) 359-366. [14] Hayes, M., and Norman, J.M., "Dynamic programming in orienteering: Route choice and the siting of controls", Journal of the Operational Research Society 35/9 (1984) 791-796. [15] Held, M., and Karp, R.M., "The traveling salesman problem and minimum spanning trees", Operations Research 18 (1970) 1138-1182.
[16] Held, M., and Karp, R.M., "The traveling salesman problem and minimum spanning trees: Part 2", Mathematical Programming 1 (1971) 6-25. [ 17] Kataoka, S., and Morito, S., "Single constraint maximum collection problem", Journal of the Operations Research Society of Japan 31/4 (1988) 515-531. [18] Kataoka, S., and Morito, S., "'Selection of relaxation problem for asymmetric traveling salesman problem instances", Journal of the Operations Research Society of Japan 34/3 (1991) 233-249. [19] Kataoka, S., "Structure of problems and feature of instances for routing problems", Doctor's Thesis, Department of Industrial Engineering and Management, Waseda University, 1993 (in Japanese). [20] Keller, C.P., "Algorithm to solve the orienteering problem: A comparison", European Journal of Operational Research 42 (1989) 224-231. [21] Konno, H., and Suzuki, H., "Integer programming problem and combinatorial optimization", Nikka-giren Syuppan, 1982 (in Japanese). [22] Laporte, G., and Martello, S., "The selective travelling salesman problem", Discrete Applied Mathematics 26 (1990) 193-207. [23] Lawler, E.L., Combinatorial Optimization: Networks and Matroids, Holt, Rinehart Winston, New York, 1977. [24] Lawler, E.L., Lenstra, J.K., Rinnooy Kan, A.H.G., and Shmoys, D.B., The Traveling Salesman Problem: A Guided Tour of Combinatorial Optimization Problems, Wiley, New York, 1985. [25] Leifer, A.C., and Rosenwein, M.B., "Strong linear programming relaxations for the orienteering problem", European Journal of Operational Research 73 (1994) 517-523. [26] Pekny, J.F., Miller, D.L., and McRae, G.J., "An exact parallel algorithm for scheduling when production costs depend on consecutive system states' ', Computers and Chemical Engineering 14/9 (1990), 1009-1023. [27] Ramesh, R., and Brown, K.M., "An efficient foar-phase heuristic for the generalized orienteering problem' ', Computers & Operations Research 18/2 (1991) 151-165. [28] Ramesh, R., Yoon, Y., and Karwan, M.H., "An optimal algorithm for the orienteering tour problem", ORSA Journal on Computing 4/2 (1992) 155-165. [29] Smith, TH.C., "A LIFO implicit enumeration algorithm for the asymmetric travelling salesman problem using a one-arborescence relaxation", Mathematical Programming Study 12 (1980) 108-114. [30] Tarjan, R.E., "'Finding optimum branchings", Networks 7 (1977) 25-35. [31] Tsiligirides, T., "Heuristic methods applied to orienteering", Journal of the Operational Research Society 35/9 (1984) 797-809. [32] Wong, R.T., "A dual ascent approach for Steiner tree problems on a directed graph", Mathematical Programming 23 (1984) 271-287.