- Email: [email protected]
Bounds for the single source modular capacitated plant location problem
Computers & Operations Research 33 (2006) 2991 – 3003 www.elsevier.com/locate/cor
Bounds for the single source modular capacitated plant location problem Isabel Correiaa,∗ , Maria Eugénia Captivob a Departamento de Matemática - C.M.A., Faculdade de Ciências e Tecnologia, Universidade Nova de Lisboa,
Quinta da Torre, 2829-516 Monte da Caparica, Portugal b Universidade de Lisboa, Faculdade de Ciências, Centro de Investigação Operacional Bloco C6, Piso 4,
Campo Grande, 1749-016 Lisboa, Portugal Available online 8 April 2005
Abstract In this paper, we propose a discrete location problem, which we call the Single Source Modular Capacitated Location Problem (SS-MCLP). The problem consists of finding the location and capacity of the facilities, to serve a set of customers at a minimum total cost. The demand of each customer must be satisfied by one facility only and the capacities of the open facilities must be chosen from a finite and discrete set of allowable capacities. Because the SS-MCLP is a difficult problem, a lagrangean heuristic, enhanced by tabu search or local search was developed in order to obtain good feasible solutions. When needed, the lower bounds are used in order to evaluate the quality of the feasible solutions. Our method was tested computationally on randomly generated test problems some of which are with large dimensions considering the literature related to this type of problem. The computational results obtained were compared with those provided by the commercial software Cplex. 䉷 2005 Elsevier Ltd. All rights reserved. Keywords: Capacitated location; Lagrangean heuristic; Tabu search
0. Introduction Location problems have been widely studied in the last years due to their application in many real situations. In fact, for many companies one of the key strategic decisions is to decide where to build ∗ Corresponding author. Tel.: +351 212948388; fax: +351 212948391.
E-mail address: [email protected] (I. Correia). 0305-0548/$ - see front matter 䉷 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.cor.2005.02.030
2992
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
facilities (e.g., warehouses, concentrators) and how to serve customers in order to optimize a given objective. These kinds of problems are not restricted to plant location alone, they also crop up in cases, such as design of telecommunication networks and distribution systems planning. In the SS-MCLP, a set of potential facility locations and a set of customers, each with a known demand, are given. In each potential location, there is a finite and discrete set of allowable capacities for the facilities and the demand of each customer must be satisfied by a single facility. The purpose of the SS-MCLP is to choose a subset of facilities and to fully assign one to each customer in such a way that the total cost is minimized and the capacities of the facilities are not exceeded. The total cost includes setting-up costs to open the facilities plus operating and transportation costs required to satisfy customers’s demand. As far as we know, this is the first study on the SS-MCLP although the Modular Capacitated Location Problem (MCLP), in which each customer may be partially served by more than one facility, has been studied by Correia and Captivo [1]. The MCLP is a generalization of the Staircase Cost Facility location problem presented in Holmberg [2], Holmberg and Ling [3] and Harkness and Revelle [4]. This study was prompted by a similar problem that arose during the design of an Internet access network by a telecommunications company in Portugal. This problem involves determining where to locate the equipment available in 19 different sizes corresponding to different costs and different levels of traffic. The SS-MCLP is a combinatorial optimization problem belonging to the NP-hard class. This means that, as the problem’s dimension becomes greater, heuristic methods are the only alternative to determine feasible solutions. To reinforce this idea, the computational results reported in this paper reveal problems involving 10 potential facility sites, three allowable capacities at each site and 100 customers, which cannot be solved optimally by the state-of-the-art solvers due to memory limitations. One kind of heuristics with many successful applications in the literature is of the lagrangean type. This is an iterative method that tries to build feasible solutions from the solution provided by a lagrangean relaxation of the problem. As for research in the discrete location area, Correia and Captivo [1] obtained good computational results for the MCLP using a lagrangean heuristic combined with penalty tests for problem reduction. Their heuristic is similar to the ones presented in Beasley [5,6] for the Capacitated Plant Location Problem (CPLP). Moreover, Cortinhal and Captivo [7], in a study of the Single Source Capacitated Plant Location Problem (SS-CPLP), reported good computational experiments with a lagrangean heuristic where the upper bounds were computed using Tabu search (TS). TS is a well-known meta-heuristic based on neighborhood search that tries to explore the solution space beyond local optimality (see Glover and Laguna [8] for details). Another applications of TS to location problems may be found in Klincewicz [9] and in Rolland et al.[10], for example. In this work, we adapt the formulations already developed for the MCLP to the SS-MCLP. The bounds on the problem’s optimal value are computed by a lagrangean heuristic in which the lower bounds are obtained by means of a lagrangean relaxation reinforced with an inequality corresponding to a lower bound on the total number of open facilities. The feasible solutions are calculated by a two-phase procedure. The objective of the first phase is to build an initial solution using a greedy algorithm. In the second phase, one applies a TS procedure or a local search (LS) scheme to try to improve the initial solution. Both methods are very flexible since the second phase may start from either a feasible or an infeasible solution. The lagrangean heuristic was tested on a set of randomly generated problems some of which may be considered to have large dimensions in the light of the published work in this field. This paper is structured in the following way. The mathematical formulation for the SS-MCLP is presented in Section 1. Section 2 presents the solution approach. Section 3 contains the computational experience and some conclusions are drawn in Section 4.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2993
1. Formulations for the SS-MCLP To formulate the SS-MCLP we consider the following data: a set I = {1, . . . , i, . . . , n} of potential facility sites; a set J = {1, . . . , j, . . . , m} of customers; a set Ki = {1, . . . , k, . . . , pi } of the pi possible sizes for the potential facility at location i (i ∈ I ); a demand aj of customer j (j ∈ J ); a capacity Bik for a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a lower bound bik on the shipments from a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a fixed-cost fik of establishing a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a unit cost qik of operating a size k facility at potential location i (i ∈ I ; k ∈ Ki )and a unit transportation cost dij from facility i (i ∈ I ) to customer j (j ∈ J ). The formulation involves the Wik variables representing the amount provided by a size k (k ∈ Ki ) facility located at i (i ∈ I ), the binary variables Yik indicating whether or not a facility of size k (k ∈ Ki ) is open at location i (i ∈ I ) and the binary variables Xij indicating whether customer j (j ∈ J ) is assigned to the facility open at location i (i ∈ I ). An integer formulation for the SS-MCLP is as follows: (Q1) Min
m n
dij aj Xij +
i=1 j =1
s.t.
n
pi n
qik Wik +
i=1 k=1
Xij = 1
(j ∈ J ),
pi n
fik Yik
(1)
i=1 k=1
(2)
i=1 pi
Yik 1
(i ∈ I ),
(3)
k=1
Wik Bik Yik
(i ∈ I ; k ∈ Ki ),
(4)
Wik bik Yik
(i ∈ I ; k ∈ Ki ),
(5)
m
aj Xij =
j =1
pi
Wik
(i ∈ I ),
(6)
k=1
Xij ∈ {0, 1}
(i ∈ I ; j ∈ J ),
(7)
Yik ∈ {0, 1}
(i ∈ I ; k ∈ Ki ).
(8)
The objective function (1) involves the minimization of total costs; constraints (2) specify that each customer is served by precisely one facility; constraints (3) ensure that no more than one facility can be open at each location i; constraints (4) imply both, that if no facility of size k is open at i the amount provided by this facility is zero and if a facility of that size is open then its maximum capacity cannot be exceeded; constraints (5) impose that if a facility of size k is open at i then it provides at least bik units, which means that these constraints correspond to the lower bounds on the variables Wik ; constraints (6) establish the connection between Xij and Wik variables and ensure that, for each location, the amount delivered to the customers is produced; (7) and (8) are the integrality requirements.In the remainder of this paper, and without loss of generality, we assume that all data are non-negative integers. We also consider that Bi1 < · · · < Bipi (i ∈ I ) and Bik + 1 bi(k+1) (i ∈ I ; k = 1, . . . , pi − 1).
2994
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
p i With the replacement of Xij = k=1 Zij k and Wik = m j =1 Zij k , in model Q1, that is, with the disaggregation of the Xij and Wik variables, the following model is obtained: (Q2) Min
s.t.
m n i=1 j =1 pi n
(dij + qik )aj Zij k +
pi n
(9)
fik Yik
i=1 k=1
Zij k = 1
i=1 k=1 pi
Yik 1
(j ∈ J ),
(10)
(i ∈ I ),
(11)
k=1
m
aj Zij k Bik Yik
(i ∈ I ; k ∈ Ki ),
(12)
aj Zij k bik Yik
(i ∈ I ; k ∈ Ki ),
(13)
Zij k ∈ {0, 1} (i ∈ I ; j ∈ J ; k ∈ Ki ), Yik ∈ {0, 1} (i ∈ I ; k ∈ Ki ).
(14) (15)
j =1 m j =1
Here each Zij k variable takes the value one if customer j is assigned to a size k facility in location i and takes the value zero otherwise. Q1 and Q2 are integer linear models. Despite the combinatorial nature of these problems, with a reasonable effort general integer software is able to solve some instances of the SS-MCLP as we shall see in Section 3. Unfortunately, as will be found later, some of our small testing instances, could not be solved in this way and require that we develop heuristic methods. 2. Solution approach 2.1. Lagrangean relaxation In order to develop a lagrangean heuristic one must choose some constraints of the problem to be relaxed in a lagrangean way (for details about lagrangean relaxation see Geoffrion [11]). The relaxation of constraints (6) in Q1 generates the lower bounds already obtained for the MCLP (for details about the resolution procedure see Correia and Captivo [1]). With the lagrangean relaxation of constraints (10) in Q2, the following model was obtained: pi pi pi m n n m n (Q2 ) Min (dij + qik )aj Zij k + fik Yik + j 1 − Zij k i=1 j =1 k=1
s.t.
i=1 k=1
j =1
i=1 k=1
(11), (12), (13), (14), (15).
Q2 can be decomposed into knapsack problems, one knapsack for each i and k value, with lower capacity constraints. In this situation the relaxed problem also belongs to the NP-hard class. To overcome
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2995
this difficulty we worked with the model RQ2 , which was obtained by replacing in Q2 the constraints Zij k ∈ {0, 1} by 0 Zij k 1. Since RQ2 has been studied for the MCLP, details about its solution procedure are available in Correia and Captivo [1]. As previously pointed out in Correia and Captivo [1], RQ2 may produce better lower bounds than Q1 because it does not satisfy the integrality property. That is why, in the rest of the paper, the model Q2 and the lagrangean relaxation RQ2 will be used. Maximization of RQ2 was performed using the subgradient optimization procedure proposed by Held et al. [12]. In an attempt to strengthen the lower bounds given by RQ2 , we added to it the valid inequality (16) that states a lower limit R (for details about the computation of R see Correia and Captivo [1]) on the total number of open facilities: pi n
Yik R.
(16)
i=1 k=1
2.2. Upper bounds The upper bounds were computed by means of a two-phase heuristic, namely a constructive phase and an improvement phase. In the constructive phase, a greedy-type heuristic tries to build an initial feasible solution from the set of open facilities provided by the lagrangean relaxation. The aim of the second phase is to improve the initial solution by using a tabu search meta-heuristic or by resorting to a local search procedure. The description of the two phases is given below. 2.2.1. Greedy heuristic (GH) Let us consider F = {(i, k) : a size k facility is open at i, i ∈ I, k ∈ Ki }. Hence F is a feasible set of open facilities if no more than one facility is open at each potential location and (i,k)∈F bik j ∈J aj (i,k)∈F Bik . The preceding inequalities are necessary conditions to enable a feasible assignment for the SS-MCLP. At the end of each subgradient iteration, a solution to the relaxed problem with at least R open facilities is available. If in this solution, all the customers are assigned to only one open facility, then a feasible solution to the SS-MCLP has been found. Unfortunately, this situation is unlikely to happen and more work is necessary to obtain such a feasible assignment. Let I t be the set of facilities opened during the iteration t of the subgradient. If I t does not satisfy the necessary conditions, the set I t is modified by the same procedure described in Correia and Captivo [1] until those inequalities are verified. If the procedure fails to find a feasible set of open facilities, no feasible solution is obtained. Otherwise, the next steps must be followed in order to get a feasible assignment of the customers to the open facilities. Step 1: For each open facility (i, k), (i, k) ∈ I t , define Ji ={j ∈ J : Xij k =1} and i ={Bik − j ∈Ji aj }, that is, Ji contains the customers served by the facility open in location i and i is the capacity still available in location i. Initially set Ji = ∅ and i = Bik . Step 2: For each customer j, j ∈ J , define Ij = {(i, k) ∈ I t : i − aj 0}, that is, Ij is the set of open facilities with enough capacity to serve customer j. Step 3: Sort the customers by descending order of their demand. Step 4: Following the order, for all customers j, j ∈ J , perform (a): If Ij = ∅, assign j to the facility (i∗, k∗) : di∗j + qi∗k∗ = min(i,k)∈Ij {dij + qik }. Update Ji∗ and i∗ accordingly.
2996
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
(b): If Ij = ∅, assign j to the facility (i , k ) : di j + qi k = min(i,k)∈I t {dij + qik }. Update Ji and i accordingly. The previous procedure is based on a similar approach proposed in Cortinhal and Captivo [7] for the SS-CPLP. When the procedure stops, all the customers are totally assigned to one open facility. However, if for the assignment of some customer step 4(b) is needed, this means that the Greedy heuristic will provide a solution that violates the capacity constraints, that is, an infeasible one. Even so, the solution moves on to the improvement phase because, as we shall see later, our TS and LS algorithms were designed to handle this kind of situation. 2.2.2. Tabu search (TS) One of the ways to try to improve a solution for the SS-MCLP is by making changes in the set of open facilities and/or in the assignment of the customers. Thus, our TS algorithm has two steps. In the assignment step, a TS algorithm is applied on the assignment variables in order to get a good assignment. During the location step, a TS procedure tries to improve the set of open facilities. The idea of using a two step TS approach, each related to a type of decision required by a given problem, is not new. In fact, this strategy has already been used by Tuzun and Burke [13] for the Location Routing Problem. In the remainder of this section, we explain the main features of the TS algorithm: the aspiration criteria, the move evaluation, the strategic oscillation and the intensification. Assignment step: The neighborhood is explored in this step through shift and swap moves. In a shift move the assignment of one customer is switched to another open facility, whereas a swap move interchanges the assignment between two customers originally assigned to different open facilities. A shift move has the advantage of generating solutions where the number of customers assigned to each open facility may vary. A swap move searches for a good solution involving a certain number of customers served by each open facility. In this way, these two kinds of moves give added flexibility to the search. The assignment step consists in alternately performing shift or swap moves until max1 consecutive moves are executed with no improvement in the best known upper bound. One of the purposes of TS is to encourage the exploration of parts of the solution space that have not previously been visited. This can be done by prohibiting the reversal of previous visited solutions. Hence, as shift moves are performed, the previous assignment of each customer is recorded in a tabu list. A shift move is tabu if it changes the assignment of a customer to a forbidden facility. A swap move is tabu if it assigns at least one customer to a forbidden facility. When a move is declared tabu, it stays in this state during u1 + v1 ∗ /max iterations where u1 and v1 are given parameters, is the number of times this assignment has already been performed and max is the maximum value of for all the assignments undertaken throughout the whole lagrangean heuristic (this kind of tabu tenure was also used in Cortinhal and Captivo [7]). More precisely, is a frequency matrix whose (i, j, k) element represents the number of times customer j was assigned to the size k facility at location i. With this kind of tabu tenure, the neighborhood of a given assignment configuration is not static, but rather changes with the history of the search. Location step: Once the improvement phase reaches this step, a swap move is executed and the method returns to the Assignment step unless max2 consecutive moves have already been made with no improvement in the best known upper bound. Given a certain configuration for the facilities, a swap move in the set of open facilities simultaneously closes one of the open facilities and opens a facility at another potential location. All the customers
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2997
assigned to the facility that has been closed are reassigned to the new facility, although an infeasible solution may be obtained. The main disadvantage of this kind of move is that it only explores solutions with a certain number of open facilities. However the introduction of ADD and/or DROP moves, that is, moves which can modify the number of open facilities in a solution, would be difficult to evaluate due to the major difficulties in the reassignment process. After a swap move is performed, both locations involved are declared tabu during u2 +v2 ∗ ik / max ik where u2 and v2 are given parameters, ik denotes the number of times that facility (i, k) was opened up to this point and max ik represent the maximum value for ik . We consider that a swap move is tabu if either its open facility or its close facility are tabu. Two-phase TS procedure: The previous description of this two-phase TS procedure has been condensed in the next pseudo-code, where the best known upper bound will be denoted by UB. Step 1: Initialize max1 , max2 , t1 = 0, t2 = 0 and all the other parameters. Step 2: Assignment phase—Perform a shift or a swap move. If UB does not improve then t1 = t1 + 1 else t1 = 0. Step 3: If t1 = max1 then goto step 4 else goto step 2. Step 4: Location phase - Perform a swap move between two facilities. If UB does not improve then t2 = t2 + 1 else t2 = 0. If t2 = max2 then stop else goto step 2. Aspiration criteria: The aspiration criteria is a rule used to override the tabu status of moves once better solutions are encountered. In this context, we consider that a better solution is the one that improves the best known upper bound. Move evaluation: The three types of moves described above can lead to feasible or infeasible solutions. A feasible solution is achieved if the capacity requirements of each facility involved in the move are still satisfied after the move has been performed, otherwise an infeasible solution is obtained. Nevertheless, when a move is evaluated, its cost is computed without considering the capacity constraints, which are only relevant to classify a move as feasible or infeasible. Since the moves are not costly to examine, the entire neighborhoods are explored. Strategic oscillation: With the strategic oscillation approach, the search is driven in the direction of a boundary represented, for example, by feasibility, where usually the method would stop. Instead of stopping, the search crosses that boundary. After some time spent in the infeasible region, the search turns around and the feasible region is reached again. Our TS method attains the infeasible region in the following situations: (i) if the best upper bound does not improve in max3 consecutive moves; (ii) if all the possible feasible moves are tabu and the aspiration criteria is not achieved or (iii) if during the Greedy heuristic any capacity constraints were violated. In the previous situations the best move is the one that minimizes infeasibility (ties are broken by the objective function value). Once outside the feasible region and with the purpose of returning to feasibility, the moves are performed in order to reduce the infeasibility of the current solution. One of the problems of the TS approach is that sometimes it tends to spend a significant amount of time in a restricted region of the search space. The strategic oscillation was also used with the purpose of allowing diversification. Intensification: The use of intensification strategies within a TS procedure drives the search to attractive regions in order to explore them more thoroughly. In our approach the intensification phase chooses the best feasible solution found so far. Next, % of the open facilities with the highest ik values are definitely
2998
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
fixed and the two-phase TS procedure, outlined above, is restarted. On completion of the intensification scheme all the variables become free. 2.2.3. Local search (LS) The LS procedure consists of two steps that are both descent local search methods. The aim of the first step is to try to get a better set of assignment variables while the second is devoted to improving of the open facilities set. In both steps the entire neighborhood is explored. Hence, in the first step, the best move is the one over all possible shift and swap moves that improves the current solution most. This step is repeated until no move exists that improves the current solution. In the second step, the best move is the one, over all the possible swap moves that most improves the current best solution. When such a move exists it is performed and the search returns to the first step, otherwise the LS procedure stops. As in the TS, this procedure can also start from an infeasible solution. In this situation, the best move is the one that decreases the degree of infeasibility most. 2.3. Lagrangean heuristic Let be the step length parameter of the subgradient method; maxiter the maximum number of iterations allowed; LB the maximum lower bound found; t the current subgradient iteration; LB t the lower bound obtained in iteration t. An outline of the lagrangean heuristic can be stated as follows: Step 1: Initialize t = 1, tj = 0 (j ∈ J ), U B = +∞, LB = 0, = 1.5, N = 0, maxiter, N1 and N2 . Step 2: Solve RQ2t . Step 3: If (LB t > LB) then LB = LB t , N = 0 else N = N + 1. Step 4: If (N = N1 ) then N = 0, = /2. Step 5: Apply the Greedy heuristic to try to find a feasible solution. Step 6: If the Greedy heuristic fails then goto step 9 else goto step 7. Step 7: Improve the solution obtained in the previous step by the TS (with intensification) or by the LS procedure. Step 8: Update the best upper bound UB. Step 9: If (U B − LB) < 1 then stop. Step 10: Update the lagrangean multipliers and t. Step 11: If t is a multiple of N2 then = 1.5 Step 12: If (t maxiter) and the subgradient vector is not null then return to step 2. 3. Computational experience All the computational tests were run on a PC with a Pentium IV processor, a 2.6 GHz internal clock and 512 Mb of RAM memory. The lagrangean heuristic was written in Fortran 77 language. We tried to obtain the optimal value of the test problems in order to evaluate the quality of the upper bounds in terms of their deviation from the optimum. In this context, the MIP optimizer of the Cplex 8.0 [14]was used. The input files for Cplex were generated by means of formulation Q1 together with constraint (16) because, due to the large number of binary variables in Q2, Cplex performs poorly with this later formulation. In order to evaluate the performance of the heuristic and because no standard set of test problems is available for the SS-MCLP, we worked with a set of randomly generated test problems (which is available
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2999
from the authors on request) already used for the MCLP, where the additional single source constraint was imposed. In the following, the main characteristics of the problem instances are briefly described (for more explanations about the generation of the test problems, the reader is referred to Correia and Captivo [1]) and the computational results are reported. 3.1. Test problems Let p be the number of possible sizes for each potential location. Five problems were generated for each of the following dimensions: n ∈ {10, 20, 50, 200, 500}, m ∈ {100, 200, 1000} and p ∈ {3, 5}. Thus, 110 test problems were obtained since only problems where n m were generated. The test problems with n = 500, m = 1000 and p = 5 have 502,500 binary variables for formulation Q1 but this number increases to 2,502,500 for formulation Q2. Therefore, large problems were used, regarding the literature in the discrete location area. In studying this kind of problem it is usual to analyze how the performance of the methods is influenced n m by the ratio between total supply and total demand. Let us consider q = . For / B a ip j i=1 j =1 i the problems described above, on average, q = 2.5. In order to evaluate the q influence, all the previous test problems were duplicated and in each data set the demand of each customer was increased by 20 units, giving rise to a mean value for q of 1.9. 3.2. Computational results After some tuning (by means of several preliminary experiments) with the objective of getting a good trade-off between solution quality and computation time, the following parameter values for the TS were used: max1 = 20, max2 = 7, max3 = 25, u1 = 7, v1 = m, u2 = 2, v2 = n and = 50%. For problems up to fifty potential locations the maxiter, N1 and N2 parameters were set to 500, 15 and 50, respectively. For problems where n=200 or n=500 the maxiter parameter was set to 2000, while the N1 and N2 parameters were set, respectively, to 50 and 150. A feasible solution (step 5 of the lagrangean heuristic) was computed whenever a different set of open facilities was obtained while resolving the relaxed problem. In all the experiments using the Cplex solver, the time limit was set to 10 h. Initially the Cplex ran with all the other parameters as default, but by making this choice and, in the case of many test problems, Cplex stopped due to memory problems before an optimal solution had been found. To try to overcome this situation, the strong branching rule to determine the best branching variable in each tree node was used (for details, see the Cplex user’s manual) and all the Cplex experiments were repeated. With this change the Cplex was able to find the optimal value of more test problems. The computational results obtained are summarized in Tables 1–4 wherein are recorded the mean gaps (in percentage) and mean execution times (in seconds) taken over 10 test problems (5 test problems with p = 3 and 5 with p = 5). For some algorithms an additional column has been included. The column “n.fails” shows the number of test problems in which the algorithm under consideration has failed the computation of a feasible solution for each set of 10 test problems. When a line of a table has “n.fails”=nf , this means, the mean gap and time in that line were computed for the remaining (10 − nf ) test problems. The gaps were computed by the expression [(U B −V )/V ∗100] where UB is the upper bound provided by each method and V denotes the optimal value, if known, otherwise it represents the maximum between the lower bound provided by RQ2 and the best lower bound at Cplex termination.
3000
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
Table 1 Results for test problems with q = 2.5 n
10 10 10 20 20 20 50 50 50 200 500
m
Lagrangean heuristic
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
Tabu search(TS)
Local search(LS)
Greedy heuristic(GH)
gap(%)
time(s)
gap(%)
time(s)
gap(%)
time(s)
n.fails
0.15 0.13 0.37 0.82 1.14 1.19 2.40 1.47 1.43 1.13 2.21 1.13
1.64 4.64 115.39 3.87 14.13 281.44 12.13 31.40 766.26 2957.93 6560.58 977.22
0.58 0.25 0.95 1.24 1.56 1.28 3.10 2.01 1.60 1.41 2.50 1.50
0.23 0.52 20.76 0.36 1.03 33.03 0.96 2.17 58.41 227.80 1666.11 182.85
1.94 1.85 2.92 2.51 3.21 1.96 4.21 2.90 2.58 1.92 3.23 2.66
0.17 0.29 1.38 0.30 0.53 2.70 0.70 1.29 6.38 105.60 274.14 35.77
0 0 0 0 0 0 0 0 0 0 6 0.55
Table 2 Results for test problems with q = 2.5 n
m
Cplex Best at Cplex termination
10 10 10 20 20 20 50 50 50 200 500
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
First integer(FI)
gap(%)
time(s)
n.fails
gap(%)
time(s)
n.fails
0.00 0.00 0.00 0.33 0.53 0.06 2.31 1.00 0.96 1.70 3.35 0.93
4.32 500.32 67.30 5686.29 14319.03 3801.72 30732.89 21431.60 35911.25 36006.36 36018.51 16770.87
0 0 0 0 0 0 0 0 0 0 5 0.45
11.15 7.65 1.98 5.36 4.67 6.55 5.33 4.08 3.81 2.59 3.54 5.16
0.11 0.58 11.70 0.32 0.69 10.47 5.73 3.11 49.14 516.96 2942.02 321.89
0 0 0 0 0 0 0 0 0 0 5 0.45
Tables 1 and 3 show the mean gaps and execution times obtained by the lagrangean heuristic when feasible solutions were computed by using the following algorithms: TS, LS and GH (in this situation, step 7 of the lagrangean heuristic was not executed). In both tables wee observe that the GH always produced the worst mean gaps while TS produced the best ones. GH was the fastest, followed by LS and lastly by TS.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
3001
Table 3 Results for test problems with q = 1.9 n
10 10 10 20 20 20 50 50 50 200 500
m
Lagrangean heuristic
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
Tabu search(TS)
Local search(LS)
Greedy heuristic(GH)
gap(%)
time(s)
gap(%)
time(s)
gap(%)
time(s)
n.fails
0.77 1.24 1.28 1.58 1.38 1.12 3.06 2.01 1.36 1.44 3.68 1.72
1.89 5.19 146.21 3.03 11.08 252.77 12.29 27.33 626.84 3157.01 6396.40 967.28
1.41 1.52 2.00 1.92 2.04 1.70 3.91 2.60 1.71 1.60 4.24 2.24
0.23 0.58 25.18 0.37 0.96 30.26 1.17 2.48 54.81 318.22 2038.87 224.83
4.06 5.48 3.60 5.08 4.37 3.20 6.33 4.30 3.92 2.30 — 4.26
0.18 0.31 1.50 0.31 0.69 3.05 0.76 1.35 6.90 115.45 — 13.05
0 0 0 0 0 0 1 0 0 6 10 1.55
Table 4 Results for test problems with q = 1.9 n
m
Cplex Best at Cplex termination
10 10 10 20 20 20 50 50 50 200 500
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
First integer(FI)
gap(%)
time(s)
n.fails
gap(%)
time(s)
n.fails
0.17 0.19 0.33 1.04 1.21 0.43 2.51 2.03 1.11 2.41 4.80 1.48
7267.54 7908.41 13842.93 29394.45 23131.32 19395.20 34043.98 17915.09 36002.52 36005.87 36018.14 23720.50
0 0 0 0 0 0 0 0 0 0 5 0.45
5.02 4.37 4.01 3.60 3.67 3.70 5.50 3.87 3.02 2.96 4.90 4.06
0.17 0.50 13.74 0.60 1.25 3.46 39.36 4.08 33.73 1024.56 5410.44 593.81
0 0 0 0 0 0 0 0 0 0 5 0.45
All the tables show that the instances with q = 1.9 were harder to solve than the ones with q = 2.5. This conclusion is not surprising since this kind of result was also reported in Cortinhal and Captivo [7] for the SS-CPLP. TS and LS found feasible solutions for all test problems but the GH failed in 6 test problems with q = 2.5 and in 17 of the instances with q = 1.9. As for the sizes of the instances, the mean execution times increase as the dimension of the problems increases, but this conclusion is not valid for the mean gaps.
3002
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
Tables 2 and 4 report the experiments with Cplex. As we pointed out above, for each instance, Cplex stopped operating when one of the following conditions hold: an optimal solution has been found within 10 h, the time limit of 10 h was attained or an out of memory problem occurs. Because Cplex sometimes did not find an optimal solution, we present the columns with the label “Best at Cplex termination”. The other columns show the results computed with the First Integer (FI) solution provided by Cplex. The computation of the FI by Cplex did not seem a good alternative to the lagrangean heuristic. In fact, FI is not competitive since with LS, all instances were successfully solved and it gave a better mean gap and time. The average gaps and times provided by FI were worst than the ones provided by GH although FI was able to solve one more test problem. In relation to the TS, it was found that TS clearly dominates FI in terms of the average gaps and number of instances unsolved but TS proved to be more time consuming. Regarding the attempts to optimally solve the test problems, Table 4 shows that even in the case of the problems with n = 10 the gap’s column has nonzero entries because Cplex was unable to optimally solve those instances. Despite the fact that average gaps presented in the column “Best at Cplex termination” are the best, they are not much lower than the ones provided by TS and TS was able to solve all the problems with significantly less computational effort. 4. Concluding remarks This paper addresses a Single Source Capacitated Plant Location Problem with several possible capacity levels for the facility that can be opened at each potential location. Because this is a hard problem to solve a lagrangean heuristic was developed. We propose a Tabu Search and a Local Search algorithm to improve the feasible solutions during the lagrangean heuristic. In order to check the efficiency of the methods proposed in this paper, some test problems were randomly generated and some experiments were performed using the Cplex solver. The heuristic methods provide very satisfactory results regarding the mean gaps and execution times even for large scale problems. As is usual with this kind of problem, instances with tighter capacity constraints were harder to solve. Future research may focus on the development of different meta-heuristic procedures to try to improve the quality of the feasible solutions. Acknowledgements The authors thank the referees for their helpful comments. References [1] Correia I, Captivo ME. A Lagrangean heuristics for a modular capacitated location problem. Annals of Operations Research 2003;122:141–61. [2] Holmberg K. Solving the staircase cost facility location problem with decomposition and piecewise linearization. European Journal of Operational Research 1994;75:41–61. [3] Holmberg K, Ling J. A lagrangean heuristic for the facility location problem with staircase costs. European Journal of Operational Research 1997;97:63–74.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
3003
[4] Harkness J, Revelle C. A Lagrangean heuristic for the staircase cost facility location problem. International Journal of Industrial Engineering-Theory, Applications and Practice 2002;9(1):36–44. [5] Beasley J. An algorithm for solving large capacitated warehouse location problems. European Journal of Operational Research 1988;33:314–25. [6] Beasley J. Lagrangean heuristics for location problems. European Journal of Operational Research 1993;65:383–99. [7] Cortinhal MJ, Captivo ME. Upper and lower bounds for the single source capacitated location problem. European Journal of Operational Research 2003;151(2):333–51. [8] Glover F, Laguna M. Tabu search. Dordrecht: Kluwer Academic Publishers; 1997. [9] Klincewicz J. Avoiding local optima in the p-hub location problem using Tabu Search and GRASP. Annals of Operations Research 1992;40:283–302. [10] Rolland E, Schilling D, Current J. An efficient Tabu Search procedure for the p-median problem. European Journal of Operational Research 1997;96:329–42. [11] Geoffrion A. Lagrangean relaxation for integer programming. Mathematical Programming Study 1974;2:82–114. [12] Held M, Wolfe P, Crowder H. Validation of subgradient optimization. Mathematical Programming 1974;6:62–6. [13] Tuzun D, Burke L. A two-phase tabu search approach to the location routing problem. European Journal of Operational Research 1999;116:87–99. [14] Cplex. Using the Cplex callable library, Version 8.0, ILOG, Inc., 2002.
Bounds for the single source modular capacitated plant location problem Isabel Correiaa,∗ , Maria Eugénia Captivob a Departamento de Matemática - C.M.A., Faculdade de Ciências e Tecnologia, Universidade Nova de Lisboa,
Quinta da Torre, 2829-516 Monte da Caparica, Portugal b Universidade de Lisboa, Faculdade de Ciências, Centro de Investigação Operacional Bloco C6, Piso 4,
Campo Grande, 1749-016 Lisboa, Portugal Available online 8 April 2005
Abstract In this paper, we propose a discrete location problem, which we call the Single Source Modular Capacitated Location Problem (SS-MCLP). The problem consists of finding the location and capacity of the facilities, to serve a set of customers at a minimum total cost. The demand of each customer must be satisfied by one facility only and the capacities of the open facilities must be chosen from a finite and discrete set of allowable capacities. Because the SS-MCLP is a difficult problem, a lagrangean heuristic, enhanced by tabu search or local search was developed in order to obtain good feasible solutions. When needed, the lower bounds are used in order to evaluate the quality of the feasible solutions. Our method was tested computationally on randomly generated test problems some of which are with large dimensions considering the literature related to this type of problem. The computational results obtained were compared with those provided by the commercial software Cplex. 䉷 2005 Elsevier Ltd. All rights reserved. Keywords: Capacitated location; Lagrangean heuristic; Tabu search
0. Introduction Location problems have been widely studied in the last years due to their application in many real situations. In fact, for many companies one of the key strategic decisions is to decide where to build ∗ Corresponding author. Tel.: +351 212948388; fax: +351 212948391.
E-mail address: [email protected] (I. Correia). 0305-0548/$ - see front matter 䉷 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.cor.2005.02.030
2992
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
facilities (e.g., warehouses, concentrators) and how to serve customers in order to optimize a given objective. These kinds of problems are not restricted to plant location alone, they also crop up in cases, such as design of telecommunication networks and distribution systems planning. In the SS-MCLP, a set of potential facility locations and a set of customers, each with a known demand, are given. In each potential location, there is a finite and discrete set of allowable capacities for the facilities and the demand of each customer must be satisfied by a single facility. The purpose of the SS-MCLP is to choose a subset of facilities and to fully assign one to each customer in such a way that the total cost is minimized and the capacities of the facilities are not exceeded. The total cost includes setting-up costs to open the facilities plus operating and transportation costs required to satisfy customers’s demand. As far as we know, this is the first study on the SS-MCLP although the Modular Capacitated Location Problem (MCLP), in which each customer may be partially served by more than one facility, has been studied by Correia and Captivo [1]. The MCLP is a generalization of the Staircase Cost Facility location problem presented in Holmberg [2], Holmberg and Ling [3] and Harkness and Revelle [4]. This study was prompted by a similar problem that arose during the design of an Internet access network by a telecommunications company in Portugal. This problem involves determining where to locate the equipment available in 19 different sizes corresponding to different costs and different levels of traffic. The SS-MCLP is a combinatorial optimization problem belonging to the NP-hard class. This means that, as the problem’s dimension becomes greater, heuristic methods are the only alternative to determine feasible solutions. To reinforce this idea, the computational results reported in this paper reveal problems involving 10 potential facility sites, three allowable capacities at each site and 100 customers, which cannot be solved optimally by the state-of-the-art solvers due to memory limitations. One kind of heuristics with many successful applications in the literature is of the lagrangean type. This is an iterative method that tries to build feasible solutions from the solution provided by a lagrangean relaxation of the problem. As for research in the discrete location area, Correia and Captivo [1] obtained good computational results for the MCLP using a lagrangean heuristic combined with penalty tests for problem reduction. Their heuristic is similar to the ones presented in Beasley [5,6] for the Capacitated Plant Location Problem (CPLP). Moreover, Cortinhal and Captivo [7], in a study of the Single Source Capacitated Plant Location Problem (SS-CPLP), reported good computational experiments with a lagrangean heuristic where the upper bounds were computed using Tabu search (TS). TS is a well-known meta-heuristic based on neighborhood search that tries to explore the solution space beyond local optimality (see Glover and Laguna [8] for details). Another applications of TS to location problems may be found in Klincewicz [9] and in Rolland et al.[10], for example. In this work, we adapt the formulations already developed for the MCLP to the SS-MCLP. The bounds on the problem’s optimal value are computed by a lagrangean heuristic in which the lower bounds are obtained by means of a lagrangean relaxation reinforced with an inequality corresponding to a lower bound on the total number of open facilities. The feasible solutions are calculated by a two-phase procedure. The objective of the first phase is to build an initial solution using a greedy algorithm. In the second phase, one applies a TS procedure or a local search (LS) scheme to try to improve the initial solution. Both methods are very flexible since the second phase may start from either a feasible or an infeasible solution. The lagrangean heuristic was tested on a set of randomly generated problems some of which may be considered to have large dimensions in the light of the published work in this field. This paper is structured in the following way. The mathematical formulation for the SS-MCLP is presented in Section 1. Section 2 presents the solution approach. Section 3 contains the computational experience and some conclusions are drawn in Section 4.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2993
1. Formulations for the SS-MCLP To formulate the SS-MCLP we consider the following data: a set I = {1, . . . , i, . . . , n} of potential facility sites; a set J = {1, . . . , j, . . . , m} of customers; a set Ki = {1, . . . , k, . . . , pi } of the pi possible sizes for the potential facility at location i (i ∈ I ); a demand aj of customer j (j ∈ J ); a capacity Bik for a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a lower bound bik on the shipments from a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a fixed-cost fik of establishing a size k facility at potential location i (i ∈ I ; k ∈ Ki ); a unit cost qik of operating a size k facility at potential location i (i ∈ I ; k ∈ Ki )and a unit transportation cost dij from facility i (i ∈ I ) to customer j (j ∈ J ). The formulation involves the Wik variables representing the amount provided by a size k (k ∈ Ki ) facility located at i (i ∈ I ), the binary variables Yik indicating whether or not a facility of size k (k ∈ Ki ) is open at location i (i ∈ I ) and the binary variables Xij indicating whether customer j (j ∈ J ) is assigned to the facility open at location i (i ∈ I ). An integer formulation for the SS-MCLP is as follows: (Q1) Min
m n
dij aj Xij +
i=1 j =1
s.t.
n
pi n
qik Wik +
i=1 k=1
Xij = 1
(j ∈ J ),
pi n
fik Yik
(1)
i=1 k=1
(2)
i=1 pi
Yik 1
(i ∈ I ),
(3)
k=1
Wik Bik Yik
(i ∈ I ; k ∈ Ki ),
(4)
Wik bik Yik
(i ∈ I ; k ∈ Ki ),
(5)
m
aj Xij =
j =1
pi
Wik
(i ∈ I ),
(6)
k=1
Xij ∈ {0, 1}
(i ∈ I ; j ∈ J ),
(7)
Yik ∈ {0, 1}
(i ∈ I ; k ∈ Ki ).
(8)
The objective function (1) involves the minimization of total costs; constraints (2) specify that each customer is served by precisely one facility; constraints (3) ensure that no more than one facility can be open at each location i; constraints (4) imply both, that if no facility of size k is open at i the amount provided by this facility is zero and if a facility of that size is open then its maximum capacity cannot be exceeded; constraints (5) impose that if a facility of size k is open at i then it provides at least bik units, which means that these constraints correspond to the lower bounds on the variables Wik ; constraints (6) establish the connection between Xij and Wik variables and ensure that, for each location, the amount delivered to the customers is produced; (7) and (8) are the integrality requirements.In the remainder of this paper, and without loss of generality, we assume that all data are non-negative integers. We also consider that Bi1 < · · · < Bipi (i ∈ I ) and Bik + 1 bi(k+1) (i ∈ I ; k = 1, . . . , pi − 1).
2994
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
p i With the replacement of Xij = k=1 Zij k and Wik = m j =1 Zij k , in model Q1, that is, with the disaggregation of the Xij and Wik variables, the following model is obtained: (Q2) Min
s.t.
m n i=1 j =1 pi n
(dij + qik )aj Zij k +
pi n
(9)
fik Yik
i=1 k=1
Zij k = 1
i=1 k=1 pi
Yik 1
(j ∈ J ),
(10)
(i ∈ I ),
(11)
k=1
m
aj Zij k Bik Yik
(i ∈ I ; k ∈ Ki ),
(12)
aj Zij k bik Yik
(i ∈ I ; k ∈ Ki ),
(13)
Zij k ∈ {0, 1} (i ∈ I ; j ∈ J ; k ∈ Ki ), Yik ∈ {0, 1} (i ∈ I ; k ∈ Ki ).
(14) (15)
j =1 m j =1
Here each Zij k variable takes the value one if customer j is assigned to a size k facility in location i and takes the value zero otherwise. Q1 and Q2 are integer linear models. Despite the combinatorial nature of these problems, with a reasonable effort general integer software is able to solve some instances of the SS-MCLP as we shall see in Section 3. Unfortunately, as will be found later, some of our small testing instances, could not be solved in this way and require that we develop heuristic methods. 2. Solution approach 2.1. Lagrangean relaxation In order to develop a lagrangean heuristic one must choose some constraints of the problem to be relaxed in a lagrangean way (for details about lagrangean relaxation see Geoffrion [11]). The relaxation of constraints (6) in Q1 generates the lower bounds already obtained for the MCLP (for details about the resolution procedure see Correia and Captivo [1]). With the lagrangean relaxation of constraints (10) in Q2, the following model was obtained: pi pi pi m n n m n (Q2 ) Min (dij + qik )aj Zij k + fik Yik + j 1 − Zij k i=1 j =1 k=1
s.t.
i=1 k=1
j =1
i=1 k=1
(11), (12), (13), (14), (15).
Q2 can be decomposed into knapsack problems, one knapsack for each i and k value, with lower capacity constraints. In this situation the relaxed problem also belongs to the NP-hard class. To overcome
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2995
this difficulty we worked with the model RQ2 , which was obtained by replacing in Q2 the constraints Zij k ∈ {0, 1} by 0 Zij k 1. Since RQ2 has been studied for the MCLP, details about its solution procedure are available in Correia and Captivo [1]. As previously pointed out in Correia and Captivo [1], RQ2 may produce better lower bounds than Q1 because it does not satisfy the integrality property. That is why, in the rest of the paper, the model Q2 and the lagrangean relaxation RQ2 will be used. Maximization of RQ2 was performed using the subgradient optimization procedure proposed by Held et al. [12]. In an attempt to strengthen the lower bounds given by RQ2 , we added to it the valid inequality (16) that states a lower limit R (for details about the computation of R see Correia and Captivo [1]) on the total number of open facilities: pi n
Yik R.
(16)
i=1 k=1
2.2. Upper bounds The upper bounds were computed by means of a two-phase heuristic, namely a constructive phase and an improvement phase. In the constructive phase, a greedy-type heuristic tries to build an initial feasible solution from the set of open facilities provided by the lagrangean relaxation. The aim of the second phase is to improve the initial solution by using a tabu search meta-heuristic or by resorting to a local search procedure. The description of the two phases is given below. 2.2.1. Greedy heuristic (GH) Let us consider F = {(i, k) : a size k facility is open at i, i ∈ I, k ∈ Ki }. Hence F is a feasible set of open facilities if no more than one facility is open at each potential location and (i,k)∈F bik j ∈J aj (i,k)∈F Bik . The preceding inequalities are necessary conditions to enable a feasible assignment for the SS-MCLP. At the end of each subgradient iteration, a solution to the relaxed problem with at least R open facilities is available. If in this solution, all the customers are assigned to only one open facility, then a feasible solution to the SS-MCLP has been found. Unfortunately, this situation is unlikely to happen and more work is necessary to obtain such a feasible assignment. Let I t be the set of facilities opened during the iteration t of the subgradient. If I t does not satisfy the necessary conditions, the set I t is modified by the same procedure described in Correia and Captivo [1] until those inequalities are verified. If the procedure fails to find a feasible set of open facilities, no feasible solution is obtained. Otherwise, the next steps must be followed in order to get a feasible assignment of the customers to the open facilities. Step 1: For each open facility (i, k), (i, k) ∈ I t , define Ji ={j ∈ J : Xij k =1} and i ={Bik − j ∈Ji aj }, that is, Ji contains the customers served by the facility open in location i and i is the capacity still available in location i. Initially set Ji = ∅ and i = Bik . Step 2: For each customer j, j ∈ J , define Ij = {(i, k) ∈ I t : i − aj 0}, that is, Ij is the set of open facilities with enough capacity to serve customer j. Step 3: Sort the customers by descending order of their demand. Step 4: Following the order, for all customers j, j ∈ J , perform (a): If Ij = ∅, assign j to the facility (i∗, k∗) : di∗j + qi∗k∗ = min(i,k)∈Ij {dij + qik }. Update Ji∗ and i∗ accordingly.
2996
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
(b): If Ij = ∅, assign j to the facility (i , k ) : di j + qi k = min(i,k)∈I t {dij + qik }. Update Ji and i accordingly. The previous procedure is based on a similar approach proposed in Cortinhal and Captivo [7] for the SS-CPLP. When the procedure stops, all the customers are totally assigned to one open facility. However, if for the assignment of some customer step 4(b) is needed, this means that the Greedy heuristic will provide a solution that violates the capacity constraints, that is, an infeasible one. Even so, the solution moves on to the improvement phase because, as we shall see later, our TS and LS algorithms were designed to handle this kind of situation. 2.2.2. Tabu search (TS) One of the ways to try to improve a solution for the SS-MCLP is by making changes in the set of open facilities and/or in the assignment of the customers. Thus, our TS algorithm has two steps. In the assignment step, a TS algorithm is applied on the assignment variables in order to get a good assignment. During the location step, a TS procedure tries to improve the set of open facilities. The idea of using a two step TS approach, each related to a type of decision required by a given problem, is not new. In fact, this strategy has already been used by Tuzun and Burke [13] for the Location Routing Problem. In the remainder of this section, we explain the main features of the TS algorithm: the aspiration criteria, the move evaluation, the strategic oscillation and the intensification. Assignment step: The neighborhood is explored in this step through shift and swap moves. In a shift move the assignment of one customer is switched to another open facility, whereas a swap move interchanges the assignment between two customers originally assigned to different open facilities. A shift move has the advantage of generating solutions where the number of customers assigned to each open facility may vary. A swap move searches for a good solution involving a certain number of customers served by each open facility. In this way, these two kinds of moves give added flexibility to the search. The assignment step consists in alternately performing shift or swap moves until max1 consecutive moves are executed with no improvement in the best known upper bound. One of the purposes of TS is to encourage the exploration of parts of the solution space that have not previously been visited. This can be done by prohibiting the reversal of previous visited solutions. Hence, as shift moves are performed, the previous assignment of each customer is recorded in a tabu list. A shift move is tabu if it changes the assignment of a customer to a forbidden facility. A swap move is tabu if it assigns at least one customer to a forbidden facility. When a move is declared tabu, it stays in this state during u1 + v1 ∗ /max iterations where u1 and v1 are given parameters, is the number of times this assignment has already been performed and max is the maximum value of for all the assignments undertaken throughout the whole lagrangean heuristic (this kind of tabu tenure was also used in Cortinhal and Captivo [7]). More precisely, is a frequency matrix whose (i, j, k) element represents the number of times customer j was assigned to the size k facility at location i. With this kind of tabu tenure, the neighborhood of a given assignment configuration is not static, but rather changes with the history of the search. Location step: Once the improvement phase reaches this step, a swap move is executed and the method returns to the Assignment step unless max2 consecutive moves have already been made with no improvement in the best known upper bound. Given a certain configuration for the facilities, a swap move in the set of open facilities simultaneously closes one of the open facilities and opens a facility at another potential location. All the customers
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2997
assigned to the facility that has been closed are reassigned to the new facility, although an infeasible solution may be obtained. The main disadvantage of this kind of move is that it only explores solutions with a certain number of open facilities. However the introduction of ADD and/or DROP moves, that is, moves which can modify the number of open facilities in a solution, would be difficult to evaluate due to the major difficulties in the reassignment process. After a swap move is performed, both locations involved are declared tabu during u2 +v2 ∗ ik / max ik where u2 and v2 are given parameters, ik denotes the number of times that facility (i, k) was opened up to this point and max ik represent the maximum value for ik . We consider that a swap move is tabu if either its open facility or its close facility are tabu. Two-phase TS procedure: The previous description of this two-phase TS procedure has been condensed in the next pseudo-code, where the best known upper bound will be denoted by UB. Step 1: Initialize max1 , max2 , t1 = 0, t2 = 0 and all the other parameters. Step 2: Assignment phase—Perform a shift or a swap move. If UB does not improve then t1 = t1 + 1 else t1 = 0. Step 3: If t1 = max1 then goto step 4 else goto step 2. Step 4: Location phase - Perform a swap move between two facilities. If UB does not improve then t2 = t2 + 1 else t2 = 0. If t2 = max2 then stop else goto step 2. Aspiration criteria: The aspiration criteria is a rule used to override the tabu status of moves once better solutions are encountered. In this context, we consider that a better solution is the one that improves the best known upper bound. Move evaluation: The three types of moves described above can lead to feasible or infeasible solutions. A feasible solution is achieved if the capacity requirements of each facility involved in the move are still satisfied after the move has been performed, otherwise an infeasible solution is obtained. Nevertheless, when a move is evaluated, its cost is computed without considering the capacity constraints, which are only relevant to classify a move as feasible or infeasible. Since the moves are not costly to examine, the entire neighborhoods are explored. Strategic oscillation: With the strategic oscillation approach, the search is driven in the direction of a boundary represented, for example, by feasibility, where usually the method would stop. Instead of stopping, the search crosses that boundary. After some time spent in the infeasible region, the search turns around and the feasible region is reached again. Our TS method attains the infeasible region in the following situations: (i) if the best upper bound does not improve in max3 consecutive moves; (ii) if all the possible feasible moves are tabu and the aspiration criteria is not achieved or (iii) if during the Greedy heuristic any capacity constraints were violated. In the previous situations the best move is the one that minimizes infeasibility (ties are broken by the objective function value). Once outside the feasible region and with the purpose of returning to feasibility, the moves are performed in order to reduce the infeasibility of the current solution. One of the problems of the TS approach is that sometimes it tends to spend a significant amount of time in a restricted region of the search space. The strategic oscillation was also used with the purpose of allowing diversification. Intensification: The use of intensification strategies within a TS procedure drives the search to attractive regions in order to explore them more thoroughly. In our approach the intensification phase chooses the best feasible solution found so far. Next, % of the open facilities with the highest ik values are definitely
2998
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
fixed and the two-phase TS procedure, outlined above, is restarted. On completion of the intensification scheme all the variables become free. 2.2.3. Local search (LS) The LS procedure consists of two steps that are both descent local search methods. The aim of the first step is to try to get a better set of assignment variables while the second is devoted to improving of the open facilities set. In both steps the entire neighborhood is explored. Hence, in the first step, the best move is the one over all possible shift and swap moves that improves the current solution most. This step is repeated until no move exists that improves the current solution. In the second step, the best move is the one, over all the possible swap moves that most improves the current best solution. When such a move exists it is performed and the search returns to the first step, otherwise the LS procedure stops. As in the TS, this procedure can also start from an infeasible solution. In this situation, the best move is the one that decreases the degree of infeasibility most. 2.3. Lagrangean heuristic Let be the step length parameter of the subgradient method; maxiter the maximum number of iterations allowed; LB the maximum lower bound found; t the current subgradient iteration; LB t the lower bound obtained in iteration t. An outline of the lagrangean heuristic can be stated as follows: Step 1: Initialize t = 1, tj = 0 (j ∈ J ), U B = +∞, LB = 0, = 1.5, N = 0, maxiter, N1 and N2 . Step 2: Solve RQ2t . Step 3: If (LB t > LB) then LB = LB t , N = 0 else N = N + 1. Step 4: If (N = N1 ) then N = 0, = /2. Step 5: Apply the Greedy heuristic to try to find a feasible solution. Step 6: If the Greedy heuristic fails then goto step 9 else goto step 7. Step 7: Improve the solution obtained in the previous step by the TS (with intensification) or by the LS procedure. Step 8: Update the best upper bound UB. Step 9: If (U B − LB) < 1 then stop. Step 10: Update the lagrangean multipliers and t. Step 11: If t is a multiple of N2 then = 1.5 Step 12: If (t maxiter) and the subgradient vector is not null then return to step 2. 3. Computational experience All the computational tests were run on a PC with a Pentium IV processor, a 2.6 GHz internal clock and 512 Mb of RAM memory. The lagrangean heuristic was written in Fortran 77 language. We tried to obtain the optimal value of the test problems in order to evaluate the quality of the upper bounds in terms of their deviation from the optimum. In this context, the MIP optimizer of the Cplex 8.0 [14]was used. The input files for Cplex were generated by means of formulation Q1 together with constraint (16) because, due to the large number of binary variables in Q2, Cplex performs poorly with this later formulation. In order to evaluate the performance of the heuristic and because no standard set of test problems is available for the SS-MCLP, we worked with a set of randomly generated test problems (which is available
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
2999
from the authors on request) already used for the MCLP, where the additional single source constraint was imposed. In the following, the main characteristics of the problem instances are briefly described (for more explanations about the generation of the test problems, the reader is referred to Correia and Captivo [1]) and the computational results are reported. 3.1. Test problems Let p be the number of possible sizes for each potential location. Five problems were generated for each of the following dimensions: n ∈ {10, 20, 50, 200, 500}, m ∈ {100, 200, 1000} and p ∈ {3, 5}. Thus, 110 test problems were obtained since only problems where n m were generated. The test problems with n = 500, m = 1000 and p = 5 have 502,500 binary variables for formulation Q1 but this number increases to 2,502,500 for formulation Q2. Therefore, large problems were used, regarding the literature in the discrete location area. In studying this kind of problem it is usual to analyze how the performance of the methods is influenced n m by the ratio between total supply and total demand. Let us consider q = . For / B a ip j i=1 j =1 i the problems described above, on average, q = 2.5. In order to evaluate the q influence, all the previous test problems were duplicated and in each data set the demand of each customer was increased by 20 units, giving rise to a mean value for q of 1.9. 3.2. Computational results After some tuning (by means of several preliminary experiments) with the objective of getting a good trade-off between solution quality and computation time, the following parameter values for the TS were used: max1 = 20, max2 = 7, max3 = 25, u1 = 7, v1 = m, u2 = 2, v2 = n and = 50%. For problems up to fifty potential locations the maxiter, N1 and N2 parameters were set to 500, 15 and 50, respectively. For problems where n=200 or n=500 the maxiter parameter was set to 2000, while the N1 and N2 parameters were set, respectively, to 50 and 150. A feasible solution (step 5 of the lagrangean heuristic) was computed whenever a different set of open facilities was obtained while resolving the relaxed problem. In all the experiments using the Cplex solver, the time limit was set to 10 h. Initially the Cplex ran with all the other parameters as default, but by making this choice and, in the case of many test problems, Cplex stopped due to memory problems before an optimal solution had been found. To try to overcome this situation, the strong branching rule to determine the best branching variable in each tree node was used (for details, see the Cplex user’s manual) and all the Cplex experiments were repeated. With this change the Cplex was able to find the optimal value of more test problems. The computational results obtained are summarized in Tables 1–4 wherein are recorded the mean gaps (in percentage) and mean execution times (in seconds) taken over 10 test problems (5 test problems with p = 3 and 5 with p = 5). For some algorithms an additional column has been included. The column “n.fails” shows the number of test problems in which the algorithm under consideration has failed the computation of a feasible solution for each set of 10 test problems. When a line of a table has “n.fails”=nf , this means, the mean gap and time in that line were computed for the remaining (10 − nf ) test problems. The gaps were computed by the expression [(U B −V )/V ∗100] where UB is the upper bound provided by each method and V denotes the optimal value, if known, otherwise it represents the maximum between the lower bound provided by RQ2 and the best lower bound at Cplex termination.
3000
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
Table 1 Results for test problems with q = 2.5 n
10 10 10 20 20 20 50 50 50 200 500
m
Lagrangean heuristic
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
Tabu search(TS)
Local search(LS)
Greedy heuristic(GH)
gap(%)
time(s)
gap(%)
time(s)
gap(%)
time(s)
n.fails
0.15 0.13 0.37 0.82 1.14 1.19 2.40 1.47 1.43 1.13 2.21 1.13
1.64 4.64 115.39 3.87 14.13 281.44 12.13 31.40 766.26 2957.93 6560.58 977.22
0.58 0.25 0.95 1.24 1.56 1.28 3.10 2.01 1.60 1.41 2.50 1.50
0.23 0.52 20.76 0.36 1.03 33.03 0.96 2.17 58.41 227.80 1666.11 182.85
1.94 1.85 2.92 2.51 3.21 1.96 4.21 2.90 2.58 1.92 3.23 2.66
0.17 0.29 1.38 0.30 0.53 2.70 0.70 1.29 6.38 105.60 274.14 35.77
0 0 0 0 0 0 0 0 0 0 6 0.55
Table 2 Results for test problems with q = 2.5 n
m
Cplex Best at Cplex termination
10 10 10 20 20 20 50 50 50 200 500
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
First integer(FI)
gap(%)
time(s)
n.fails
gap(%)
time(s)
n.fails
0.00 0.00 0.00 0.33 0.53 0.06 2.31 1.00 0.96 1.70 3.35 0.93
4.32 500.32 67.30 5686.29 14319.03 3801.72 30732.89 21431.60 35911.25 36006.36 36018.51 16770.87
0 0 0 0 0 0 0 0 0 0 5 0.45
11.15 7.65 1.98 5.36 4.67 6.55 5.33 4.08 3.81 2.59 3.54 5.16
0.11 0.58 11.70 0.32 0.69 10.47 5.73 3.11 49.14 516.96 2942.02 321.89
0 0 0 0 0 0 0 0 0 0 5 0.45
Tables 1 and 3 show the mean gaps and execution times obtained by the lagrangean heuristic when feasible solutions were computed by using the following algorithms: TS, LS and GH (in this situation, step 7 of the lagrangean heuristic was not executed). In both tables wee observe that the GH always produced the worst mean gaps while TS produced the best ones. GH was the fastest, followed by LS and lastly by TS.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
3001
Table 3 Results for test problems with q = 1.9 n
10 10 10 20 20 20 50 50 50 200 500
m
Lagrangean heuristic
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
Tabu search(TS)
Local search(LS)
Greedy heuristic(GH)
gap(%)
time(s)
gap(%)
time(s)
gap(%)
time(s)
n.fails
0.77 1.24 1.28 1.58 1.38 1.12 3.06 2.01 1.36 1.44 3.68 1.72
1.89 5.19 146.21 3.03 11.08 252.77 12.29 27.33 626.84 3157.01 6396.40 967.28
1.41 1.52 2.00 1.92 2.04 1.70 3.91 2.60 1.71 1.60 4.24 2.24
0.23 0.58 25.18 0.37 0.96 30.26 1.17 2.48 54.81 318.22 2038.87 224.83
4.06 5.48 3.60 5.08 4.37 3.20 6.33 4.30 3.92 2.30 — 4.26
0.18 0.31 1.50 0.31 0.69 3.05 0.76 1.35 6.90 115.45 — 13.05
0 0 0 0 0 0 1 0 0 6 10 1.55
Table 4 Results for test problems with q = 1.9 n
m
Cplex Best at Cplex termination
10 10 10 20 20 20 50 50 50 200 500
100 200 1000 100 200 1000 100 200 1000 1000 1000 Average
First integer(FI)
gap(%)
time(s)
n.fails
gap(%)
time(s)
n.fails
0.17 0.19 0.33 1.04 1.21 0.43 2.51 2.03 1.11 2.41 4.80 1.48
7267.54 7908.41 13842.93 29394.45 23131.32 19395.20 34043.98 17915.09 36002.52 36005.87 36018.14 23720.50
0 0 0 0 0 0 0 0 0 0 5 0.45
5.02 4.37 4.01 3.60 3.67 3.70 5.50 3.87 3.02 2.96 4.90 4.06
0.17 0.50 13.74 0.60 1.25 3.46 39.36 4.08 33.73 1024.56 5410.44 593.81
0 0 0 0 0 0 0 0 0 0 5 0.45
All the tables show that the instances with q = 1.9 were harder to solve than the ones with q = 2.5. This conclusion is not surprising since this kind of result was also reported in Cortinhal and Captivo [7] for the SS-CPLP. TS and LS found feasible solutions for all test problems but the GH failed in 6 test problems with q = 2.5 and in 17 of the instances with q = 1.9. As for the sizes of the instances, the mean execution times increase as the dimension of the problems increases, but this conclusion is not valid for the mean gaps.
3002
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
Tables 2 and 4 report the experiments with Cplex. As we pointed out above, for each instance, Cplex stopped operating when one of the following conditions hold: an optimal solution has been found within 10 h, the time limit of 10 h was attained or an out of memory problem occurs. Because Cplex sometimes did not find an optimal solution, we present the columns with the label “Best at Cplex termination”. The other columns show the results computed with the First Integer (FI) solution provided by Cplex. The computation of the FI by Cplex did not seem a good alternative to the lagrangean heuristic. In fact, FI is not competitive since with LS, all instances were successfully solved and it gave a better mean gap and time. The average gaps and times provided by FI were worst than the ones provided by GH although FI was able to solve one more test problem. In relation to the TS, it was found that TS clearly dominates FI in terms of the average gaps and number of instances unsolved but TS proved to be more time consuming. Regarding the attempts to optimally solve the test problems, Table 4 shows that even in the case of the problems with n = 10 the gap’s column has nonzero entries because Cplex was unable to optimally solve those instances. Despite the fact that average gaps presented in the column “Best at Cplex termination” are the best, they are not much lower than the ones provided by TS and TS was able to solve all the problems with significantly less computational effort. 4. Concluding remarks This paper addresses a Single Source Capacitated Plant Location Problem with several possible capacity levels for the facility that can be opened at each potential location. Because this is a hard problem to solve a lagrangean heuristic was developed. We propose a Tabu Search and a Local Search algorithm to improve the feasible solutions during the lagrangean heuristic. In order to check the efficiency of the methods proposed in this paper, some test problems were randomly generated and some experiments were performed using the Cplex solver. The heuristic methods provide very satisfactory results regarding the mean gaps and execution times even for large scale problems. As is usual with this kind of problem, instances with tighter capacity constraints were harder to solve. Future research may focus on the development of different meta-heuristic procedures to try to improve the quality of the feasible solutions. Acknowledgements The authors thank the referees for their helpful comments. References [1] Correia I, Captivo ME. A Lagrangean heuristics for a modular capacitated location problem. Annals of Operations Research 2003;122:141–61. [2] Holmberg K. Solving the staircase cost facility location problem with decomposition and piecewise linearization. European Journal of Operational Research 1994;75:41–61. [3] Holmberg K, Ling J. A lagrangean heuristic for the facility location problem with staircase costs. European Journal of Operational Research 1997;97:63–74.
I. Correia, M. Eugénia Captivo / Computers & Operations Research 33 (2006) 2991 – 3003
3003
[4] Harkness J, Revelle C. A Lagrangean heuristic for the staircase cost facility location problem. International Journal of Industrial Engineering-Theory, Applications and Practice 2002;9(1):36–44. [5] Beasley J. An algorithm for solving large capacitated warehouse location problems. European Journal of Operational Research 1988;33:314–25. [6] Beasley J. Lagrangean heuristics for location problems. European Journal of Operational Research 1993;65:383–99. [7] Cortinhal MJ, Captivo ME. Upper and lower bounds for the single source capacitated location problem. European Journal of Operational Research 2003;151(2):333–51. [8] Glover F, Laguna M. Tabu search. Dordrecht: Kluwer Academic Publishers; 1997. [9] Klincewicz J. Avoiding local optima in the p-hub location problem using Tabu Search and GRASP. Annals of Operations Research 1992;40:283–302. [10] Rolland E, Schilling D, Current J. An efficient Tabu Search procedure for the p-median problem. European Journal of Operational Research 1997;96:329–42. [11] Geoffrion A. Lagrangean relaxation for integer programming. Mathematical Programming Study 1974;2:82–114. [12] Held M, Wolfe P, Crowder H. Validation of subgradient optimization. Mathematical Programming 1974;6:62–6. [13] Tuzun D, Burke L. A two-phase tabu search approach to the location routing problem. European Journal of Operational Research 1999;116:87–99. [14] Cplex. Using the Cplex callable library, Version 8.0, ILOG, Inc., 2002.