Computers & Operations Research 33 (2006) 3363 – 3383 www.elsevier.com/locate/cor
Lower and upper bounds for the mixed capacitated arc routing problem José-Manuel Belenguera , Enrique Benaventa , Philippe Lacommeb , Christian Prinsc,∗ a Dept d’Estadística i Investigació Operativa, Universitat de València, C / Dr. Moliner, 50-46100 Burjassot (València), Spain b LIMOS, Université Blaise Pascal, BP 10125, 63173 Aubière Cedex, France c ISTIT, Université de Technologie de Troyes, BP 2060, 10010 Troyes Cedex, France
Available online 30 March 2005
Abstract This paper presents a linear formulation, valid inequalities, and a lower bounding procedure for the mixed capacitated arc routing problem (MCARP). Moreover, three constructive heuristics and a memetic algorithm are described. Lower and upper bounds have been compared on two sets of randomly generated instances. Computational results show that the average gaps between lower and upper bounds are 0.51% and 0.33%, respectively. 䉷 2005 Elsevier Ltd. All rights reserved. Keywords: Capacitated arc routing problem; Mixed graph; Lower bound; Cutting plane; Heuristic; Memetic algorithm; Waste collection
1. Introduction and background The capacitated arc routing problem (CARP), introduced by Golden and Wong [1], considers an undirected network in which each edge has a traversal cost. A fleet of identical vehicles of limited capacity is based at a depot node. A subset of required edges with known demands must be serviced by a vehicle. The CARP consists of determining a set of vehicle trips of minimum total cost, such that each trip starts and ends at the depot, each required edge is serviced by one single trip and during one traversal, ∗ Corresponding author. Tel.: +33 325 71 56 41; fax: +33 325 71 56 49.
E-mail addresses:
[email protected] (J.-M. Belenguer),
[email protected] (E. Benavent),
[email protected] (P. Lacomme),
[email protected] (C. Prins). 0305-0548/$ - see front matter 䉷 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.cor.2005.02.009
3364
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
and the total demand for each trip does not exceed vehicle capacity. Assad and Golden [2] list various applications for the CARP, like municipal waste collection, winter gritting and inspection of power lines. The CARP is NP-hard. Exact methods [3] are still limited to 20 required edges, explaining the importance of heuristics. Path-scanning [4], augment-merge [1] and Ulusoy’s heuristic [5] are good examples of fast heuristics to compute one acceptable solution. Much better results can be obtained by metaheuristics. They include tabu search algorithms designed by Eglese and Li [6] and by Hertz et al. [7], memetic algorithms proposed by Lacomme et al. [8,9], Greistorfer’s tabu-scatter search approach [10] and a guided local search from Beullens et al. [11]. Comparisons with tight lower bounds computed by Belenguer and Benavent [12] show that a majority of instances with up to 100 required edges can be solved to optimality by these algorithms. In fact, the CARP is too simplistic to model real street networks encountered for instance in urban refuse collection. Indeed, its undirected graph can only model 2-way streets whose both sides are collected in parallel and in any direction. Even if this practice is common in low-traffic areas, actual networks also include two-way streets with independent sides and one-way streets. Clearly, such extensions can be tackled by extending the CARP with a mixed network, leading to a new problem called Mixed CARP or MCARP. The MCARP is defined on a mixed network G = (V , E ∪ A) with a set V of n nodes, an edge-set E and an arc-set A. A fleet of K identical vehicles with limited capacity Q is based at a depot node s. K is either fixed or left as a decision variable. The m elements of E ∪ A are called links. Each time a link u is traversed without service by a vehicle, a cost c(u) called deadheading cost is incurred. A subset of required links must be serviced by the vehicles. It comprises a subset of required edges ER ⊆ E and a subset of required arcs AR ⊆ A. Each required link u has a demand r(u) and a service cost w(u). All costs and demands are non-negative integers. The goal is to determine a least-cost set of vehicle trips, such that each trip starts and ends at the depot, each required link is serviced by one single trip and during one traversal (no split service), and the total demand processed by a trip fits vehicle capacity. The MCARP is NP-hard, since it includes NP-hard particular cases: the CARP or UCARP (undirected graph), the Directed CARP or DCARP (directed graph), the mixed rural postman problem or MRPP (if Q = ∞) and the mixed Chinese postman problem or MCPP (if Q=∞ and all links are required). To the best of our knowledge, no published paper addresses the MCARP. The paper is organized as follows. Section 2 deals with a linear programming approach leading to a cutting plane algorithm that computes tight lower bounds for the MCARP. Fast constructive heuristics are presented in Section 3. Section 4 introduces a memetic algorithm. Section 5 is devoted to a computational evaluation, in which lower and upper bounds are compared on two sets of randomly generated MCARP instances. 2. A lower bound for the MCARP 2.1. A linear formulation This section presents a cutting plane algorithm to obtain a lower bound for the MCARP. The method is based on the formulations proposed by Belenguer and Benavent [12] for the UCARP, by Nobert and Picard [13] for the MCPP and by Corberán et al. [14] for the general routing problem on a mixed graph. Nobert and Picard introduced the so-called balance constraints, that allow formulating the MCPP by
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
3365
considering only one variable per edge and arc. This makes the formulation more compact and, probably, tighter. They used the following characterization of an eulerian mixed graph, due to Ford and Fulkerson [15]: • Evenness condition. The degree of each node (number of links incident with it) is even. • Balanced sets condition. For every proper subset of nodes S, the number of arcs entering S minus the number of arcs leaving S is less than or equal to the number of edges between S and V \S. Let us introduce some notation. Given a subset of nodes S ⊆ V \{s}, S = ∅, we denote by + (S) the set of arcs leaving S, by − (S) the set of arcs entering S, and by (S) the set of edges between S and V \S. The cutset induced by S is defined as L (S) = + (S) ∪ − (S) ∪ (S). Let L(S) be the set of links having both extreme nodes in S. A subscript R in any of these sets indicates the corresponding subset of required edges/arcs/links: for instance + R (S) denotes the set of required arcs leaving S. r(S) denotes the total demand of the required links of L R (S) ∪ LR (S) and k(S) = r(S)/Q. Note that at least k(S) vehicles are needed to service the demands of the required links in L R (S) ∪ LR (S). The unbalance of the + (S)| − | (S)| − | (S)|. cutset L (S) is defined as b(S) = |− R R R Given an MCARP solution (a set of trips for the vehicles), let xu be the total number of times the link u ∈ E ∪ A is traversed without service by the vehicles. We may construct an eulerian supergraph that corresponds to this MCARP solution, by adding to graph GR =(V , ER ∪AR ) xu copies (called deadheading links) of every link u ∈ E ∪A. For any subsets E ⊆ E, and A ⊆ A, we denote x(E ∪A )= u∈E ∪A xu . This leads to the following relaxed formulation for the MCARP: Min
c(u)xu ,
(2.1)
x(+ (S)) − x(− (S)) + x((S)) b(S) ∀S ⊆ V , S = ∅, x(L (S)) 2k(S) − |L R (S)| ∀S ⊆ V \{s}, S = ∅, L x( (S)) 1 ∀S ⊆ V \{s}, S = ∅, |L R (S)| is odd, xu 0 ∀u ∈ E ∪ A, xu integer ∀u ∈ E ∪ A.
(2.2) (2.3) (2.4) (2.5) (2.6)
u∈ER ∪AR
s.t.:
w(u) +
u∈E∪A
Constraints (2.2) or balance constraints express the balanced sets condition to be satisfied by an eulerian mixed graph. Constraints (2.3) are called capacity constraints and are explained as follows: as noted above, since at least k(S) vehicles are needed to service the total demand of the required links L R (S)∪LR (S), the L total number of deadheading traversals of the links in the cutset (S) must be at least 2k(S) − |L R (S)|. L Obviously, these constraints are only interesting when 2k(S) − |R (S)| > 0. Constraints (2.4) are called odd cutset constraints and directly derive from the evenness condition. Note that they can be dominated, in some cases, by the capacity constraints. A directed version of the capacity constraints is also possible. Note that any vehicle servicing a link L in L R (S) ∪ LR (S) has to cross the cutset (S) twice, once from V \S to S and another time from S to V \S. Thus, the total number of times that the cutset has to be crossed from V \S to S is at least k(S), so we have the following valid constraints: x(− (S)) + x((S)) k(S) − |− R (S)| − |R (S)|
∀S ⊆ V \{s}, S = ∅.
(2.7)
3366
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
Nevertheless, these constraints, as well as the corresponding ones in the reverse direction, have not been used because, as it can be easily shown, they can be generated from constraints (2.2) and (2.3). The above formulation is not complete for the MCARP because it may contain integer solutions that do not correspond to any MCARP solution. Furthermore, no complete formulation using the sole variables xu is known for the MCARP. This can be explained by the fact that variables xu provide information on the deadheading traversals of the set of routes in the MCARP solution, but nothing about the individual trips (see Belenguer and Benavent [12] for a more detailed discussion on this topic for the UCARP). Nevertheless, the formulation can be improved by adding new sets of valid constraints. Namely, the disjoint path (DP) inequalities, that were introduced by Belenguer and Benavent [12] for the UCARP, can be easily adapted for the MCARP. In what follows, we show that any valid inequality in the UCARP formulation of Belenguer and Benavent [12] can be translated into a valid inequality for the MCARP. The UCARP formulation uses the same variables xu , but they are associated with edges only, obviously. Consider the MCARP defined on the mixed graph G = (V , A ∪ E), and let GU = (V , E ∪ U ) be the undirected graph obtained from G by substituting every arc a ∈ A by an edge with the same cost and demand. Let us denote by f (u) ∈ A the arc that has originated the edge u ∈ U . Note that any feasible solution for the MCARP is also a feasible solution for the UCARP defined on GU (as an edge can be E∪U traversed in any direction). Now, let b ∈ R and a0 ∈ R be such that the inequality u∈E∪U bu xu a0 is a valid inequality for the UCARP defined on GU . Then, the inequality u∈E bu xu + u∈U bu xf (u) a0 is obviously valid for the MCARP defined on G. 2.2. A cutting plane algorithm for the MCARP We propose a cutting plane algorithm (CPA) to obtain a lower bound for the MCARP which is based on the above formulation. This procedure can briefly be described as follows. At each iteration, we solve a linear program (LP) that contains: the objective function (2.1), the non-negativity constraints (2.5) and a subset of valid constraints (2.2)–(2.4) and DP inequalities. Note that the optimal cost of this LP is a lower bound for the MCARP. Then, we look for a set of valid inequalities violated by the optimal LP solution, we add them to the LP and proceed as before. The algorithm stops when no violated inequality is found. Note that even if the LP solution is integer, it may not correspond to an optimal solution for the MCARP because our formulation is not complete. The set of constraints for the first LP is built as follows. Starting at a node v ∈ V \{s}, we generate a sequence of node sets (Wi ) such that W0 = {v} and, for i 1, Wi+1 is generated from Wi by adding to it all the nodes that are incident in G with at least one node of Wi (excluding the depot, if necessary). The procedure stops when the current node set contains all the nodes except the depot. For each node set Wi , constraints (2.2) corresponding to Wi and V \Wi are added to the LP if their right-hand sides are positive; furthermore, constraints (2.3) corresponding to Wi and V \(Wi ∪ {s}) are added to the LP if their right-hand sides are positive, and otherwise, the constraint that is added is (2.4) in the case that |L R (Wi )| is odd. At each iteration of the cutting plane algorithm several identification procedures have been used to generate violated inequalities. Let x¯ ∈ RE∪A be the optimal LP solution at the current iteration of the cutting plane algorithm. Violated balance constraints can be identified as follows. Let ze = x¯e + 1 for each required link e ∈ ER ∪ AR , ze = x¯e for a non-required link, and consider the weighted graph G (z) = (V , E ∪ A) with weights ze . It is easy to see that the balance constraint (2.2), corresponding to set S ⊆ V , is violated by x¯
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
3367
iff f (S) = z(+ (S)) − z(− (S)) + z((S)) < 0. A set S for which f (S) < 0 is called an unbalanced set. Nobert and Picard [13] showed that the most unbalanced set in G (z) can be found by solving a maximum flow problem on a graph with |V | + 2 nodes (see also Benavent et al. [16] for a short description of the maximum flow problem to be solved). Violated constraints (2.3), (2.4) and DP inequalities can be identified on a graph where the arcs are considered as (non-directed) edges. The identification procedures described by Belenguer and Benavent in [12] are briefly recalled below. Let G
(x) ¯ be the undirected weighted graph induced by the links e ∈ E ∪ A such that x¯e > 0. The identification of violated odd cutset constraints (2.4) can be done in polynomial time using the odd minimum cut algorithm of Padberg and Rao [17]. Unfortunately, no exact polynomial algorithm is known to identify violated capacity constraints (2.3). We have used two heuristics that generate several node sets S ⊆ V \{s} to see whether (2.3) is violated. Furthermore, for every node set that is generated, constraints (2.2) and (2.4) are also tested. The first heuristic consists of finding the connected components of G
(x). ¯ The second heuristic consists of increasing the edge demands by a given percentage and finding a minimum cut in an auxiliary graph. The procedure is inspired by an algorithm [12] that identifies a relaxed version of the capacity constraints, called fractional capacity constraints, and is applied with 10 different values of increasing percentage. The identification of violated DP constraints is much more complex than the previous ones. One such constraint may involve one or several cutsets, so we use the strategy of storing in a pool the node sets generated by the above heuristics and try to combine them to build valid DP inequalities that are violated by the LP solution. The preceding identification algorithms were combined using the following strategy: the heuristics that generate node sets S ⊆ V \{s} are called first; in case that these heuristics do not find any violated constraint of types (2.2), (2.3) or (2.4), the exact identification procedures for (2.2) and (2.4) are called. Finally, when all the previous procedures fail, the identification of violated DP inequalities is invoked. 3. Constructive heuristics for the MCARP 3.1. Introduction In this section, three CARP heuristics cited in introduction are generalized for the MCARP: pathscanning [4], augment-merge [1] and Ulusoy’s heuristic [5]. Several possible improvements are also suggested, before their validation in the computational evaluation of Section 5. In fact, the heuristics are extended for an enriched MCARP with five additional features: • • • • •
a mixed multigraph, allowing parallel links between two nodes, prohibited turns and turn penalties, a common trip cost limit L (also called range constraint), windy edges, whose deadheading and/or service costs depend on the direction, several dumping sites.
These additional features play a key role in applications. For instance, in solid waste collection, a service road parallel to a main road is easily tackled with a multigraph. Prohibited turns are widespread and must be taken into account. Some undesirable turns like U-turns can be allowed, but penalized with an
3368
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
extra-cost (turn penalty). A cost limit on the trips is useful to handle a vehicle range or a maximum driving time. A sloping street can be handled as a windy edge. Finally, several dumping sites are in general available to unload a vehicle before its return to the depot. 3.2. Additional notations for heuristics For all heuristics of Sections 3 and 4, the mixed multigraph G is coded as a fully directed multigraph
G∗ = (V , A∗ ), by replacing each edge by two opposite arcs and by adding one dummy loop for the depot
and for each dumping site (the role of these loops is explained below). A∗ is a set of m∗ arcs, identified by indexes from 1 to m∗ instead of pairs of nodes, to distinguish between parallel arcs. Each arc u in A∗ has a head h(u), a tail t (u) and a deadheading cost c(u). The = + required links of G correspond in G∗ to a subset R ⊆ A∗ with = 2 · + arcs. Each arc u ∈ R has a demand r(u), a service cost w(u) and a pointer inv(u). Each arc of AR is coded in R as one arc u with inv(u) = 0. Each edge from ER gives in R two opposite arcs u and v, with inv(u) = v, inv(v) = u, r(u) = r(v), possibly with distinct costs if the edge is windy. G∗ can include also a set of permitted turns (u, v), with associated turn penalties pen(u, v). Following Benavent and Soler [18] for the RPP with turn penalties, we define a feasible path from arc u to arc v as a sequence of arcs = (u = u1 , u2 , . . . , uk = v), such that, for i = 1 to k − 1, h(ui ) = t (ui+1 ) and the turn (ui , ui+1 ) is allowed. The cost of is defined by Eq. (3.1). Using an adaptation of Dijkstra’s shortest path algorithm [19], forbidden turns are made transparent by pre-computing an m∗ × m∗ arc-to-arc distance matrix D, in which d(u, v) is the cost of shortest feasible paths from u to v. Note that c(u) and c(v) are not counted in d(u, v), to ease operations like inserting one required arc z between two others u and v in a trip: the cost variation simply is d(u, z) + w(z) + d(z, v) − d(u, v) c() = pen(u, u2 ) +
k−1
(c(ui ) + pen(ui , ui+1 )).
(3.1)
i=2
Let P be the set of allowed dumping nodes. Distances between one arc and the depot s or one dumping node are simply tackled by including in A∗ one fictitious loop = (s, s) for the depot and one similar loop for each dumping node. Note that |A∗ | = m∗ = 2 · |E| + |A| + |P | + 1. The best dumping loop dump(u) for a trip whose last required arc is u and the minimal cost term(u) to finish the trip after u can be pre-computed by Eqs. (3.2–3.3), in,in which is the set of dumping loops and the cost of one dump operation. dump(u) = arg min{d(u, v) + + d(v, ) : v ∈ },
(3.2)
term(u) = d(u, dump(u)) + + d(dump(u), ).
(3.3)
Solutions and trips are coded by lists. In the sequel, for any list , || denotes the cardinal and (k) or k the element at rank k. A trip is a list of arcs from R representing distinct required links. Implicit shortest deadheading paths are assumed between these arcs. The load of the trip, load( ), is the sum of its demands and must not exceed Q (Eq. (3.4)). The cost of the trip, cost( ), includes the service costs of its required links, the costs of its connecting paths, and the cheapest cost for completing the trip. It must not exceed L (Eq. (3.5)). An MCARP solution S simply is a list of trips (S1 , S2 , . . . , S|S| ). Its cost
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
6
3369
d
u 8 3
3
c
a
6
b
9
depot
Fig. 1. Path-scanning iteration (with the first criterion).
is the sum of its trip costs. r( i ) Q, load( ) =
(3.4)
i=1,| |
cost( ) = d(, 1 ) +
(w( i ) + d( i , i+1 )) + w( | | ) + term( | | ) L.
(3.5)
i=1,| |−1
3.3. Path-scanning heuristic Path-scanning (PS) is a sequential heuristic building one trip at a time, inspired by the nearest-neighbour heuristic for the TSP. The emerging trip is extended recursively by adding at the end the closest required link. However, in arc routing, the distance from the last link serviced u to the nearest link v not yet serviced is often 0, e.g., when u and v are adjacent streets. This is why PS breaks ties by minimizing one of the five following criteria: • • • • •
F1 (v) = term(v), the cost for completing the trip after v; F2 (v) = −term(v), equivalent to maximizing F1 (v); F3 (v) = r(v)/w(v), a kind of “productivity” for link v; F4 (v) = −r(v)/w(v), equivalent to maximizing F3 (v); F5 (v) = if load( ) W/2 then F5 (v) := F2 (v) else F5 (v) := F1 (v) endif
In F5 , the trip completion cost is maximized as long as the vehicle is at most half-full. Afterward, it is minimized. Fig. 1 gives an example, containing required edges only, with one dumping site coinciding with the depot and rule F1 . The required edges are in boldface. Thin lines indicate shortest deadheading paths. The current trip ends at required edge u. The minimal distance to the edges not yet serviced is 3, with a set of two candidates Z = {a, b}. Edge a is chosen to extend the trip, since it minimizes the return cost to the depot. For one given criterion F ∈ {F1 , F2 , F3 , F4 , F5 }, a solution can be computed in O(2 ), i.e., in O(n2 ) for real street networks with m∗ ≈ 4n. The PS heuristic consists of constructing five solutions
3370
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383 Trip S(i ) : link u u v Trip S( j): link v
S(i ) absorbs S( j )
Concatenation of S(i ) and S( j ) u v 10 s(i ) s( j) 8
6
Fig. 2. Augment (left) and Merge (right) operations in augment-merge.
(for F = F1 .F5 ) and returning the best one. By randomizing the computation of v, it is possible to construct more than five solutions and to improve solution costs on average. We suggest two variants PSRC (random criterion) and PSRL (random link) that are compared with PS in Section 5.4. PSRC randomly draws F in {F1 , F2 , F3 , F4 , F5 } before computing v. PSRL ignores the criteria, it randomly draws v among the closest links (set Z), instead of choosing the link minimizing F on Z. 3.4. Augment-merge Fig. 2 illustrates augment-merge (AM) for the UCARP. AM starts with a trivial solution S with one trip per required link. These trips are sorted in decreasing cost order. For each trip S(i), i = 1, 2, . . . , − 1, a first phase Augment scans each smaller trip S(j ), j = i + 1, i + 2, . . . , . If the unique required link v of S(j ) is on a deadheading path of S(i) and load(S(i)) + r(v) Q, then trip S(j ) is absorbed. The cost of S(i) does not change because deadheading and service costs are equal in the UCARP, but the total cost decreases by cost(S(j )). A second phase Merge evaluates the concatenation of any two trips {S(i), S(j )}, subject to Q: e.g., concatenating S(i) then S(j ) in Fig. 2 yields a saving 8 + 6 − 10 = 4. Merge concatenates the two trips with the largest positive saving and repeats this process until no such concatenation is possible. There are 8 ways of concatenating two trips: S(i) then S(j ) or S(j ) then S(i), with each trip inverted or not. In fact, only 4 ways need to be tested thanks to symmetries: S(i) then S(j ) gives the same saving as the inverse of S(j ) then the inverse of S(i). In the MCARP, cost variations are affected by the multiple dumping sites, the distinction between service and deadheading costs, and the fact that a trip is no longer equivalent to its inverse, due to the asymmetry of D. To illustrate these changes, assume that several dumping sites are available in Fig. 2. The best initial “daisy” is obtained by selecting the best direction for each required edge u, i.e., the trip for u costs min(d(, u) + w(u) + term(u), d(, inv(u)) + w(inv(u)) + term(inv(u))). The saving of the Augment operation is now cost(S(j )) + term(u) − d(u, v) − w(v) − term(v). But, since v is on the path from u to the closest dumping site, term(u) = d(u, v) + c(v) + term(v) and the saving can be simplified as cost(S(j )) + c(v) − w(v). It is still non-negative, since w(v) is contained in cost(S(j )). However, in addition to the capacity constraint, the absorption must now satisfy the range constraint cost(S(i)) + w(v) − c(v) L. In Merge, the 8 ways of merging two trips must now be evaluated. Let last( ) and rev( ) denote, respectively, the last required link and the inverse of a trip . Note that all required edges of are flipped in rev( ), e.g., if = (u, v), then rev( ) = (inv(v), inv(u)). The merger in Fig. 2 yields a saving of term(u) + d(, v) − d(u, v), while merging S(i) and rev(S(j )) would save cost(S(j )) + term(u) − cost(rev(S(j )))+d(, inv(last(S(j ))))−d(u, inv(last(S(j )))). In practice, hopefully, many mergers can be discarded because the trips containing required arcs are not invertible. We do not detail the algorithm of
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383 c(5)
20
b(3) 5 10
11
d(1) 15
20
5
5
5 9
13
10 7 e(6)
12 20
depot
(a)
10
20
5
15
11
14
11 a(4)
5 3 Trip
20
5
13
9
Trip 1
5
20
3371
7 16
16 (c) bcd(80)
Node index: Node label: (b)
bc(56) a(37) 0 1 b(27) 2 c(40) 0 37 51 ab(51)
de(50) 3 d(32) 4 e(33) 5 91 115 141 cd(64)
Fig. 3. Principle of split: (a) initial giant tour S = (a, b, c, d, e), (b) auxiliary graph and shortest path (labels in each node), (c) resulting trips.
AM, which is rather tedious because of the eight kinds of mergers. It can be implemented in a non-trivial way in O(2 . log ), i.e., O(n2 . log n) for real street networks with m∗ ≈ 4n. Two improvements are proposed. Firstly, the augment phase is disputable. If all links are required, it works well because the links absorbed by a trip are often contiguous. If not, the internal deadheads created between absorbed links cannot be recovered during the merge phase, leading to degraded results. This suggests a version of AM without augment, called M. Secondly, several pairs of trips may give the same saving in the merge phase. Some testing with M has shown an improvement if ties are broken by merging the pair {S(i), S(j )} with the largest discrepancy in load, i.e., such that |load(S(i)−load(S(j ))| is maximal. This version is called IM (Improved M). M and IM can be preceded or not by an Augment phase, giving four versions AM, M, AIM and IM which are evaluated in Section 5.5. 3.5. Ulusoy’s heuristic The principle of the original UCARP heuristic is to compute a giant tour covering the required links (by relaxing vehicle capacity) and to split it into capacity feasible trips, using a shortest path procedure Split described below. Since an optimal giant tour does not necessarily provide an optimal solution after splitting, the first step may be solved by any heuristic. Fig. 3 depicts Split on a simple undirected example for a giant tour = (a, b, c, d, e) of five required links, with demands in brackets and deadheading costs. Assuming here Q = 9, Split builds an acyclic auxiliary graph H with + 1 nodes indexed from 0 onward. Each sub-list ( (i), . . . , (j )) corresponding to a feasible trip is modeled in H as one arc (i − 1, j ), weighted by the trip cost. The sub-list is feasible if Q and L are not violated. The load and cost of the corresponding trip are computed using Eqs. (3.4) and (3.5). A shortest path from node 0 to node in H (boldface) indicates the best splitting, here 3 trips and a total cost 141. Note that Split returns an optimal MCARP solution, subject to the sequence fixed by the giant trip. The shortest path is computed with Bellman’s algorithm for directed acyclic graphs [19], whose complexity is linear in the number of arcs. For H, this number is at most · ( − 1)/2. This O(2 ) becomes
3372
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
O(n2 ) for real street networks with m∗ ≈ 4n. It can also be improved when the smallest demand rmin is large enough: since a trip contains at most = Q/rmin required links, H contains O( · ) arcs. Our implementation runs in O() space, without generating H explicitly: for i = 1, 2, . . . , , all feasible sub-lists ( (i), . . . , (j )) are scanned and evaluated, and their costs are used to directly update node labels in Bellman’s algorithm. Ulusoy’s heuristic can be improved by splitting several giant tours and by returning the best MCARP solution obtained. We suggest to relax vehicle capacity Q and trip cost limit L and to call path-scanning or its variants PSRC and PSRL to get giant tours. Recall that PS builds five giant tours (one per criterion) while PSRC and PSRL generate a specified number of tours (Section 3.3). The resulting heuristics are called UPS, URC and URL, respectively. A last kind of improvement can be obtained by playing on the auxiliary graph H. For instance, in Fig. 3, arc (1, 4) in H corresponds to trip (b, c, d) with a cost 80. Assume a distance of 5 from link d to link b. Starting the trip by c and scanning the same sublist of links circularly leads to a better trip (c, d, b) with cost 71. Better results can be achieved by assigning to each arc of H the cost corresponding to the best rotation. This technique can be applied to UPS, URC and URL, giving improved versions respectively called IUPS, IURC and IURL. All these variants are compared in Section 5.6. 4. Memetic algorithm 4.1. The existing memetic algorithm A memetic algorithm (MA) is a genetic algorithm hybridised with a local search, see Moscato [20]. Lacomme et al. [8] proposed an MA for a mixed CARP with forbidden turns. In [9], this algorithm is extended to tackle turn penalties, parallel links, a range constraint and other objective functions like the duration of the longest trip. The average deviation to lower bounds is also improved by performing short restarts based on a partial replacement of the population. Both algorithms were evaluated on classical UCARP instances only. The last version [9], which is now one of the most efficient solution methods for the CARP, is briefly recalled. A chromosome is an ordered list of the required links, in which each required edge may appear as one of its two directions. Implicit shortest paths are assumed between successive required links. The chromosome does not include trip delimiters and can be viewed as a giant tour for a vehicle with infinite capacity. The fitness (CARP objective) is evaluated with the Split procedure (Section 3.5), which optimally partitions the giant tour into feasible routes. The initial population (nc chromosomes) comprises three good solutions computed by the basic versions of path-scanning, augment-merge and Ulusoy’s heuristic, completed by random permutations. Clones (identical chromosomes) are forbidden. At each iteration, two parents are selected by binary tournament and reproduce according to a modified version of the classical order crossover (OX). One child is randomly selected, the other is discarded. The mutation is replaced by a local search procedure, called with probability pm , which works on the individual routes instead of the giant tour. The moves include the removal of one or two consecutive required links from a route, with reinsertion at another position, the exchange of two required links, and 2-opt moves. All moves may involve one or two routes and the two traversal directions of a required edge are tested in the reinsertions. Each iteration scans all these moves in O(2 ) to find the first improving move. The
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
required arc
u
3373
y
trip S(i) v
x
trip S( j)
required edge
depot
Fig. 4. Example of 2-opt move on two trips.
whole local search stops when no more improvement can be found. The trips are then concatenated into a chromosome (giant tour), which is re-evaluated by Split. The resulting child replaces one of the nc /2 worst chromosomes, randomly chosen in the population and such that no clone is created. The MA stops after a given number of iterations ni , after a stabilization period of ns iterations without improving the best solution or when reaching a lower bound. When the lower bound is not reached, nr short restarts are performed with an intensified local search (rate pmr > pm ). Each restart begins with a special procedure which replaces ncr chromosomes by new ones (ncr < nc ), while preserving the best solution. It stops when the lower bound is achieved or after a given number of iterations nir (nir < ni ). 4.2. Adaptation for the extended MCARP The standard versions of PS, AM and Ulusoy’s heuristic that provide the initial population with three good solutions are replaced by better variants described in the previous subsections. The best variants are selected in the computational testing (Section 5). Multiple dumping sites and windy edges are easily added, they affect only the Split procedure and the local search. The MA of [9] tests the 2-opt moves only for the trips equivalent to their inverses (e.g., no required arc), to keep a low complexity in O(2 ). The main modification in the new MA is to scan the 2-opt moves in all cases, while keeping the O(2 ) complexity. Fig. 4 gives an example for two trips S(i) and S(j ). Required links are in boldface. Arrows indicate required arcs. Thin lines represent shortest paths. The path between two required links u and v is denoted u → v. We call 2-opt move the replacement of two deadheading paths by two others, like the replacement of u → v and x → y by u → inv(x) and inv(v) → y in the figure. Note that two sections are inverted: the end of S(i), from v onward, and the beginning of S(j ), up to x included. This is possible because they have no required arc. In general, the inverted sections have different costs, due to the windy edges and multiple dumping sites accepted in the extended MCARP of Section 3.1. For any trip and any required link u in , let rev( ) denote the inverted trip, arcsto( , u) the number of required arcs counted from the beginning of up to u included and arcsfrom( , u) the number of required arcs counted from u (included) up to the end of the trip. Similarly, let loadto( , u) denote the total demand up to u, loadfrom( , u) the cumulated demand from u onward, costto( , u) the cost of
3374
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
up to u, and costfrom( , u) the cost from u onward. The load and cost of u are included in the four last functions. Using this notation, the move in Fig. 4 is possible if arcsfrom(S(i), v)=arcsto(S(j ), x)=0 (no required arc in the sections to be inverted) and if loadto(S(i), u) + loadto(S(j ), x) Q and loadfrom(S(i), v) + loadfrom(S(j ), y) Q (vehicle capacity is respected). In that case, the cost variation is d(u, inv(x)) + d(inv(v), y) − d(u, v) − d(x, y) + costto(rev(S(i)), inv(v)) − costfrom(S(i), v) + costfrom(rev(S(j )), inv(x)) − costto(S(j ), x). Compared to the UCARP, these computations are much more involved. Nevertheless, it is still possible to search in O(2 ) an improving 2-opt move for the incumbent solution. Roughly speaking, our implementation scans all pairs of trips (S(i), S(j )), each required link u in S(i) and each required link x in S(j ), giving O(2 ) pairs (u, x) in total. At the beginning of the search, the values of arcsfrom, arcsto, loadto, loadfrom, costto and costfrom are initialized for each trip. The complexity of this pre-processing is O() only. During the search, when u (or v) is replaced by its successor in its trip, these values are incrementally updated, in O(1). 5. Computational evaluation 5.1. General The heuristics were programmed in the Pascal-like language Delphi 7. The lower bound was implemented in Visual C++ because it calls the linear programming library CPLEX. All algorithms were benchmarked on a 2 GHz Pentium IV PC with Windows 2000. Sections 5.2 and 5.3 propose two sets of instances for the basic MCARP, without any extension, to allow comparisons with the lower bound. Our goal is also to provide the other researchers with standard MCARP instances to benchmark their own solution methods. The two sets of instances can be downloaded from the Internet [21] or requested by e-mail. Sections 5.4–5.6 compare the different MCARP versions of PS, AM and Ulusoy’s heuristic (described in Section 3) on the two sets of instances. The rankings are confirmed on three sets of classical UCARP instances, see [22] to get the files: 23 gdb files from Golden et al. [4], 34 val files used by Belenguer and Benavent [12], and 24 egl instances derived by the same authors from data supplied by Eglese. The goal is to determine the best variants to provide the MA with three good initial solutions. Section 5.7 describes the parameters used by the memetic algorithm and by the CPA for the two sets of MCARP instances. The solution costs and lower bound values are detailed for each instance, together with the associated running times. 5.2. The mval instances These instances are derived from 34 CARP instances, called val files, designed by Benavent et al. [23] and used for instance by Belenguer and Benavent [12], Hertz et al. [7], Beullens et al. [11] and Lacomme et al. [9]. The original networks contain 24–50 nodes and 34–97 edges, all required. Each of them was converted into an mval one (“m” like mixed), with the same total demand, as follows. Each edge was kept with probability 0.4, replaced by two opposite arcs (same cost, but demands halved) with probability 0.4, and replaced by one single arc (same cost and demand) with probability 0.2, while preserving the strong connectivity. The resulting instances contain 24–50 nodes and 43–138 links, all required.
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
(a) 1
Vertical
Horizontal
1st diagonal
2nd diagonal
3375
(b) 2
3
4
5 1
6
7
8
9
2
10
4
3
6
8
8
9
15 11
12
16
17
18
19
13
12
11
20
21
18
19
17
10
13 14
16
14
5
15
20 21
(c)
(d)
Fig. 5. Random generation of lpr instances: (a) initial grid with basic squares, (b) random modifications of squares, (c) planar result with node numbers, (d) random distortions.
5.3. The lpr instances These large, sparse and planar networks are intended to mimic the shape of real street networks met in municipal waste collection. The costs correspond to deadheading and collecting times. The random generation starts with a rectangular grid made of basic squares defined by a given side length side, see (a) in Fig. 5. Four kinds of modifications are then applied with fixed probabilities pv, ph, pd1, pd2 to each square, see (b): split vertically (V), horizontally (H), along the first diagonal (D1) and along the second one (D2). Note that V and H create two new nodes. A square can undergo up to four modifications. As from two, a central node is created to preserve planarity, like node 7 in the figure. This step gives a planar network of 2-way streets, see (c), which can be already encoded as a symmetric directed graph G∗ like in Section 3.2. Since this graph is still too regular, each node is randomly moved in a small circle specified by a radius radius, see (d). Provisional lengths in meters can be computed from node coordinates. To simulate curved streets without drawing them, a second perturbation consists of applying a random growth factor growth to each length (between 0 and 10% for instance). Each 2-way street is then converted with a given probability pone into a 1-way street, by deleting at random one of its two arcs while checking strong connectivity. For each remaining arc u, we compute a deadheading cost c(u) from the length and a given vehicle speed speed, we draw the status (required or not) with a probability preq and assign to each required arc a non-zero random demand r(u) drawn in a given interval [qmin, qmax]. Finally, for each 2-way street, we randomly decide with probability pedge
3376
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
Table 1 Average deviations of path-scanning and its variants to lower bounds in % Files
PS
PSRC(20)
PSRC(50)
Break-even
PSRL(20)
PSRL(50)
Break-even
gdb val egl mval lpr
10.98 16.19 26.02 25.59 3.57
7.49 13.96 22.85 23.04 3.57
5.45 11.78 21.66 21.09 3.29
5 solns 11 solns 4 solns 8 solns 21 solns
8.77 14.65 22.26 25.29 3.77
7.23 12.60 20.98 22.17 3.48
7 solns 10 solns 4 solns 16 solns 42 solns
if it must be considered as one edge (bilateral waste collection). If yes, the two arcs are linked with the inv pointer (see Section 3.2) and the edge demand is defined as the total demand of its two arcs. The service cost w(u) comprises the deadheading time and a duration resulting from bin manipulations. It is computed as c(u) + r(u)/binload × bintime. The constant binload defines the average amount of waste per bin, the quotient corresponds to the number of bins in the street, and the constant bintime is the average time spent by the crew to unload a bin into the truck. The random generation ends by drawing the vehicle capacity W, the depot-node s and a non-zero dumping time . This dumping time (see Section 3.2) is the only extension to the basic MCARP accepted in the lpr files. The other extensions like forbidden turns, windy edges or parallel links are not included and the dumping site coincides with the depot. Let rtot denote the total demand. The dumping time is counted in the lower bound by computing the minimum number of vehicles required k ∗ = rtot /Q and by adding
· k ∗ to the value obtained by the cutting plane algorithm. Fifteen files were generated with the following parameters: side = 200 m, ph = pv = 0.2.0.7, pd1 = pd2=0.1.0.2, radius=20 m, growth=0.1, pone=0.0.85, speed =20 km/h, preq=0.6.0.9, qmin=10 kg, qmax =400 kg, pedge=0.1.0.30, binload =10 kg, bintime=10 s, W =10 tons and =300 s. The time per bin seems small: in fact, the two workers behind the truck need 20 s each, giving an apparent throughput of one bin per 10 s. The resulting graphs contain 28–401 nodes, 52–1056 links and 50–806 required links. Hence, they include large scale instances, much larger than the mval networks and classical UCARP instances. The files are partitioned into three groups of five called a, b and c. In each group, the files are numbered in increasing order of size. The lpr-a files model modern towns with a majority of large 2-way streets and few edges: their graphs are rather directed and symmetric (pone=0.0.15, pedge 0.1). Lpr-b files correspond to old town centres, with narrow 1-way streets: their graphs are rather directed and anti-symmetric (pone 0.85, pedge 0.3). Low-traffic residential areas with a majority of streets with bilateral service (edges) are imitated by the lpr-c problems: their graphs are rather undirected (pone 0.25, pedge 0.85). The features of each instance are listed in the Appendix (Table 6). 5.4. Evaluation of path-scanning and its variants Table 1 compares the average deviations of PS and its variants PSRC and PSRL to good lower bounds. Recall (see Section 3.3) that PS builds 5 solutions (one per criterion), while PSRC and PSRL can generate any specified number of solutions. In each iteration, PSRC randomly selects the criterion while PSRL randomly takes a required link among the closest neighbours. PSRC and PSRL are tested with 20 and 50 solutions. The bounds are taken from [12] for all UCARP instances. They are computed by the cutting
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
3377
Table 2 Average deviations of augment-merge and its variants to lower bounds in % Files
AM
M
AIM
IM
gdb val egl mval lpr
7.37 10.73 23.11 16.04 2.40
6.54 12.54 6.77 16.96 2.20
7.11 9.95 23.03 14.39 2.47
5.61 6.01 5.89 9.93 1.62
Table 3 Average deviations of splitting heuristics to lower bounds in % Files
UPS
URC(20)
URL(20)
IURL(20)
gdb val egl mval lpr
8.36 14.53 17.43 21.27 3.14
5.28 11.56 14.86 17.73 2.79
4.56 11.39 14.65 17.34 2.77
3.99 9.60 13.35 14.30 2.46
plane algorithm and listed in Appendices 1 and 2 for the mval and lpr files. The Break-even columns indicate the average number of solutions that must be examined by these heuristics to outperform PS. The results show that randomised versions can easily do better than PS, provided they inspect enough solutions, and that the gain increases with the number of solutions. Of course, the price to pay is a longer running time, roughly proportional to the number of solutions computed. For the same number of solutions, PSRC outperforms PSRL, except on egl files. 5.5. Evaluation of augment-merge and its variants Table 2 indicates the average deviations of AM, M, AIM and IM (in %) to the lower bounds mentioned in Section 5.4 and reveals two surprising conclusions: (a) the Augment phase is in general useless (except for the val and mval files, where AM outperforms M) and (b) the improved version of merge without augment (IM) strongly improves M and gives by far the best results. The augment phase is disastrous on the egl files, in which some edges are not required, contrary to the other UCARP benchmarks. As discussed in Section 3.4, this favours the creation of internal deadheads in the trips during the Augment phase. Moreover, the average demand per required edge in these files is large compared to vehicle capacity (around 1/7). Therefore, the total number of Augment or Merge operations is limited, and running the Augment phase lets very few feasible mergers to the Merge phase. 5.6. Evaluation of Ulusoy’s heuristic and its variants Table 3 compares UPS (5 giant tours) to URC and URL with 20 giant tours. Recall (see Section 2.5) that these heuristics apply the splitting procedure to giant tours obtained by calling PS, PSRC and PSRL with an infinite vehicle capacity. In Table 1, PSRC gives better solutions than PSRL, but here, surprisingly, split produces better results on the giant tours supplied by PSRL (heuristic URL).
3378
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
The best version URL is improved further in the fifth column (IURL), by allowing the rotations of sub-sequences when building the auxiliary graph, as explained at the end of Section 3.5. Comparisons with Tables 1 and 2 indicate that IURL is the best heuristic on the gdb files. It is outperformed only by IM on the other benchmarks, and also by M on the egl and lpr instances.
5.7. The memetic algorithm and the lower bound The heuristics PSRC(20), IM and IURL(20) are called to include three good solutions in the MA initial population. The standard setting for the MA parameters of Section 4.1 is as follows: population size nc = 30, local search rate pm = 0.1, stop after ni = 20, 000 crossovers, ns = 6000 crossovers without improving the best solution, or when the lower bound is reached. If the bound is not achieved, the MA performs up to nr = 20 restarts with an augmented local search rate pmr = 0.2. Each restart begins by replacing ncr = 7 chromosomes in the population and ends after nir = 2000 crossovers or when reaching the bound (in that case, the sequence of restarts is interrupted). A time limit of 1 h is imposed to the whole algorithm. The lower bound is computed with the CPA described in Section 2.2. Nevertheless, the identification of DP inequalities is inhibited for the largest lpr instances (more than 500 required links) because it is very time-consuming and it has been observed that it improves the bound on a small number of instances only. The detailed results are given in Appendices 1 (mval files) and 2 (lpr files). The first columns respectively indicate the file name, then the numbers of nodes, links, edges and arcs, then the numbers of required links, edges and arcs. The numbers of links, edges and arcs are suppressed in the mval table because all links are required. The next columns provide the cost and the running time for the lower bound LB, for the three heuristics PSRC(20), IM and IURL(20), and for the MA with the standard setting of parameters. The last column give the best result found for the MA during the testing, using various settings. The simple heuristics are very fast in all cases. The Merge heuristic IM is by far the best one on average (9.93% on the mval and 1.62% on the lpr), but also the slowest one on the largest instances (nearly 9 s for lpr-c-05), simply because it is implemented in O(3 ) instead of the better complexity O(2 . log ) mentioned in Section 3.4. The average deviation to the lower bound is smaller for the lpr instances. This can be explained by the relatively small deadheading times, compared to the uncompressible collecting times. The comparison between the MA and the lower bound LB proves their effectiveness in practice: using the standard setting of parameters, the average deviation is 0.72% for the mval files, with 19 proven optima (when the bound is reached), and 0.42% for the lpr files, with 6 proven optima. Using several settings (MA best column), it was possible to decrease the deviations down to 0.51% (with 23 optima) and 0.33% (with 6 optima), respectively. All proven optima are found very quickly. When the bound is not reached, the execution times are reasonable on the mval files (5 min maximum) but still excessive on the biggest lpr (numbered 4 and 5), for which the MA reaches its maximum time limit. However, these networks are enormous compared to the largest graphs used in the arc routing literature (up to 1200 links instead of 190), and our implementations could be accelerated, using for instance neighbourhood reduction techniques like the lists of neighbours described in [11]. Similar comments apply to the computation of LB with the cutting plane algorithm: execution times are negligible on the mval files but too large on the biggest lpr.
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
3379
Table 4 Average deviations of the MA to lower bounds for different time limits in % Files
mval lpr
Best initial soln
1 min
5 min
Deviation
Optima
Deviation
Optima
Deviation
Optima
9.21 1.43
0 0
0.81 0.85
16 6
0.72 0.66
19 out of 34 6 out of 15
In fact, the MA still provides excellent solution values for small time limits, as shown in Table 4 for the standard setting of parameters. The best solution in the initial population is always one of the three heuristic solutions, never a random one. After 1 min, the average deviation to the lower bounds drops below 1% for the two sets of instances, and the MA has already discovered the six proven optima for the lpr files, and 16 out of 20 for the mval. 6. Conclusion This paper proposes a relaxed linear model for the Mixed CARP and a cutting plane algorithm to solve it. Several simple heuristics, derived from classical UCARP heuristics, and one memetic algorithm are also described. They are designed to tackle complications like windy edges or forbidden turns. Since it seems very difficult to include these complications in linear programming models, the lower bounds and upper bounds are compared on two sets of basic MCARP instances, i.e., without the complications. The simple heuristics are very fast, even on the largest instances, but the average deviation to the lower bound may reach 23%. The study also validates the improvements suggested for path-scanning, augment-merge and Ulusoy’s heuristic, even on standard UCARP instances. In particular, it reveals the bad influence of the Augment phase in augment-merge and suggests the use of the variants PSRC, M and IURL as the best heuristics to initialize the MA. The MA reaches the CPA lower bound on two thirds of mval files and almost one half of lpr. The average deviation to the bound is at most 0.5% but the lpr instances with more than 500 links require 1 h of computation time. Nevertheless, 1 min is enough to obtain an average deviation below 1%. Moreover, the results obtained on the large instances (lpr files) are very encouraging to envisage practical applications. Acknowledgements The contribution by J.M. Belenguer has been partially supported by the Ministerio de Ciencia y Tecnologia of Spain through project HF2001-0071. The contribution by E. Benavent has been partially supported by the Ministerio de Ciencia y Tecnologia of Spain through projects HF2001-0071 and TIC2000-C06-01. P. Lacomme and C. Prins were partially funded by the France–Spain Integrated Action 04283QB from the Ministère des Affaires Etrangères of France. Appendix Tables 5 and 6.
n
24 24 24 24 24 24 24 24 24 41 41 41 41 34 34 34 34 31 31 31 40 40 40 30 30 30 50 50 50 50 50 50
File
mval1a mval1b mval1c mval2a mval2b mval2c mval3a mval3b mval3c mval4a mval4b mval4c mval4d mval5a mval5b mval5c mval5d mval6a mval6b mval6c mval7a mval7b mval7c mval8a mval8b mval8c mval9a mval9b mval9c mval9d mval10a mval10b
55 51 53 44 52 49 48 45 43 95 102 103 104 96 91 98 92 69 66 68 86 91 90 96 91 83 132 120 125 131 138 134
m=
20 13 17 16 12 14 15 16 18 26 19 21 21 22 35 17 29 22 22 23 36 25 28 20 27 28 32 44 42 38 32 33
35 38 36 28 40 35 33 29 25 69 83 82 83 74 56 81 63 47 44 45 50 66 62 76 64 55 100 76 83 93 106 101
Table 5 Computational results for mval instances
230 261 309 324 395 521 115 141 166 580 650 630 746 597 613 697 719 326 317 365 364 412 424 581 531 617 458 453 428 514 634 661
LB cost 0.11 0.11 0.11 0.44 0.05 0.11 0.33 0.16 0.17 0.87 0.27 4.12 0.82 0.22 0.39 0.38 0.49 0.11 0.22 0.44 0.99 0.44 1.53 0.16 0.22 0.39 0.38 0.71 24.17 9.99 3.30 1.93
LB time 243 314 427 409 471 644 133 162 191 699 775 828 1015 733 718 809 883 392 406 526 439 507 578 666 619 842 529 552 529 695 735 753
PSRC(20) cost 0.016 0.013 0.015 0.012 0.014 0.013 0.013 0.012 0.011 0.031 0.033 0.034 0.033 0.031 0.030 0.031 0.030 0.020 0.019 0.020 0.027 0.029 0.028 0.031 0.030 0.026 0.049 0.044 0.046 0.048 0.052 0.050
PSRC(20) time 243 276 352 360 407 560 119 163 174 653 693 702 810 686 677 743 821 370 346 402 381 470 451 639 568 718 500 534 479 575 710 717
IM cost 0.004 0.003 0.004 0.002 0.003 0.003 0.003 0.003 0.002 0.018 0.020 0.020 0.022 0.017 0.017 0.017 0.016 0.007 0.007 0.007 0.015 0.015 0.015 0.017 0.016 0.013 0.043 0.037 0.040 0.044 0.047 0.043
IM time 231 292 357 374 434 601 128 150 192 684 737 740 905 683 677 811 855 367 354 444 390 491 504 651 611 762 514 502 498 622 705 714
IURL(20) cost 0.004 0.005 0.004 0.006 0.004 0.004 0.004 0.005 0.002 0.012 0.011 0.013 0.012 0.012 0.011 0.013 0.008 0.006 0.004 0.006 0.011 0.009 0.006 0.011 0.012 0.008 0.018 0.016 0.013 0.011 0.025 0.017
IURL(20) time
MA time 0.10 0.27 34.08 0.20 0.31 31.57 0.45 53.36 6.65 17.82 2.23 147.30 134.25 19.40 141.14 81.07 113.30 11.65 15.95 58.18 0.75 7.47 113.40 76.78 6.02 77.84 24.00 17.95 246.28 253.22 128.60 313.70
MA cost 230∗ 261∗ 315 324∗ 395∗ 526 115∗ 142 166∗ 580∗ 650∗ 631 776 597∗ 615 697∗ 757 326∗ 317∗ 375 364∗ 412∗ 428 581∗ 531∗ 638 458∗ 453∗ 434 520 634∗ 662
230∗ 261∗ 315 324∗ 395∗ 526 115∗ 142 166∗ 580∗ 650∗ 630∗ 770 597∗ 613∗ 697∗ 739 326∗ 317∗ 371 364∗ 412∗ 426 581∗ 531∗ 638 458∗ 453∗ 429 520 634∗ 661∗
MA best
3380 J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
100 87
623 643
2.74 2.31
0
44.11%
23.04%
751 847
0.052s
0.029s
0.051 0.048
0
17.88%
9.93%
680 706
0.047s
0.019s
0.047 0.044
0
23.50%
14.30%
714 760
0.025s
0.010s
0.018 0.015
19
5.29%
0.72%
624 650
314.05s
79.89s
314.05 267.05
23
3.40%
0.51%
623∗ 649
Notes: All links required. Cost columns: average and worst deviations to LB in %. Time columns: seconds on a 2 GHz PC, Windows 2000. Asterisks denote proven optima.
Optima
24.17s
36 42
Worst
136 129 1.74s
50 50
Average
mval10c mval10d
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383 3381
LB cost
LB time 13597 28377 77331 128566 207597 14918 29285 80591 129449 215883 18744 36485 112462 171823 262089
0
1.028s 3.54%
0.330s 1.62%
0.019 0.044 0.216 0.424 0.966 0.016 0.038 0.183 0.382 0.798 0.018 0.043 0.217 0.549 1.028
PSRC IM (20) cost time 13537 28586 78151 131884 212167 14868 28947 79910 132241 219702 18706 36763 114539 173161 266058
IURL (20) cost
0
8.945s 4.72%
1.687s 2.46%
0.002 0.012 0.308 1.288 5.201 0.002 0.012 0.296 1.146 4.910 0.003 0.021 0.567 2.595 8.945
IM time
1.27 1.48 1824.25 3600.00 3600.00 0.06 1.14 1654.17 3600.00 3600.00 0.31 10.05 3600.00 3600.00 3600.00
13484∗ 28052∗ 76155 127930 206086 14835∗ 28654∗ 77878 127454 212279 18639∗ 36339∗ 111632 169487 260538
6
0.311s 1.65%
13484∗ 28052∗ 76155 127352 205499 14835∗ 28654∗ 77878 127454 211771 18639∗ 36339∗ 111632 169254 259937
MA best
6
3600.0s 1.36%
1912.8s 0.33%
MA time
MA cost
0.104s 0.42%
0.010 0.020 0.070 0.140 0.301 0.000 0.010 0.070 0.120 0.270 0.010 0.010 0.060 0.151 0.311
IURL (20) time
Notes: Cost columns: average and worst deviations to LB in %. Time columns: in seconds on a 2 GHz PC with Windows 2000. Asterisks denote proven optima.
0
3557.80s 6.52%
Worst
Optima
727.67s 3.57%
13600 29094 79083 133055 215153 15047 29522 80017 133954 223473 18897 36929 115763 174416 268368
PSRC (20) cost
Average
28 94 0 94 52 0 52 13484 0.11 53 169 5 164 104 5 99 28052 0.93 146 469 33 436 304 33 271 76108 58.44 195 651 34 617 503 34 469 126941 339.28 321 1056 58 998 806 58 748 202735 2424.46 28 63 5 58 50 5 45 14835 0.11 53 117 9 108 101 9 92 28654 0.88 163 361 26 335 305 26 279 77837 56.80 248 582 8 574 501 8 493 126932 412.66 401 876 37 839 801 37 764 209791 3032.00 28 52 39 13 50 39 11 18639 0.11 53 101 77 24 100 77 23 36339 1.48 163 316 241 75 302 241 61 111117 181.59 277 604 362 242 504 362 142 168441 848.38 369 841 387 454 803 387 416 257890 3557.80
|A|
lpr-a-01 lpr-a-02 lpr-a-03 lpr-a-04 lpr-a-05 lpr-b-01 lpr-b-02 lpr-b-03 lpr-b-04 lpr-b-05 lpr-c-01 lpr-c-02 lpr-c-03 lpr-c-04 lpr-c-05
|E|
n
File
m
Table 6 Computational results for lpr instances
3382 J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
J.-M. Belenguer et al. / Computers & Operations Research 33 (2006) 3363 – 3383
3383
References [1] Golden BL, Wong RT. Capacitated arc routing problems. Networks 1981;11:305–15. [2] Assad AA, Golden BL. Arc routing methods and applications. In: Ball MO., et al., editors. Handbooks in OR and MS, vol. 8. Amsterdam: Elsevier; 1995. p. 375–483. [3] Hirabayashi R, Saruwatari Y, Nishida N. Tour construction algorithm for the capacitated arc routing problem. Asia-Pacific Journal of Operational Research 1992;9:155–75. [4] Golden BL, Dearmon JS, Baker EK. Computational experiments with algorithms for a class of routing problems. Computers and Operations Research 1983;10(1):47–59. [5] Ulusoy G. The fleet size and mix problem for capacitated arc routing. European Journal of Operational Research 1985;22:329–37. [6] Eglese RW, Li LYO. A tabu search based heuristic for arc routing with a capacity constraint and time deadline. In: Osman IH, Kelly JP, editors. Metaheuristics: theory and applications. Kluwer; 1996. p. 633–50. [7] Hertz A, Laporte G, Mittaz M. A tabu search heuristic for the Capacitated Arc Routing Problem. Operations Research 2000;48(1):129–35. [8] Lacomme P, Prins C, Ramdane-Chérif W. A genetic algorithm for the capacitated arc routing problem and its extensions. In: Boers EJW, et al., editors. Applications of evolutionary computing. Lecture notes in computer science, vol. 2037. Berlin: Springer; 2001. p. 473–83. [9] Lacomme P, Prins C, Ramdane-Chérif W. Competitive memetic algorithms for arc routing problems, Research Report LOSI-2001-01, University of Technology of Troyes, France, 2001. Annals of Operations Research 2004;131:159–85. [10] Greistorfer P. A tabu scatter search metaheuristic for the arc routing problem. Computers and Industrial Engineering 2003;44:249–66. [11] Beullens P, Muyldermans L, Cattrysse D, Van Oudheusden D. A guided local search heuristic for the capacitated arc routing problem. European Journal of Operational Research 2003;147(3):629–43. [12] Belenguer JM, Benavent E. A cutting plane algorithm for the capacitated arc routing problem. Computers and Operations Research 2003;30(5):705–28. [13] Nobert Y, Picard JC. An optimal algorithm for the mixed Chinese postman problem. Networks 1996;27:95–108. [14] Corberán A, Romero A, Sanchis JM. The mixed general routing polyhedron. Mathematical Programming 2003;96: 103–37. [15] Ford LR, Fulkerson DR. Flows in networks. Princeton, NJ: Princeton University Press; 1962. [16] Benavent E, Corberán A, Sanchis JM. Linear programming based methods for solving arc routing problems. In: Dror M, editor. Arc routing, theory, solutions and applications. Boston: Kluwer Academic Publishers; 2000. p. 231–75. [17] Padberg MW, Rao MR. Odd minimum cut-sets and b-matchings. Mathematics of Operations Research 1982;7:67–80. [18] Benavent E, Soler D. The directed Rural Postman Problem with turn penalties. Transportation Science 1999;33(4): 408–18. [19] Cormen TH, Leiserson CE, Rivest RL, Stein C. Introduction to algorithms. 2nd ed., Cambridge, MA: The MIT Press; 2001. [20] Moscato P. Memetic algorithms: a short introduction. In: Corne D, Dorigo M, Glover F, editors. New ideas in optimization. McGraw-Hill; 1999. p. 219–34. [21] Belenguer JM, Benavent E, Lacomme P, Prins C, Ramdane-Chérif W. Directory of MCARP instances, 2003. http://www.uv.es/∼belengue/mcarp/index.html. [22] Belenguer JM, Benavent E. Directory of UCARP instances, 1997. http://www.uv.es/∼belengue/carp.html. [23] Benavent E, Campos V, Corberán A, Mota E. The Capacitated Arc Routing Problem: lower bounds. Networks 1992;22: 669–90.