The traveling purchaser problem with budget constraint

Computers & Operations Research 36 (2009) 2263 -- 2274

Contents lists available at ScienceDirect

Computers & Operations Research journal homepage: w w w . e l s e v i e r . c o m / l o c a t e / c o r

The traveling purchaser problem with budget constraint Renata Mansini, Barbara Tocchella ∗ Department of Electronics for Automation, University of Brescia, Italy

A R T I C L E

I N F O

Available online 11 September 2008 Keywords: Traveling purchaser problem Budget constraint Local search Multi-start variable neighborhood search

A B S T R A C T

Let us consider a set of markets plus a depot and a set of products for each of which a positive demand is specified. Each product is made available in a subset of markets in each of which only a given quantity, less than or equal to the required one, can be purchased at a given unit price. The distance between each couple of markets and between each market and the depot is known. The Traveling Purchaser Problem with Budget constraint (TPP-B) looks for a simple cycle starting at and ending to the depot and visiting a subset of markets at a minimum traveling cost and such that the demand for each product is satisfied and the cost globally spent for purchasing the products does not exceed a defined budget threshold. As the TPP also this problem arises in several application domains, but while the former has been largely studied, very few contributions exist in the literature for the TPP-B. We propose and compare two solution algorithms, an enhanced local-search heuristic and a variable neighborhood search (VNS) approach also tested in a multi-start variant. The proposed algorithms have been used to solve both the capacitated and the uncapacitated version of the problem. Test problems have been obtained by adding a budget constraint to known benchmark instances for the TPP. © 2008 Elsevier Ltd. All rights reserved.

0. Introduction The capacitated Traveling Purchaser Problem with Budget constraint (capacitated TPP-B) can be described as follows. Consider a set of markets M := {1, . . . , m} plus a depot (indexed 0) and a set of products K := {1, . . . , n}. Each product k, k ∈ K, can be purchased at a subset Mk ⊆ M of markets at a non-negative unit price fki depending on the market i, i ∈ Mk . A demand dk of units is specified for each product k ∈ K and only qki units of product k are available  at market i, i ∈ Mk . We assume that 0 < qki  dk and i∈M qki  dk k

for all k ∈ K and i ∈ Mk . A traveling cost between each couple of markets and between each market and the depot is known and given by cij , ∀i, j ∈ M ∪ {0}. A budget B is also specified and represents a threshold on the total amount which can be spent on purchasing products. The problem consists in selecting a subset of markets and visiting them with a tour starting at and ending to the depot such that the demand dk for each product k ∈ K is satisfied at a global purchasing cost which does not exceed the given budget threshold B while minimizing the total traveling cost. The special case where dk = 1 ∀k ∈ K and qki = 1 for each i ∈ Mk , defines the uncapacitated version of the problem (uncapacitated TPP-B).

∗

Corresponding author. Tel.: +390 30 3715512; fax: +390 30 380014. E-mail addresses: [email protected] (R. Mansini), [email protected], [email protected] (B. Tocchella). 0305-0548/$ - see front matter © 2008 Elsevier Ltd. All rights reserved. doi:10.1016/j.cor.2008.09.001

The TPP-B is a problem involving two types of costs. The first one is the traveling cost which is minimized in the objective function. The second one is the purchasing cost, the value of which is bounded in the constraints. Most of the previous studies on similar problems have restricted the analysis to the case where the two cost measures were both included in the objective function. This single-objective combinatorial optimization problem is known as TPP and has found several applications in scheduling and routing contexts (see Ramesh [1] and Singh and van Oudheusden [2]) as well as in the telecommunication network design (see Ravi and Salman [3]). The TPP is known to be NP-hard in the strong sense, reducing to the traveling salesman problem (TSP) when each market offers a product which is not provided by the remaining ones. It also contains the uncapacitated facility location problem and the set covering problem as special cases. Due to its computational hardness several heuristic procedures have been proposed to solve it (see, for instance, Golden et al. [4], Ong [5], Pearn and Chien [6], Vo [7], Renaud et al. [8] and Mansini et al. [9]) while recently a branch-and-cut procedure has allowed the optimal solution of problems with up to 200 markets and 150 products (see Laporte et al. [10] and references therein for other exact methods previously proposed in literature for the same problem). An evident drawback of the TPP is the possible different nature of the two costs merged in the objective function: the traveling cost is usually represented as distance or traveling time whereas the purchasing cost is typically measured as currency. Moreover, a trade-off frequently characterizes the two objectives since reducing the purchasing cost may imply an increase in the distance traveled. To avoid

2264

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

such disadvantage Riera-Ledesma and Salazar-González (see [11]) introduce the bi-objective TPP: this is a bi-criteria version of the TPP where minimizing the purchasing cost and the traveling distances are two separate objectives. The bi-objective is a generalization of the TPP whose solution provides insight into the trade-off between the two costs. The authors tackle the problem by generating the set of all supported and non-supported efficient points in the objective space. For each efficient point a Pareto optimal solution in the decision space is computed by solving a single-objective problem while bounding the remaining cost function. In the case the single objective is represented by traveling costs while purchasing costs are bounded in the constraints, the problem solved is the TPP-B. The authors provide exact solutions by means of a branch-and-cut approach for instances with up to 100 markets and 200 different products but only for the uncapacitated version of the problem. The uncapacitated TPP-B is also analyzed in Ravi and Salman [3] as an application of telecommunication network design. In that paper, the authors introduce an algorithm for the bi-criteria version of the UTPP with metric distances based on the rounding of an LP relaxation solution whose worst-case ratio is poly-logarithmic. They also show that, for a special case of the UTPP which models the ring-star network design problem with proportional costs, a constant factor approximation algorithm exists. In the present paper, we propose two different algorithms to solve the TPP-B. To the best of our knowledge these procedures represent the first heuristics proposed in the literature to solve this bounded version of the TPP. The first procedure is based on an enhanced localsearch scheme where a neighbor solution is obtained by the current one by removing l-consecutive markets and inserting as many markets as required to restore solution feasibility. Since the value of the parameter l is changed during the search, the procedure is more than a simple local search. The idea of removing l-consecutive vertices is not new in the literature. One of the most effective heuristic algorithm for the TPP (see Riera-Ledesma and Salazar-González [12]) makes use of such idea obtained as generalization of the procedures proposed by Voss [7] and Keller [13]. With respect to this heuristic, we introduce a tabu structure characterized by a variable length tabu list to help the search to escape from local optimal solutions and a problem-specific procedure to improve the objective function value without losing feasibility. Moreover, in our algorithm the maximum number of consecutive markets which can be removed is not predefined as in [12] but changes during the search. Extensive computational results have shown that our heuristic is very efficient and, in those instances where budget constraint is not too tight (i.e. not too close to the minimum purchasing cost), it is also able to provide very effective solutions. It represents the immediate choice for a decision maker interested in a practical tool to obtain good solutions in few seconds. The second procedure is based on a variable neighborhood search (VNS) scheme (see Mladenovic and Hansen [14] and references therein). We generate a sequence of neighborhoods Nq (·), q = 1, . . . , qmax . For a given q and a given solution s, a solution s belonging to the neighborhood Nq (s) is generated from s by inserting q random markets. Then a local search procedure is applied to s . At this aim we use a variant of the proposed enhanced local search algorithm. If no solution improvement is obtained, the parameter q is increased allowing the search to explore farthest and possibly more promising parts of the solution space. We also implement a variant of this VNS scheme where, if no solution improvement is obtained, a maximum number tmax of random solutions is generated in the same neighborhood Nq (s) before the value of q can be changed (Multi-start VNS). The proposed VNS procedure and its variant are straightforward and produce very effective solutions. To test all the proposed algorithms we have introduced a broad class of benchmark problems obtained by inserting budget thresholds to the TPP instances proposed in [10]. We have solved instances with both tight

and untight budget constraints and tested both the capacitated and the uncapacitated version of the problem. Our experimental analysis provides the first computational results available for the TPP with budget constraint. Both optimal (when available) and heuristic solutions of the solved benchmark problems can be downloaded at the web site http://www.ing.unibs.it/or-group/testproblems.htm. The paper is organized as follows. In Section 1, we present the mathematical formulation of the TPP with budget constraint. The proposed solution methods are introduced in Section 2 where we also provide the description of a lower bound for the problem used to evaluate our algorithms performance. Section 3 is devoted to computational results. For the uncapacitated case we have used the instances solved in [11], where the authors deal with the UTPP with budget constraint to compute Pareto optimal solutions of a bi-objective TPP. For the capacitated version of the problem, since no optimal solutions are available, we compare the algorithms each other and with respect to the problem lower bound in terms of efficiency and effectiveness. Finally conclusions and future developments are drawn in Section 4.

1. Problem formulation The capacitated TPP-B can be formally defined on an undirected complete graph G = (V, E), where V := {0} ∪ M is the vertex set with vertex 0 representing the depot and E := {(i, j) : i, j ∈ V, i < j} is the edge set. The problem looks for a simple cycle in G, starting at and ending to vertex 0, which visits a subset of vertices at a minimum traveling cost while purchasing all products according to predefined demands at a global purchasing cost which does not exceed a budget threshold B. To avoid explicit consideration of trivial cases it is assumed that none of the markets can provide each of the n products at the cheapest cost. To formulate the capacitated version of the problem the following three sets of decision variables are introduced:  xe :=  yi :=

1 0

if the edge e belongs to the optimal cycle otherwise

1 0

if market i belongs to the optimal cycle otherwise

∀e ∈ E

∀i ∈ V

zki  0 the amount of product k purchased at market i ∀k ∈ K and ∀i ∈ Mk . For any S ⊂ V, we define (S) as the set {(i, j) ∈ E : i ∈ S, j ∈ V\S}. Moreover, we indicate as M∗ , the set M∗ := {0} ∪ {i ∈ M :  ∃k ∈ K such that j∈M \{i} qkj < dk }, representing the markets which k necessarily have to make part of any feasible solution (see also [10]). The capacitated TPP with budget constraint can be formulated as follows: (TPP-B)

min

v :=



ce xe

(1)

e∈E



xe = 2yi

∀i ∈ V

(2)

e∈({i})



xe  2yi

∀S ⊆ M and ∀i ∈ S

(3)

e∈(S)

  k∈K i∈Mk

fki zki  B

(4)

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274



zki = dk

∀k ∈ K

(5)

i∈Mk

zki  qki yi yi = 1

∀ k ∈ K and ∀i ∈ Mk

∀i ∈ M∗

(6) (7)

xe ∈ {0, 1} ∀e ∈ E

(8)

yi ∈ {0, 1} ∀i ∈ V\M∗

(9)

zki  0 ∀k ∈ K and ∀i ∈ Mk

(10)

Objective function (1) establishes the minimization of the routing costs. The set of constraints (2) ensures that the degree of each market i belonging to the solution (i.e. yi = 1) has to be equal to two, i.e. each selected market is visited only once. Inequalities (3) ensure connectivity. Constraint (4) imposes a budget threshold B on the total purchasing cost. Equalities (5) guarantee that the exact amount dk for each product k is purchased, whereas for inequalities (6) it is not possible to purchase a product k at an unvisited market or to buy an amount larger than the offered quantity qki at a visited one. Equalities (7) impose that all the markets belonging to set M∗ are necessarily selected in any feasible solution. Constraints (8)–(10) are binary and non-negativity conditions on the variables. Notice that the uncapacitated TPP-B formulation can be straightforwardly obtained by setting dk =1 ∀k ∈ K and qki =1∀k ∈ K, ∀i ∈ Mk while variables zki now range between 0 and 1. 2. Solution algorithms Both proposed solution algorithms receive as input the same initial feasible solution sI , use the same routines to restore feasibility and apply the same procedure to possibly improve the objective function value of a feasible solution as described below. From now on, given a current solution s, we define as V s the set of markets in the solution and as external the markets which do not belong to it. Finally, we refer to Tabu as to the set of markets which cannot be removed from the solution in the current iteration. Feasibility restoring procedures. A solution s may result to be infeasible if the required demand for some products is not satisfied or if the total purchasing cost exceeds the threshold B. In the first case the GetAll(s) routine is called, whereas the BudgetRestore(s) procedure is applied in the second one. • Procedure GetAll(s). The procedure receives as input a solution s which is infeasible with respect to products demand. It iteratively inserts external markets one after the other until all products are bought in the required quantities. Markets are selected and inserted according to the cheapest insertion rule. At this aim, for each external market, we first evaluate its insertion cost between any two consecutive markets of the tour, then we insert it between the two markets corresponding to the minimum of the computed insertion costs. External market inducing the minimum increase of traveling cost is selected first. A formal description of the cheapest insertion rule and a proof of its time complexity equal to O(m2 log(m)) for the construction of a tour with m vertices can be found in Rosenkrantz et al. [15]. Since the procedure aims at satisfying products demand while minimizing routing costs, when selecting the markets the possible saving in the purchasing cost is not taken into account. All the inserted markets become tabu with a tabu tenure value r depending on the iteration. Every rmax iterations the tabu tenure r reduces to zero and the set Tabu becomes

2265

void. This means that, when entering a solution, a market will remain tabu for a number of iterations equal to the current value of parameter r, where r is less than or equal to rmax . Let V s be the set of markets belonging to the input solution s, then the time complexity of the procedure is given by the insertion of |V\V s | external markets by cheapest insertion. The output solution may not be feasible with respect to the budget threshold constraint. In such case the following BudgetRestore(s) procedure is applied. • Procedure BudgetRestore(s). The procedure receives as input a solution s which is infeasible with respect to the budget. As a first step, it determines all the subsets of external markets with cardinality one and two such that, if singularly inserted in the solution s, they guarantee the budget to be satisfied. In particular, the subset which allows the budget to be restored with the smallest increase in the traveling cost is selected. The routing cost increase is computed using the cheapest insertion rule (see [15]): if a subset with two markets is inserted the traveling cost increase is estimated as the sum of the costs obtained by singularly inserting each market. If no subset of cardinality less than or equal to two is enough to restore budget feasibility, in a second step the procedure iteratively inserts one by one the external markets, previously sorted according to their routing costs. The procedure stops when the budget constraint is satisfied. In case of tie, the market which allows the largest cost saving is selected first. The cost saving is optimally computed as the difference between the total purchasing cost before and after the market insertion. To efficiently evaluate the total purchasing cost of a solution after markets insertion we have implemented the data structure proposed in [12]. According to such data structure the purchasing cost saving produced by the insertion of a new market in a solution s can be directly computed as the sum of the savings obtained for each of the products offered by the new market requiring a computational time equal to O(|K| log(|Mk ∩ V s |)). Finally, markets inserted in a solution become tabu for a number of iterations equal to r, with r  rmax , as described in the previous routine. Given the input solution s, the computational time complexity of BudgetRestore(s) is O(|K| log(|Mk ∩ V s |)|V\V s |2 ) which is the computational time required in step one to find all subsets of two external markets and to verify if the purchasing cost saving, allowed by each of them, is enough to restore feasibility. In practice, such complexity is reduced by considering only the more promising external markets, i.e. those markets that are more likely to shrink the purchasing cost if inserted. At this aim, we identify the products for which their average price with respect to the markets in s is higher than the corresponding value computed over all the markets. We limit the search to the external markets offering such products. Initial solution. An initial solution sI is obtained by a two steps procedure. In the first step, the markets are added one by one to a solution until all products have been purchased (construction phase). We have used the nearest neighborhood rule to select the market which has to be inserted next and we have started insertion with markets belonging to set M∗ . If the constructed solution results to be infeasible with respect to the budget threshold B, the BudgetRestore() routine is applied. In the second step, the previously obtained solution is possibly improved by recomputing its traveling cost with a TSP-like heuristic (improving phase). At this aim, we have applied GENIUS algorithm. This procedure consists in a construction phase (GENI) the main feature of which is that the insertion of a vertex in a tour does not necessarily take place between two consecutive vertices and a post-optimization routine (US) based on different rules for removing a vertex from a tour and inserting it back (see Gendreau et al. [16] for details).

2266

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

Fig. 1. Pseudo-code of EJEMO algorithm.

Procedure DRF(s) . The idea that underlies the procedure Drop Retaining Feasibility (DRF(s)) is to try and drop as many markets as possible from a feasible solution s by improving the objective function value without losing feasibility. Provided with an input feasible solution s, DRF(s) procedure finds the subset U of markets belonging to s which are singularly not necessary for feasibility, i.e. each of these markets can be eliminated without affecting products demand and the budget threshold B satisfaction. Then, if U ∅, the algorithm evaluates for each market i ∈ U the routing cost decrease produced by its elimination from the current solution s and drops the market which entails the largest traveling cost reduction. Before recomputing set U, the routing cost of the new solution is possibly improved by applying GENIUS algorithm. The procedure stops as far as U = ∅. Notice that finding set U for a solution s requires a computational time complexity of O(|V s K| log(|Mk ∩ V s |)). 2.1. Algorithm EJEct and MOve (EJEMO) The first algorithm EJEMO (EJEct and MOve) consists in an enhanced local search procedure where a neighbor solution is obtained from a current one by removing a chain of l markets and inserting as many markets as required to restore feasibility. The procedure is more than a simple local search since parameter l and its maximum value lmax are changed according to the status of the search. The detailed description of the algorithm is provided in Fig. 1. Procedure ChainEjection(s, l), that allows to explore the neighborhood of a current solution s, consists of two main steps. In the first step the length of the current cycle s is reduced by removing the chain of l markets that yields the largest reduction in the traveling costs. At this aim the procedure looks for the set P(s) of all possible

chains with l markets belonging to the current solution s and which do not contain markets in M∗ or markets which are tabu. If set P(s) is void, an aspiration rule is applied which allows to consider tabu markets for chains construction. For each chain p ∈ P(s), the procedure evaluates its routing cost reduction as the difference between the traveling costs of the removed chain and that of directly joining the chain extreme points. The chain yielding the largest reduction is removed first. Given s, the number of possible l-chains is O(|V s |). The second step controls if the new cycle is feasible. If this is not the case solution feasibility is restored by consecutively applying, when necessary, GetAll(s) and BudgetRestore(s) procedures. Inserted markets (set S) become tabu for a maximum of rmax iterations during which they cannot be removed from the current solution. If the new solution s has a better value than s, the parameter l (initially set to a maximum value lmax ) remains unchanged (internal while loop) and the algorithm moves from s to s , otherwise l is reduced. As far as l reaches a minimum value lmin (for loop), the algorithm tries and improves the current solution s with procedure DRF(s). Then the search is diversified by generating a new solution s through the RandomShaking(s , h) routine. Such procedure randomly inserts in the current solution s as many markets as necessary to increase its traveling cost by at least a predefined percentage h (initially set to a maximum value hmax ). The procedure aims at possibly escaping from a local minimum by leading the search towards unexplored regions of the solution space. After the insertion of the new markets, the heuristic starts again with l set to a new initial value lmax computed as a fixed percentage of the markets making part of the current solution, whereas the diversification parameter h is reduced by a constant percentage h0 . The algorithm's stopping rule is given by h = hmin (main while loop). EJEMO is an efficient procedure depending on two main parameters: the maximum length of the chain to be removed (parameter lmax ) and the parameter controlling diversification and imposing a maximum number of restarts (parameter hmax ). In some preliminary experiments, we have analyzed the variant of the algorithm in which parameter l is initialized to its minimum value lmin and then increased up to lmax . Since this version has provided a worse performance we have decided to abandon it. We have also tested the introduction of an improvement routine for the traveling costs inside the Chain-Ejection(s, l) procedure. The choice has come out to be, on average, effective but the improvement has not been so significant to justify the additional computational time. We have thus decided to avoid it and to introduce an improving heuristic for the traveling cost only in the DRF(s) procedure just before the perturbation scheme (see Fig. 1). Finally, the idea of removing chains of consecutive markets is not new. We have decided to use it as providing extremely good results when applied to the TPP (see local search algorithm proposed in [12]). The main differences of EJEMO with respect to such heuristic are given by the idea of changing parameter lmax during the search and by the introduction of a dynamic tabu structure: markets inserted into a solution become tabu with a tabu tenure which is reduced by one at each iteration. As soon as the tabu tenure is zero it is reset at its maximum value rmax . The introduction of this tabu structure avoids cycling and facilitates to escape from local optimal solutions. 2.2. A VNS scheme In a general VNS scheme a sequence q = 1, . . . , qmax of neighborhoods Nq (·) is generated, each one characterized by a different value of a parameter q (also called the radius). Given a feasible solution s, if intensification of the search is required, the radius q is initialized to its minimum value so that the closest neighborhood is selected first. Otherwise q is set to its maximum value and then decreased

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

2267

many external markets as necessary to restore the feasibility by applying GetAll() and BudgetRestore() procedures. Moreover, with respect to EJEMO, no tabu structure is implemented (i.e. rmax = 0). In fact, some preliminary experiments have put in evidence how, in most of the solved instances, the tabu structure is made useless by the VNS scheme. As soon as an improving solution is obtained the evaluation is interrupted (first improvement rule), the algorithm moves to the new solution, whereas the parameter l is set to a new maximum value lmax (computed as a fixed percentage of the number of markets in the current solution) and the local search is restarted. If no improvement is achieved after evaluating all feasible solutions obtainable by the elimination of l consecutive markets, the parameter l is decreased by one unit and the procedure is repeated. When l = lmin the local search stops providing the best solution s found and the algorithm tries to improve its value by applying the DRF(s ) routine. If this solution is better than the incumbent one (sVNS ), the algorithm moves to this new solution (flag move is set to true) and the MVNS procedure restarts with q = qmin . Otherwise, if t  tmax , a new random solution s in the neighborhood Nq (sVNS ) is generated (multi-starts). As far as tmax random solutions (restarts) have been generated in the same neighborhood and the corresponding applications of the local search have not improved the incumbent solution value the parameter q is increased by a value qstep and the algorithm moves to a new neighborhood. The MVNS stops when q takes its maximum value qmax . Parameter q controls search diversification. When no better solutions can be found, q is increased allowing the search to explore farthest and possibly more promising regions of the solution space. The proposed VNS and its multi-start variant are not as efficient as EJEMO algorithm but, on average, provide more effective solutions. 2.3. Problem lower bound Fig. 2. Pseudo-code of the multi-start variable neighborhood search algorithm.

during the search. The VNS scheme was first proposed by Mladenovic and Hansen in [14] and has been effectively applied to many routing problems. In a pure VNS scheme the value of the parameter q is changed as far as no solution improvement can be obtained. In this paper we also consider a simple variant, called Multi-start VNS (MVNS), where the algorithm generates tmax random solutions in each neighborhood Nq (·) before the value of q is changed and a new neighborhood is taken into account. More precisely, the implemented VNS and its multi-start variant analyze increasing neighborhoods and use, as local search scheme to explore a neighborhood, a modification of the EJEMO algorithm. Since EJEMO is more than a simple local search algorithm the resulting VNS scheme has embedded a meta-heuristic approach as improving routine. Fig. 2 shows the pseudo-code of the MVNS approach. The basic VNS scheme can be directly obtained by setting tmax = 1. The procedure MVNS starts by setting the parameter q to its minimum value qmin . Given q and a solution sVNS a new solution s , belonging to the neighborhood Nq (sVNS ), is generated from sVNS by inserting q external random markets. Markets insertion is made according to the cheapest insertion rule (procedure RandomIns(sVNS , q)). Then procedure Improvement(s , l) is called. As procedure ChainEjection() (cf. EJEMO algorithm), this routine is based on the drop of l consecutive markets. However, while in the first case, the chain to be ejected is selected a priori as the one which yields the highest travel costs reduction, in Improvement() the procedure iteratively evaluates each feasible solution which can be obtained by removing l consecutive markets and by inserting as

Since no integer optimal solution for the capacitated version of the problem is available in the literature, we have removed connectivity constraints (3) from the original formulation (TPP-B) and solved the following relaxed integer problem: (RTPP-B)

min

vR :=



ce xe

(11)

e∈E

subject to: (2), (4)–(10)

(12)

Avoiding connectivity constraints usually results in a poor relaxation. Although such constraints are in exponential number, they can be added by a separation problem solvable in polynomial time. Let (x¯ , z¯ , y¯ ) be the optimal solution of problem (RTPP-B), to determine the most violated connectivity constraint a minimum-cut problem is solved on the capacitated undirected graph G = (V, E) with edge capacity equal to x¯ e , e ∈ E. In particular, for each vertex i, such that y¯ i > 0, the most violated constraint corresponds to the minimum capacity cut (S, V\S), with i ∈ S and vertex 0 in V\S, the value of which is less than 2. Since minimum cut problem is solvable in polynomial time O(m3 ) and the separation problem has to be solved at maximum m times such inequalities can be added to (RTPP-B) formulation in polynomial time O(m4 ). More precisely, we have solved problem (RTPP-B) with CPLEX 10 branch-and-cut routine and iteratively added to it the most violated inequalities belonging to the class of connectivity constraints up to a maximum running time is elapsed. If, within such time, all possible inequalities are added, the final integer solution is the optimal one, otherwise the best bound found by the branch-and-cut routine is provided as lower bound. The

2268

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

algorithm by Padberg-Rinaldi [17] has been used to solve the mincut problem. 3. Computational results 3.1. Testing environment Since no instances for the TPP with budget constraint are available in the literature, we have tested our algorithms on the Euclidean instances proposed for the TPP by Laporte et al. [10]. These instances can be downloaded at the web site http:// webpages.ull.es/users/jjsalaza. More precisely, we have used the instances of the Class 3 for the uncapacitated version of the problem and the instances of the Class 4 for the capacitated one. In both classes, markets are uniformly distributed with integer coordinates in a [0, 1000] × [0, 1000] square. Routing costs are represented by Euclidean distances. Each product k is associated with |Mk | randomly selected markets, where |Mk | is an integer number uniformly distributed in [1, |V| − 1]. Product prices are generated in the interval [1, 500] according to a discrete uniform distribution. Moreover, in the capacitated instances, for each product k and each market i, integer quantities qki are randomly generated in [1, 15], whereas product  demands dk are computed as  maxi∈M {qki } + (1 − ) i∈M qki with k

k

0 <  < 1. For the uncapacitated case we have set the budget threshold B equal to the optimal purchasing costs provided by Riera-Ledesma and Salazar-González [11] when solving the bi-objective formulation of the uncapacitated TPP on the instances of the Class 3. They have implemented an iterative procedure to find the efficient solutions set of the bi-criteria problem in the objective space. The first step of this general method determines, for each instance, two initial efficient points by optimizing hierarchically the two objective functions represented by the purchasing and the routing costs. More precisely, an efficient point is computed by finding the minimum purchasing cost over all markets and then by setting this value as a constraint threshold in a problem which minimizes the routing cost. The other one is obtained by first solving the problem that looks for the cycle at minimum cost while purchasing all products and then by finding the minimum purchasing cost over the selected markets. For each efficient point obtained by solving a given instance of the bi-objective problem we have generated a corresponding instance for the UTTPB by setting its budget threshold B equal to the optimal purchasing cost of the efficient point. This means that for each instance of Class 3 solved in [11], two different instances for the UTPP-B are created. Notice that while the instance corresponding to the first efficient point will have a tight budget threshold (the optimal purchasing cost corresponds to the minimum one), this is not the case for that corresponding to the second efficient point. Moreover, it is immediate to see that the optimal solutions for these UTPP-B instances correspond to the optimal traveling costs of the corresponding efficient points. The instances corresponding to the first efficient point and characterized by a budget equal to the minimum purchasing cost are identified as tight budget instances, whereas those associated to the second efficient point as untight budget instances. In the following we analyze the computational results for the instances of Class 3 with m = 50 and 100 markets and n = 50, 100, 150 and 200 products, for which the authors in [11] provide optimal solutions. Since for each combination of m and n five instances are generated, this means that we have generated and solved 80 uncapacitated instances, 40 out of which are tight budget instances and 40 untight budget ones. For the capacitated case, the optimal solutions for the TPP instances in [10] could be used as optimal solutions for our problem by setting the threshold B equal to the optimal purchasing cost. Unfortunately such information is not provided by the authors where the objective function value is given as the cumulative value of routing

and purchasing costs. Thus, for each capacitated instance of Class 4 we have introduced three different values of the budget threshold B obtained by increasing the corresponding minimum purchasing cost over all markets by a given percentage  set equal to 0, 50% and 100%. In this case tight budget instances are defined for  =0, whereas budget untightness increases with  =50% and 100%. We generate the budget threshold by using the minimum purchasing cost since such value can be easily computed. On the contrary, the maximum purchasing cost requires the solution of a minmax problem looking for that subset of markets which satisfies the products demand and for which the minimum purchasing cost is maximum. We have solved instances of Class 4 with m = 50 and 100 markets, n = 50, 100 and 150 products and  = 0.8, 0.9 and 0.99. This means five instances for each combination of n, m and  each one considered with three different values of the budget threshold for 270 instances altogether. For such instances we have compared our algorithms performance to the solutions value found by solving problem (RTPP-B) (cf. Section 2.3) with the addition of connectivity inequalities within a computation time threshold of 2 h. In order to set parameter values for each algorithm we have conducted some preliminary experiments. Details on this parameters tuning for the case of the uncapacitated instances can be found in [18]. In particular, for such instances the following choices have been implemented. For heuristic EJEMO we have set the parameter hmax = 35%. This means that as many markets as necessary to increase the traveling costs by at least 35% will be added to the current solution when h = hmax . Parameters hmin and h0 have been set to 1%. The decision has been taken after trying alternative values for the maximum percentage hmax equal to 20%, 25%, 40%, and 50%, respectively. The maximum value for the tabu tenure rmax has been set to 4, parameter lmin has been set to 1 and lmax equal to 10% of the number of the markets making part of the current solution (values are truncated to integer numbers). For algorithm VNS we have easily set qmax = 20 and qmin = qstep = 1 as the best parameter values, while some testing has been required in order to choose the initial value for parameter l which is resulted to be more crucial. At this aim we have solved several instances with different values of parameter lmax and found out that the best algorithm performance can be obtained by setting lmax equal to 5% of the number of markets in the current solution. Parameter lmin has been always set to 1. The same parameters setting has been used for MVNS algorithm where the number of restarts tmax has been set to 5. For the capacitated instances the tuning of the parameters for EJEMO has provided almost the same results found for the uncapacitated case, so no relevant changes have been implemented. For VNS and MVNS approaches the parameters tuning has required some additional analysis as strictly depending on the combination of  and budget B values. In particular, for tight budget instances (i.e.  = 0) the best trade-off between quality of the solution and computational time can be obtained by setting lmax equal to the 5% of the number of markets in the current solution. For the untight budget instances (i.e.  = 50% and 100%) lmax is set equal to 10% and to 15% of the number of markets in the current solution for  = 0.8, 0.9 and 0.99, respectively. For all instances, lmin = 1. Moreover, qmax has been set according to the size of the instance equal to 15 m for the VNS, while for the MVNS the best trade-off between computational time and quality of the solutions has led to set qmax equal to m/10. In both approaches qstep = 2 and qmin = 1. The number tmax of random solutions generated for each neighborhood has been set to 5. Algorithms have been coded in ANSI C and run on a PC Pentium with 3.5 GHz and 2 GB of RAM. Problem RTPP-B has been solved by using Concert Technology 22 and CPLEX 10.0. Finally, given the random nature of the proposed solution algorithms, each of them has been run 5 times (5 trials) over each instance (either capacitated or uncapacitated).

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

3.2. Uncapacitated TPP-B: computational results The following two tables provide the main findings for the UTPPB out of the 40 tight budget instances (Table 1) and out of the 40 untight budget ones (Table 2), respectively. Each table is divided into three parts, the first one provides the computational results for heuristic EJEMO, whereas the second and the third one those for the VNS procedure and its multi-start variant (MVNS). The column headings in each table have the following meaning: m and n identifies the size of the instances solved; avg, min and max refers to the average, the minimum and the maximum percentage error of the heuristic solution value from the optimal one out of the five different instances tested and the five trials; column referred as # returns how many times the heuristic reaches the optimal value over the five instances tested (considering the five trials); finally, column time provides the corresponding average computational time in seconds. The last line of each table provides the average values out of the 40 solved instances but for column # which indicates the total number of instances solved to optimality by each solution algorithm. By comparing the two tables it is evident how tight budget instances are the most difficult ones. In such instances EJEMO algorithm behaves rather poorly providing an average error equal to 6.11% out of the 40 solved instances (given the five trials). Moreover, it is never able to find the optimal solution value. Tight budget instances are characterized by a limited number of feasible solutions. While analyzing the results, we have noticed that, in some cases, the bad performances of the algorithm on such instances are due to ChainEjection() procedure which fails in selecting the right chain of markets to be removed. The simple selection rule used by such procedure behaves acceptably well only when the number of good quality solutions is quite large as for untight budget instances. However, there are several instances, especially for m=100, for which the errors are also due to the bad route construction. In all such cases GENIUS algorithm has failed to improve the routing costs of the solution received as input. The VNS procedure has come out to be quite effective yielding an average error equal to 0.52% and finding the optimal solution value in 32 out of the 40 tight budget instances solved (considering the five trials). In the MVNS the number of optimally solved instances increases to 34 even if such better performance is obtained with an increase in the computation time more than proportional. From the efficiency point of view, the computational time required by EJEMO is almost constant for all instances (its average value out of all instances is equal to 5.75 s), whereas the average computational time for the VNS algorithm and the multi-start variant grows quickly with the problem size reaching for instances with 100 markets and 200 products average values of 277 and 3363.1 s, respectively. All the heuristics performance improves when untight budget instances are taken into account. Table 2 shows that EJEMO provides an average percentage error out of all instances equal to 0.37% with a maximum error never larger than 3.84%. Moreover, given the five trials, it is able to find the optimal solution values in 28 out of 40 solved instances. The results obtained by the VNS approach are even better: this procedure yields an average percentage error equal to 0.12% out of the 40 solved instances and gets the optimal solution value (given the five trials) in all but three instances. The average computational times are close to 7 s for both the heuristics. The MVNS algorithm finds the optimal solutions in all trials and all instances solved with m = 50 (notice that this is not the case for the VNS approach), while it increases by 2 the number of optimally solved instances with respect to VNS. Its computational time is definitely larger than those of the other heuristics but still acceptable being no more than few minutes. To conclude, if we consider tight budget instances the ranking of the proposed algorithms in terms of solution values suggests VNS as

2269

a first choice followed by EJEMO. Nevertheless, if running time is a crucial issue EJEMO may be used at a price of a lower performance. On the contrary for untight budget instances procedure VNS is the dominant one providing solutions which are more effective and as efficient as those yielded by EJEMO. 3.3. Capacitated TPP-B: computational results Table 3 contains average results grouped according to the value of n, B and , for all the heuristics out of the 270 instances tested for the capacitated problem. The table is organized as followed. The first and the second columns represent the pair of parameters identifying the instances with respect to which the average results have been obtained. For instance, the first line shows the average results out of all possible values of  and of the budget B obtained when considering instances with 50 markets and 50 products. Notice that each entry in the table provides an average value computed out of 45 instances. The third column of the table shows the average computation time required by CPLEX 10.0 to solve problem (RTPP-B) with the addition of connectivity constraints when a maximum computation time of 2 h is assigned to the exact routine. Each entry in the fourth column (symbol #) provides the number of times out of the 45 solved instances, CPLEX has found the integer optimal solution within such time threshold. The remaining columns are divided into three parts, one for each heuristic. Each part provides the average and the minimum percentage errors with respect to the best bound found by CPLEX within 2 h (gap), the average and minimum computation time (in seconds) required by the algorithm (time) and the total number of instances in which each algorithm finds the optimal solution (column #). By comparing the results, it is evident how, on average, EJEMO behaves quite poorly but is extremely efficient never requiring an average time larger than 15 s. The results obtained by the VNS and the MVNS are comparable in terms of percentage error: out of all 270 solved instances the average gap is equal to 0.93% and to 0.84% for the VNS and the MVNS, respectively. However, while VNS has an average computation time around 1 min for instances with 50 markets and equal to 946 s out of all 270 instances, the same values for the MVNS increase to 2 min and 1600 s, respectively. Since in our problem the main parameter is the budget threshold, in the following we further analyze the section of Table 3 corresponding to the results grouped according to the budget parameter  (i.e. lines 7–12). Tables 4 and 5 provide such details. In particular, Table 4 refers to the instances with 50 markets and Table 5 to those with 100 markets, respectively. Both tables are organized as followed. The first three columns identify the five instances with respect to which the average results are computed. For instance, the first line refers to the instances where the value of the budget B is set to its minimum ( = 0), the number of products n is equal to 50 and the value of  is 0.8. The meaning of the remaining columns is as in Table 3. Since for m = 50 CPLEX has been always able to find the optimal solution within the 2 h of computation time threshold (cf. Section 2.3), some interesting conclusions can be drawn by analyzing the percentage errors reported in Table 4. Most of the instances with  = 0.8, 0.9 and  = 0 result to be TSP instances. In such case our algorithms are ineffective and the objective values provided correspond to the initial solution values found by applying GENIUS algorithm. On the contrary if we take into account instances with larger values of  and  (larger budgets) the resulting problems are real TPP-B instances. In general, our VNS algorithm is able to find the optimal solution (128 instances out of 135) in a low computation time (on average less than 1 min). The quality of the solutions found by MVNS is even better (130 instances solved to optimality). Nevertheless, being the results for the VNS almost optimal, MVNS improvement does

2270

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

Table 1 Computational results for the UTPP-B: tight budget instances. m

n

EJEMO

VNS min

max

MVNS min

max

50

50 100 150 200

3.54 6.93 6.09 5.43

2.49 4.07 4.48 4.05

5.07 8.97 8.60 7.31

0 0 0 0

1.9 2.8 6.0 9.9

0.47 0.15 0.79 0.34

0.12 0.00 0.41 0.00

1.32 0.42 1.17 0.94

4 5 3 5

3.5 9.1 17.2 27.4

0.12 0.02 0.91 0.00

0.00 0.00 0.91 0.00

0.28 0.06 0.93 0.00

5 5 2 5

134.2 257.0 335.8 442.8

100

50 100 150 200

3.78 6.09 9.25 7.76

1.76 3.33 4.78 4.76

6.00 8.87 12.63 10.03

0 0 0 0

13.8 3.9 4.8 3.0

0.21 0.54 1.08 0.55

0.01 0.01 0.12 0.06

0.59 2.48 3.38 0.91

4 4 3 4

29.9 157.4 288.3 277.0

0.02 0.01 0.30 0.33

0.00 0.00 0.03 0.31

0.02 0.03 0.69 0.35

5 5 4 3

318.8 968.4 2379.6 3363.1

6.11

3.72

8.44

0

5.7

0.52

0.09

1.40

32

101.2

0.21

0.16

0.30

34

1025.0

time

avg

min

max

#

time

avg

min

max

#

time

Average

avg

#

time

avg

min

max

#

time

avg

#

time

Table 2 Computational results for the UTPP-B: untight budget instances. m

n

EJEMO

VNS

avg

min

max

#

MVNS

50

50 100 150 200

0.00 0.53 0.11 0.05

0.00 0.37 0.11 0.02

0.00 0.84 0.11 0.10

5 4 3 4

1.26 1.61 13.92 14.43

0.01 0.19 0.09 0.09

0.00 0.00 0.00 0.00

0.00 0.32 0.24 0.31

5 5 5 5

0.54 1.06 1.34 2.41

0.00 0.00 0.00 0.00

0.00 0.00 0.00 0.00

0.00 0.00 0.00 0.00

5 5 5 5

43.60 69.62 89.47 100.94

100

50 100 150 200

0.37 0.83 0.14 0.92

0.24 0.23 0.00 0.87

0.56 1.16 0.47 0.96

4 1 5 2

8.29 2.32 2.36 3.20

0.07 0.08 0.00 0.44

0.00 0.05 0.00 0.02

0.11 0.19 0.01 0.71

5 2 5 5

4.08 13.83 14.86 18.26

0.00 0.04 0.00 0.02

0.00 0.04 0.00 0.00

0.00 0.05 0.00 0.05

5 4 5 5

99.07 209.03 164.80 247.16

0.37

0.23

0.53

28

5.92

0.12

0.01

0.24

37

7.05

0.01

0.00

0.01

39

127.96

Average

Table 3 Average values for the CTTP-B instances. m

CPLEX time

50

100

products 50 100 150 50 100 150

EJEMO #

VNS

gap

time

avg

min

avg

# min

MVNS

gap

time

avg

min

avg

# min

gap

time

avg

min

avg

# min

94.66 87.81 163.59

45 45 45

1.63 1.64 4.73

0.85 0.99 4.09

2.01 2.74 3.61

1.89 2.68 3.47

19 15 13

0.13 0.10 0.13

0.02 0.05 0.02

30.06 55.86 73.14

24.48 45.16 57.47

43 41 44

0.06 0.06 0.08

0.01 0.04 0.02

60.13 108.66 154.11

50.33 92.43 129.41

44 42 44

3802.78 4687.91 4775.60

34 25 19

5.10 4.68 5.44

3.17 3.28 4.13

5.49 9.27 14.39

5.28 8.92 13.96

1 1 0

1.63 1.48 2.08

1.24 1.31 1.86

852.79 1599.60 3063.91

627.41 1380.98 2179.83

20 18 10

1.36 1.43 2.05

1.20 1.29 1.88

2189.24 2404.03 4923.47

1819.25 2151.55 3933.86

23 18 11

budget () 50

0 50 100

12.76 133.34 190.13

45 45 45

2.56 1.03 0.98

2.00 0.28 0.21

2.82 2.66 2.88

2.72 2.57 2.74

14 11 22

0.07 0.15 0.14

0.03 0.03 0.03

22.10 76.33 58.51

18.69 60.41 46.50

43 42 43

0.06 0.08 0.06

0.03 0.01 0.02

56.23 156.32 107.00

50.29 132.87 86.74

43 43 44

100

0 50 100

1529.44 5462.62 5855.13

41 20 17

5.42 5.68 5.91

4.55 3.66 4.12

8.18 9.17 11.81

8.03 8.83 11.32

0 0 2

0.89 1.67 2.63

0.78 1.39 2.24

230.45 2344.03 2941.81

160.81 1751.37 2276.04

17 15 16

0.84 1.59 2.40

0.76 1.35 2.26

820.01 3932.89 4531.23

574.50 3337.61 3800.44

19 17 16

lambda 50

0.8 0.9 0.99

123.55 134.20 78.48

45 45 45

4.54 1.17 2.30

4.21 0.47 1.25

3.02 2.81 2.53

2.91 2.70 2.42

12 12 23

0.16 0.07 0.14

0.05 0.01 0.03

62.24 59.81 34.89

50.17 48.14 27.29

42 44 42

0.10 0.03 0.07

0.04 0.01 0.02

137.51 121.44 60.60

119.77 101.00 49.12

44 44 42

100

0.8 0.9 0.99

3722.58 4092.22 5032.40

30 29 19

2.63 5.25 9.13

1.67 3.86 6.81

13.37 10.52 5.26

13.00 10.14 5.04

1 1 0

0.91 1.26 3.03

0.78 0.91 2.72

2993.90 2078.07 444.32

2224.54 1612.29 351.38

14 17 17

0.87 1.13 2.83

0.76 0.88 2.72

3402.50 4096.74 1784.88

2788.73 3451.50 1472.31

17 16 19

2197.24

213

4.17

3.04

6.25

6.03

49

0.93

0.75

945.54

718.97

176

0.84

0.74

1600.61

1330.41

182

Average

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

2271

Table 4 Computational results for the CTTP-B instances: average values with 50 markets.





n

CPLEX time

0

50

100

150

50

50

100

150

100

50

100

150

EJEMO #

VNS

gap

time

#

gap

MVNS time min

gap

time

#

avg

min

avg

min

avg

min

avg

min

avg

min

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

6.16 8.03 20.76 12.22 11.03 11.00 8.43 10.27 26.94 12.76

5 5 5 5 5 5 5 5 5 45

1.92 3.35 2.31 3.10 0.09 5.73 0.24 0.01 6.25 2.56

1.92 2.93 1.50 3.10 0.09 3.98 0.24 0.01 4.20 2.00

2.34 2.33 1.86 2.92 2.72 2.61 3.63 3.51 3.47 2.82

2.22 2.16 1.74 2.91 2.70 2.48 3.55 3.37 3.35 2.72

1 2 3 1 1 3 0 2 1 14

0.06 0.07 0.25 0.14 0.00 0.05 0.01 0.01 0.08 0.07

0.06 0.07 0.00 0.14 0.00 0.00 0.01 0.01 0.00 0.03

10.91 23.20 28.14 26.70 27.89 26.07 0.93 0.98 54.12 22.10

10.90 21.68 22.07 25.47 26.74 19.21 0.93 0.98 40.25 18.69

5 4 5 4 5 5 5 5 5 43

0.06 0.07 0.15 0.14 0.00 0.05 0.01 0.01 0.03 0.06

0.06 0.07 0.00 0.14 0.00 0.01 0.01 0.01 0.00 0.03

29.92 58.17 65.64 69.89 74.54 62.10 0.93 0.98 143.87 56.23

29.86 55.39 53.36 66.80 71.29 50.58 0.93 0.98 123.44 50.29

5 4 5 5 5 4 5 5 5 43

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

84.39 116.33 72.61 109.87 105.67 103.27 340.76 174.35 92.80 133.34

5 5 5 5 5 5 5 5 5 45

0.70 1.54 1.34 1.11 1.35 0.93 1.27 0.13 0.87 1.03

0.17 0.17 0.41 0.41 0.02 0.42 0.67 0.11 0.14 0.28

2.12 1.99 1.61 2.96 2.57 2.29 3.69 3.51 3.23 2.66

1.98 1.87 1.55 2.89 2.53 2.23 3.54 3.39 3.12 2.57

0 1 1 2 0 2 0 2 3 11

0.00 0.11 0.22 0.17 0.03 0.08 0.64 0.04 0.08 0.15

0.00 0.00 0.08 0.05 0.00 0.00 0.17 0.00 0.00 0.03

58.03 24.53 26.77 94.28 99.71 24.60 158.09 141.26 59.70 76.33

46.94 20.48 23.32 75.95 79.38 18.21 122.63 108.13 48.64 60.41

5 5 4 4 5 5 5 5 4 42

0.00 0.05 0.11 0.02 0.00 0.02 0.47 0.01 0.05 0.08

0.00 0.00 0.00 0.00 0.00 0.00 0.13 0.00 0.00 0.01

65.87 55.67 27.67 259.54 115.57 56.59 388.57 369.98 67.41 156.32

55.93 50.95 22.53 213.13 96.38 48.67 351.17 307.12 49.98 132.87

5 5 5 4 5 5 5 5 4 43

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

136.63 137.72 180.80 136.96 198.07 102.19 276.56 446.28 95.95 190.13

5 5 5 5 5 5 5 5 5 45

0.49 1.31 1.71 0.72 0.97 0.79 0.29 1.77 0.77 0.98

0.06 0.19 0.32 0.31 0.35 0.20 0.10 0.35 0.03 0.21

2.20 2.03 1.59 3.31 2.79 2.53 4.05 3.81 3.57 2.88

2.04 1.89 1.52 3.20 2.70 2.44 3.83 3.68 3.37 2.74

4 3 4 4 0 2 0 1 4 22

0.25 0.08 0.18 0.08 0.06 0.25 0.08 0.20 0.04 0.14

0.00 0.00 0.00 0.06 0.00 0.18 0.00 0.00 0.00 0.03

50.95 18.51 10.36 80.37 80.42 42.74 79.88 121.85 41.49 58.51

38.67 15.49 7.23 63.04 63.04 35.34 66.98 97.36 31.32 46.50

5 5 5 4 5 4 5 5 5 43

0.09 0.04 0.00 0.05 0.04 0.19 0.08 0.03 0.01 0.06

0.00 0.00 0.00 0.00 0.00 0.19 0.00 0.00 0.00 0.02

137.59 45.01 25.41 192.41 97.97 49.34 92.89 275.05 47.33 107.00

108.14 39.03 17.34 171.35 74.53 39.13 80.66 213.35 37.10 86.74

5 5 5 5 5 4 5 5 5 44

112.08

135

1.52

0.83

2.79

2.68

47

0.12

0.03

52.31

41.87

128

0.07

0.02

106.51

89.97

130

Average

avg

#

Table 5 Computational results for the CTTP-B instances: average values with 100 markets.



0



n

50

100

150

50

50

100

150

100

50

100

150

Average

CPLEX

EJEMO

VNS

time

#

gap

time

avg

min

avg

min

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

31.00 737.60 2240.20 752.60 1112.20 4135.80 1028.20 917.00 3170.20 1569.42

5 5 4 5 5 3 5 5 4 41

1.63 10.87 9.18 1.38 1.62 10.28 1.16 3.62 9.06 5.42

1.63 10.37 6.29 1.38 1.62 7.74 1.16 3.62 7.16 4.55

7.24 8.14 2.78 8.22 8.42 5.43 12.75 11.68 8.98 8.18

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

3654.80 4430.80 4931.80 5880.00 6964.20 5215.40 6365.80 5286.00 6434.80 5462.62

4 4 3 3 1 2 1 2 0 20

3.56 8.15 9.63 2.89 3.18 6.70 3.11 4.85 9.05 5.68

2.11 3.92 6.52 0.83 2.05 5.05 2.05 2.85 7.57 3.66

0.8 0.9 0.99 0.8 0.9 0.99 0.8 0.9 0.99 avg

3516.60 4874.80 6395.40 6055.20 6150.60 5925.20 6219.00 6356.80 7202.60 5855.13

5 3 1 1 3 2 1 1 0 17

2.74 4.99 11.26 3.15 4.68 8.24 4.02 5.32 8.74 5.91

4295.73

78

5.67

VNS multi-start

#

gap

time

#

gap

avg

min

7.77 7.58 3.79 7.76 7.72 4.76 7.74 7.72 5.07 6.66

0 0 0 0 0 0 0 0 0 0

0.71 1.02 0.52 1.03 0.76 1.07 1.75 0.61 0.57 0.89

0.71 0.77 0.23 1.03 0.61 0.94 1.75 0.48 0.50 0.78

33.63 249.25 272.08 4.40 150.82 221.60 4.58 694.59 443.11 230.45

avg

min

23.30 184.75 181.48 4.40 121.12 221.60 4.58 406.56 299.51 160.81

1 1 2 1 2 3 1 2 4 17

0.72 0.91 0.23 1.03 0.76 1.04 1.75 0.61 0.50 0.84

0.72 0.57 0.20 1.03 0.61 0.94 1.75 0.48 0.50 0.76

95.70 774.84 637.92 4.40 359.99 433.15 4.58 1921.33 3148.22 820.01

90.21 594.77 515.73 4.40 323.50 433.15 4.58 885.45 2318.68 574.50

1 1 4 1 2 3 1 2 4 19

7.21 4.85 2.55 12.54 9.50 4.63 18.67 14.62 7.93 9.17

4.87 5.27 3.19 4.91 5.73 3.56 4.94 5.67 3.98 4.68

0 0 0 0 0 0 0 0 0 0

0.49 0.92 2.77 0.19 1.08 2.95 0.99 1.19 4.43 1.67

0.29 0.31 2.58 0.13 0.73 2.76 0.71 0.73 4.25 1.39

2640.65 786.19 208.22 2492.55 3486.14 566.46 7974.07 2389.30 552.68 2344.03

1903.26 558.58 155.16 2044.49 3080.19 566.46 5418.94 1688.81 346.41 1751.37

2 3 3 3 1 2 0 1 0 15

0.26 0.86 2.60 0.20 0.76 3.22 0.98 1.03 4.42 1.59

0.11 0.41 2.58 0.13 0.54 2.76 0.60 0.80 4.20 1.35

7287.56 967.18 477.27 948.54 10481.43 769.28 2750.24 8529.99 3184.50 3932.89

6349.39 653.67 353.89 739.37 9600.85 769.28 2096.18 6921.20 2554.62 3337.61

4 2 3 3 1 2 1 1 0 17

1.82 3.86 7.85 1.79 2.62 6.43 2.25 3.79 6.69 4.12

8.81 5.37 2.45 17.62 12.47 4.61 27.30 19.64 7.98 11.81

3.78 4.12 2.96 3.83 4.46 3.22 3.88 4.33 3.47 3.78

1 0 0 0 1 0 0 0 0 2

0.30 1.85 6.11 1.18 1.29 3.74 1.51 2.60 5.11 2.63

0.01 1.23 5.01 1.02 1.16 3.38 1.33 2.17 4.83 2.24

2330.11 726.39 428.58 4246.32 2803.30 424.78 7218.85 7416.63 881.38 2941.81

1733.55 574.12 332.53 3894.18 2071.58 424.78 4994.17 5824.94 634.50 2276.04

4 3 1 1 3 2 1 1 0 16

0.12 1.43 5.09 1.23 1.23 3.38 1.57 2.57 5.02 2.40

0.01 1.15 5.01 1.05 1.15 3.38 1.44 2.19 4.93 2.26

5209.35 1101.29 1058.53 3808.79 3870.16 960.55 10513.38 8864.50 5394.51 4531.23

4348.35 881.30 856.91 3147.99 3384.85 960.55 8318.10 7817.92 4488.00 3800.44

4 3 1 1 3 2 1 1 0 16

4.11

9.72

5.04

2

1.73

1.47

1838.76

1396.07

48

1.61

1.45

3094.71

2570.85

52

avg

min

time avg

# min

2272

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

1.00% 0.90% 0.80%

min

VNS

0.70%

avg

MVNS

0.60%

max

0.50% 0.40% 0.30% 0.20% 0.10% 0.00%

0

50

100

0

0.8

50

100

0

50

0.9

100

0.99

Fig. 3. Percentage errors with respect to CPLEX (m = 50, n = 50).

1.00% 0.90% 0.80%

min

VNS

0.70%

avg

MVNS

0.60%

max

0.50% 0.40% 0.30% 0.20% 0.10% 0.00%

0

50

100

0

0.8

50

100

0

0.9

50

100

0.99

Fig. 4. Percentage errors with respect to CPLEX (m = 50, n = 100).

1.00% 0.90% 0.80%

min

0.70%

avg

0.60%

max

0.50%

VNS

0.40%

MVNS

0.30% 0.20% 0.10% 0.00%

0

50

100

0

0.8

50

100

0

0.9

50

100

0.99

Fig. 5. Percentage errors with respect to CPLEX (m = 50, n = 150).

not seem to be consistent enough to justify the corresponding increase in the computation time. Figs. 3–5 help to synthesize the results reported in Table 4. They show the minimum, the average and the maximum percentage errors found by the VNS and MVNS when the number of products n is set equal to 50, 100 and 150, respectively. More precisely, each figure is

divided into three sections each one associated to a different value of the parameter . For each section the percentage errors found by the two procedures (a pair of columns) when the budget threshold is computed for  = 0%, 50% and 100% are shown. In each pair of columns the first one is associated to VNS heuristic (stripped column) the other one to the MVNS approach. The height of each column

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

2273

Table 6 Computational results for the CTTP-B instances: budget with  = 100%. n



2h

4h

CPLEX time

VNS #

CPLEX

gap

#

avg

min

time

VNS #

gap

#

avg

min

50

0.8 0.9 0.99

3516.60 4874.80 6395.40

5 3 1

0.30 1.85 6.11

0.01 1.23 5.01

4 3 1

3516.60 7754.80 12155.40

5 3 1

0.30 1.31 3.90

0.01 0.68 2.83

4 3 1

100

0.8 0.9 0.99

6055.20 6150.60 5925.20

1 3 2

1.18 1.29 3.74

1.02 1.16 3.38

1 3 2

7493.60 9028.20 10242.60

4 3 2

0.47 0.94 2.62

0.28 0.80 2.26

4 3 2

150

0.8 0.9 0.99

6219.00 6356.80 7202.60

1 1 0

1.51 2.60 5.11

1.33 2.17 4.83

1 1 0

10535.60 11963.40 11521.40

2 2 2

0.86 1.83 3.58

0.68 1.40 3.31

2 2 2

5855.13

17

2.63

2.24

16

9356.84

24

1.76

1.36

23

Average

represents the maximum percentage error (upper edge of the light gray part) whereas the upper edge of the black part and that of the gray one provide the minimum and the average percentage errors computed out of the five instances identified by the values of  and . For instance, in Fig. 3 with  = 0.9 and  = 50% the VNS and MVNS procedures have the same maximum error and the minimum one equal to zero, but the average one for the VNS is higher than that obtained by the MVNS. Each figure can be read in two ways. By comparing the percentage errors within a same section or by selecting a given level of the budget threshold and analyzing how errors change through the sections. In the first case, we can study how the budget threshold affects the solutions quality. In the second case the same analysis but for the  value can be provided. The analysis of the three figures allow to conclude that whereas VNS approach is more efficient, it is, however, more unstable providing average errors comparable to those found by the MVNS but larger maximum errors (larger deviation). For example, in the third section of Fig. 3, the maximum percentage error found by the VNS tends to grow when the budget threshold becomes untight (larger  values), whereas the corresponding average error tends to decrease. For the MVNS both values decrease. Notice again that, in the instances with  = 0 and  = 0.8 and 0.9 the two algorithms show the same behaviors being TSP instances. Table 5 provides the results for the instances with m = 100. In such case CPLEX 10.0 is able to find the optimal solution within 2 h of computation time in 41 out of the 45 instances with tight budget (=0). This can be explained by the fact that a tight budget is a helpful cut in the branch-and-cut approach. For larger budget thresholds the number of optimally solved instances drastically shrinks to 20 and 17 out of 45 for instances with  = 50% and 100%, respectively. In many instances of Table 5 the errors found by algorithms are quite high. For example, with  = 50% and n = 100 the percentage error found in some instances is close to 10%. Errors are even worse for instances with  = 100%. In order to evaluate how much of these errors can be explained by a poor value (far from the integer optimal value) of the bound provided by CPLEX, we have increased the time made available to the exact routine from 2 to 4 h. Results are shown in Table 6. It can be noticed that by increasing the computation time, the percentage errors obtained by the VNS decrease, on average, by one percentage point and the number of instances solved to optimality increases from 17 to 24. Similar conclusions can be drawn for the MVNS approach. Finally, we have selected a set of eight instances where VNS heuristic has found errors larger than 9% and tried to exactly solve them by further extending the computation time threshold made available to CPLEX up to 20 h. In only four instances CPLEX was able

to find the optimal solution and in all these cases we have proved that the solutions found by our algorithm were optimal. 4. Conclusions In this paper we have introduced a variant of the traveling purchaser problem (TPP) where the objective function looks for a minimum cost cycle visiting a subset of markets such that all the required products are bought and the total purchasing cost does not exceed a given budget threshold (TPP with budget constraint). To solve the problem we have proposed an efficient enhanced local-search procedure and a more effective variable neighborhood search approach also implemented in a multi-start variant. A large set of benchmark problems with both tight and untight budget constraints and for both the capacitated and the uncapacitated version of the problem has been introduced to test the algorithms. For many of these instances we have been able to provide optimal solution values. This experimental analysis provides the first computational results made available for the TPP with budget constraint. Acknowledgments The authors wish to thank anonymous reviewers for their useful comments which helped to improve the presentation of the work. References [1] Ramesh T. Traveling purchaser problem. Opsearch 1981;18:78–91. [2] Singh KN, van Oudheusden DL. A branch and bound algorithm for the traveling purchaser problem. European Journal of Operational Research 1997;97:571–9. [3] Ravi R, Salman S. Approximation algorithms for the traveling purchaser problem and its variants in network design. In: Proceedings of the 7th annual European symposium on algorithms. Lecture notes in computer science, vol. 1643. Berlin: Springer; 1999. p. 29–40. [4] Golden BL, Levy L, Dahl R. Two generalizations of the travelling salesman problem. Omega 1981;9:439–45. [5] Ong HL. Approximate algorithms for the traveling purchaser problem. Operations Research Letters 1982;1:201–5. [6] Pearn WL, Chien RC. Improved solutions for the traveling purchaser problem. Computers & Operations Research 1998;25:879–85. [7] Vo S. Dynamic tabu search strategies for the traveling purchaser problem. Annals of Operations Research 1996;63:253–75. [8] Boctor FF, Laporte G, Renaud J. Heuristics for the traveling purchaser problem. Computers & Operations Research 2003;30(4):491–504. [9] Mansini R, Pelizzari M, Saccomandi R. An effective tabu search algorithm for the capacitated traveling purchaser problem. Technical Report 10–49, Dipartimento di Elettronica per l'Automazione, University of Brescia; 2005. [10] Laporte G, Riera-Ledesma J, Salazar-González JJ. A branch-and-cut algorithm for the undirected traveling purchaser problem. Operations Research 2003;51(6):940–51.

2274

R. Mansini, B. Tocchella / Computers & Operations Research 36 (2009) 2263 -- 2274

[11] Riera-Ledesma J, Salazar-González JJ. The biobjective traveling purchaser problem. European Journal of Operational Research 2005;160:599–613. [12] Riera-Ledesma J, Salazar-González JJ. A heuristic approach for the traveling purchaser problem. European Journal of Operational Research 2005;162(1): 142–52. [13] Keller CP. Algorithms to solve the orienteering problem: a comparison. European Journal of Operational Research 1989;41:224–31. [14] Mladenovic N, Hansen P. Variable neighbourhood search. Computers & Operations Research 1997;24:1097–100. [15] Rosenkrantz DJ, Stearns RE, Lewis II PM. An analysis of several heuristics for the traveling salesman problem. SIAM Journal of Computing 1977;6:563–81.

[16] Gendreau M, Hertz A, Laporte G. New insertion and postoptimization procedures for traveling salesman problem. Operations Research 1992;40(6):1086–94. [17] Padberg M, Rinaldi G. An efficient algorithm for the minimum capacity cut problem. Mathematical Programming 1990;47:19–36. [18] Mansini R, Tocchella B. Local-search heuristics for a bounded version of the uncapacitated TPP. Technical report 01–51, Dipartimento di Elettronica per l'Automazione, University of Brescia; 2006 (to appear in Bertazzi L, Van Nunen J, Speranza MG, editors. Innovations in distribution logistics. Berlin: Springer).