A goal programming approach to vehicle routing problems with soft time windows

European Journal of Operational Research 177 (2007) 1720–1733 www.elsevier.com/locate/ejor

A goal programming approach to vehicle routing problems with soft time windows q Herminia I. Calvete a,*, Carmen Gale´ b, Marı´a-Jose´ Oliveros c, Bele´n Sa´nchez-Valverde b a

Dpto. de Me´todos Estadı´sticos, Universidad de Zaragoza, Pedro Cerbuna 12, 50009 Zaragoza, Spain Dpto. de Me´todos Estadı´sticos, Universidad de Zaragoza, Marı´a de Luna 3, 50018 Zaragoza, Spain Dpto. de Ingenierı´a de Disen˜o y Fabricacio´n, Universidad de Zaragoza, Marı´a de Luna 3, 50018 Zaragoza, Spain b

c

Available online 18 November 2005

Abstract The classical vehicle routing problem involves designing a set of routes for a fleet of vehicles based at one central depot that is required to serve a number of geographically dispersed customers, while minimizing the total travel distance or the total distribution cost. Each route originates and terminates at the central depot and customers demands are known. In many practical distribution problems, besides a hard time window associated with each customer, defining a time interval in which the customer should be served, managers establish multiple objectives to be considered, like avoiding underutilization of labor and vehicle capacity, while meeting the preferences of customers regarding the time of the day in which they would like to be served (soft time windows). This work investigates the use of goal programming to model these problems. To solve the model, an enumeration-followed-by-optimization approach is proposed which first computes feasible routes and then selects the set of best ones. Computational results show that this approach is adequate for medium-sized delivery problems.  2005 Elsevier B.V. All rights reserved. Keywords: Vehicle routing; Logistics; Time windows; Goal programming; Enumeration

1. Introduction In such a competitive world as the current one, it is evident that companies should make operational and strategic decisions in order to optimize and manage the processes involved in their supply chains more effiq

This research work has been supported by Spanish Contract IBE2002-CIEN-03. Corresponding author. Fax: +34 976 761115. E-mail addresses: [email protected] (H.I. Calvete), [email protected] (C. Gale´), [email protected] (M.-J. Oliveros), [email protected] (B. Sa´nchez-Valverde). *

0377-2217/$ - see front matter  2005 Elsevier B.V. All rights reserved. doi:10.1016/j.ejor.2005.10.010

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1721

ciently. Amongst them, the delivery of commodities at low cost, with high quality service. This guarantees not only a good service to customers but also savings in warehousing and distribution costs. The problem of physical distribution of goods to customer locations is specially relevant since it accounts for a large proportion of the overall operational costs of a producer. Hence, effective and efficient management of transportation and distribution of goods is becoming increasingly important both from the point of view of theoretical research and from the point of view of practical applications. These problems are generically known as Vehicle Routing Problems (VRP). The classical VRP consists of determining the best set of routes for a fleet of vehicles originating and terminating at a single central depot to distribute goods to a set of customers geographically dispersed, while minimizing the total travel distance/time or the total distribution cost. In addition, some problems have pre-set time constraints on the periods of the day in which deliveries should take place. They are known as Vehicle Routing Problems with Time Windows (VRPTW). In these problems customers have time windows associated which are characterized by the earliest and the latest allowable times within which the service should begin. The goal is to determine the optimal set of routes and the optimal sequence of customers visited by each vehicle, taking into account vehicle capacities, service times and time windows. These time windows are hard constraints when a route is not feasible if the service of any customer does not start within the limits established by the time window. These are Vehicle Routing Problems with Hard Time Windows (VRPHTW). In other cases, both lower and upper bounds of the time window need not be met, but can be violated at a penalty. These are Vehicle Routing Problems with Soft Time Windows (VRPSTW). Due to its importance, there is a vast literature that addresses modeling and solution aspects of both VRP and VRPTW. Comprehensive reviews of these problems are presented in [2,6,10,11,24,28]. All these combinatorial optimization problems have been proven to be NP-hard and only relatively small instances can be solved to optimality due to their huge computational requirements. Methodology involves the use of branch and bound methods [3,15], column generation [8] and Lagrangian relaxation [14]. For larger problems, researchers usually focus on heuristic and meta-heuristic (tabu search, genetic algorithms, simulated annealing) methods to derive approximate solutions of acceptable quality in reasonable computational time [7,12,13,23]. The number of references in VRPSTW is relatively small. Sexton and Choi [22] consider a single-vehicle pickup-and-delivery routing problem where each load has its own origin and its own destination. They minimize a linear combination of total vehicle operating time and total customer penalty due to missing any of the time windows, and solve the problem by BendersÕ decomposition. Min [17] considers a single-vehicle library distribution system with a single depot where the daily distribution activity starts and ends. He incorporates customerÕs specific time window preferences into the modeling process and propose a goal approach to solve the problem. The objectives are to minimize the total travel time and the deviations from time window preferences. Koskosidis et al. [16] propose an heuristic approach for the VRPSTW which iterates between a generalized assignment problem to assign customers to vehicles and a series of routing and scheduling problems defined by traveling salesman problems with soft time windows constraints. Balakrishnan [1] proposes three heuristics for the VRPSTW with a homogeneous fleet based on the nearest neighbor, Clarke-Wright savings and space–time rules. The cost of violating time windows is assumed to be a linear function of the amount of time window violation. Taillard et al. [26] propose a tabu search heuristic to solve a VRPSTW in which total distance traveled plus total penalty incurred for not serving in the time window is minimized. Fagerholt [9] deals with a multi-ship pickup and delivery problem with soft time windows in which inconvenience costs for serving customers outside their time windows are imposed. The sum of transportation costs and inconvenience costs is minimized by solving a set partitioning problem whose variables are feasible schedules. On the other hand, it is worth pointing out that most real-life delivery problems have other objectives in addition to minimizing the total travel distance/time or the total distribution cost. For instance, to avoid underutilization of labor or vehicle capacity, while meeting the preferences of customers regarding the time

1722

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

of the day in which they would like to be served. This results in a multiobjective VRP. Past years, great attention has been paid to the development of methodologies for dealing with multiobjective programming problems [5,25,29]. Amongst the most popular of these is the goal programming (GP) approach [19,21,27]. Underlying GP is a concept of satisfying of objectives, i.e. to seek a solution that comes as close as possible to the goals. For this purpose a target level is established which represents an acceptable level of achievement for each objective. Deviational variables are defined in order to represent the amount by which we numerically exceed or are under the corresponding target. Usually, GP models are classified into two different categories [27]. In weighted goal programming, weights are assigned to deviations from targets levels which reflect the decision makerÕs rating of the importance of these deviations. Then, the total penalty is minimized. In lexicographic or preemptive goal programming priorities are assigned to goals. Next, goals are minimized one at a time in order of priority. At each successive minimization, higher priority level optimum values are maintained as determined in previous minimizations. In this paper, we consider a general medium-sized delivery problem with soft and hard time windows, a heterogeneous fleet of vehicles and multiple objectives. Its purpose is to investigate the use of goal programming to model the problem and propose an approach to solve it providing an optimal solution in a reasonable computational time in terms of real applications. As we will see in the paper, on the one hand, the inclusion of multiple objectives and specific constraints in the general setting of VRP makes the model more complicated and highlights the difficulties regarding the time needed to solve it optimally with standard optimization software. Nevertheless, on the other hand, these features can also be exploited to come up with solving procedures that, at first sight, might not seem very efficient. Keeping this in mind, we propose not to directly solve the resulting goal programming model with a general purpose software, but to carry out a two-step approach that we call Ôenumeration-followed-by-optimization (EFBO)Õ approach. This procedure consists of firstly enumerating all feasible routes and computing its deviations from targets. Secondly, the best set of routes according to the goals previously established are selected by solving a set partitioning problem. A similar approach has been taken in [9] to deal with a multi-ship pickup and delivery problem with soft time windows. Also this approach has been adopted in [20] to get an optimal solution to a real feed compound distribution problem. The special characteristics of vehicles, which are compartmentalized, reduces the number of feasible routes. As noted in both papers, the efficiency of the procedure lies on the number of feasible routes which have to be dealt with. The remainder of the paper is organized as follows. We start in Section 2 by presenting a comprehensive description of the VRPTW tackled and formulating a goal programming model. Section 3 goes on developing the EFBO approach. Finally, in Section 4 an experiment is carried out to evaluate the computational efficiency of the procedure and we conclude with a summary.

2. Model formulation We consider the problem of daily delivery of goods from a single central depot to a number of geographically dispersed customers. The aim is to design the routes for serving customers, i.e. to determine the ordered sequence of locations on each vehicle route and the arrival and departure times for each customer location. Let G ¼ ½N; A be a directed network associated to the system, where N ¼ f1; . . . ; ng is the set of nodes (each representing the central depot or a customer location) and A ¼ fði; jÞ : i; j 2 Ng is the set of directed arcs (each representing a direct connection). Index 1 refers to the central depot, while indices 2 to n refer to the customers. We assume that customers demands are known when routes are scheduled. Let qi denote the demand of customer i. Moreover, customer i requires a service time si. Notice that, usually, service time is a function of demand composed of a fixed term sfi and a variable term svi , which in general is assumed to be a linear function of demand. Hence, si ¼ si ðqi Þ ¼ sfi þ svi ðqi Þ.

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1723

Because of town council regulations or requirements of customers, each customerÕs delivery should begin and finish (i.e. has to take place) within a pre-specified time window. Let ½ehi ; lhi  define the hard time window associated with customer i. The lower bound ehi indicates the earliest time in which the service to customer i should start. Similarly, the upper bound lhi indicates the latest time in which the service to customer i should end. In addition to the hard time window, there exists a soft time window ½esi ; lsi  which lays out the preferences of the customer regarding the period of the day within which he/she would like to be served. Certainly, ehi 6 esi and lsi 6 lhi . Travel time between the central depot and each customer location, as well as between each pair of customer locations is also known. Let tij denote the time taken to travel from node i to node j through arc (i, j). We assume that delivery starts as soon as the vehicle arrives at the customer location and that the vehicle leaves the customer as soon as the delivery has been completed (once the service starts, it is carried out until completion), that is to say vehicles are not permitted to wait to begin the service when they arrive too early at a customer location. A heterogeneous fleet of vehicles is available. They are located at the depot, so they must leave from and return to the central depot. We assume that there is no limitation on the number of vehicles of each type that can be used, but in order to facilitate the model formulation we denote by V the maximum possible size of the fleet. The actual number of vehicles that will be used is one of the variables of the model. Let Uk denote the fixed capacity of vehicle k. We assume that Uk is such that at least one customer may be assigned a vehicle. Vehicles are assumed not to be overloaded, hence total demand of all customers on one particular route must not exceed the capacity of the vehicle assigned to this route. Moreover, deliveries cannot be split up amongst vehicles, that is to say, each customer is served by a single vehicle. Let ck and ckij , respectively, denote the activation cost and the cost of traveling from i to j, associated with vehicle k. Certainly, the capacity and cost of all the same type of vehicles are equal. In order to formulate the model we define the following variables:  1 if vehicle k travels directly from node i to node j; k xij ¼ ði; jÞ 2 A. 0 otherwise;  1 if vehicle k is actually used; zk ¼ k ¼ 1; . . . ; V . 0 otherwise; ai = arrival time (service starting) to customer i. pi = departure time (service ending) from customer i. pk1 ¼ departure time from the central depot of vehicle k. In most distribution companies, vehicles and drivers are not waiting to leave the depot at any moment of the day. Generally, there exists a discrete set of departure times at which a vehicle can begin its assigned route. In order to simplify notations, we assume that all vehicles have the same possible set of departure times, {d1, . . . , dm}. Hence, pk1 2 fd 1 ; . . . ; d m g; k ¼ 1; . . . ; V . Since delivery starts as soon as the vehicle arrives to the customer, pi = ai + si, i = 2, . . . , n. In the classical VRPHTW the single objective is to minimize operational costs while holding hard time window constraints. As mentioned in the previous section, managers usually have more objectives to be considered than just cost, like avoiding underutilization of labor and vehicles capacity, while meeting the preferences of customers regarding the time of the day in which they would like to be served (soft time windows). Hence, to formulate a goal programming model, we set the following goals: • • • •

Goal Goal Goal Goal

1: 2: 3: 4:

Minimize total operational cost Satisfy time window preferences of customers (soft time windows) Avoid underutilization of vehicles capacity Avoid underutilization of labor

1724

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

and establish target levels which represent acceptable levels of achievement for each objective. Let Z be a lower bound on the total operational cost of delivery. It can be obtained, for instance, by solving a classical VRP, that is to say, all time window constraints are removed and only operational costs are minimized. Let tmax stand for the total daily time that a driver can be working due to labor regulations. Let ~tmax denote the total daily time that a driver can be driving due to legal regulations. Given the above-defined goals and variables, the problem can be formulated as follows: MIGP: min

€ g€þ þ W

n X

bi g^ þ W i

n X

i¼2

s.t.

V X X

ckij xkij þ

V X

n X

qi

i¼2

X

k¼1

V X

k g ; W k

ð1Þ

k¼1

ð2Þ

k¼1



n X

W_ k g_  k þ

ck zk  g€þ ¼ Z;

s ^þ pi þ g^ i g i ¼ e i þ si ;

pi þ

V X

i¼2

k¼1 ði;jÞ2A

g~ i

ei g~þ þ W i

g~þ i

¼

lsi ;

i ¼ 2; . . . ; n;

i ¼ 2; . . . ; n;

xkij  U k zk þ g_  k ¼ 0;

ð3Þ ð4Þ

k ¼ 1; . . . ; V ;

ð5Þ

j¼1

ðsi þ tij Þxkij  tmax zk þ g k ¼ 0;

k ¼ 1; . . . ; V ;

ð6Þ

ði;jÞ2A n X V X i¼1

n X V X i¼1

xkij ¼ 1;

j ¼ 2; . . . ; n;

ð7Þ

xkji ¼ 1;

j ¼ 2; . . . ; n;

ð8Þ

k¼1

k¼1

xkij  zk 6 0; ði; jÞ 2 A; k ¼ 1; . . . ; V ; X zk  xkij 6 0; k ¼ 1; . . . ; V ;

ð9Þ ð10Þ

ði;jÞ2A n X

xk1j 6 1;

k ¼ 1; . . . ; V ;

ð11Þ

xki1 6 1;

k ¼ 1; . . . ; V ;

ð12Þ

j¼2 n X i¼2 n X i¼1

xkij 

n X

xkji ¼ 0;

j ¼ 1; . . . ; n; k ¼ 1; . . . ; V ;

ð13Þ

i¼1

pj  pk1 þ ð1  xk1j ÞM P sj þ t1j ;

j ¼ 2; . . . ; n; k ¼ 1; . . . ; V ;

ð14Þ

pj  pk1  ð1  xk1j ÞM 6 sj þ t1j ;

j ¼ 2; . . . ; n; k ¼ 1; . . . ; V ;

ð15Þ

pj  pi þ ð1  xkij ÞM P sj þ tij ;

i; j ¼ 2; . . . ; n; k ¼ 1; . . . ; V ;

ð16Þ

i; j ¼ 2; . . . ; n; k ¼ 1; . . . ; V ;

ð17Þ

pj  pi  ð1  pi P ehi þ si ;

xkij ÞM

6 sj þ tij ;

i ¼ 2; . . . ; n;

ð18Þ

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

pi 6 lhi ; i ¼ 2; . . . ; n; X tij xkij 6 ~tmax ; k ¼ 1; . . . ; V ;

1725

ð19Þ ð20Þ

ði;jÞ2A

xkij 2 f0; 1g; zk 2 f0; 1g; pk1 2 fd 1 ; . . . ; d m g; pi P 0;

i; j ¼ 1; . . . ; n; k ¼ 1; . . . ; V ;

k ¼ 1; . . . ; V ;

i ¼ 2; . . . ; n;

^þ ~ ~þ _  g€þ ; g^ i ;g i ;g k ;g k P 0; i ;g i ;g

i ¼ 2; . . . ; n; k ¼ 1; . . . ; V ;

~þ where M is a large number; g€þ , g^þ i , and g i are positive deviational variables representing the amount by ~ _  which we numerically exceed the corresponding goal; g^ i , g i , g k , and g k are negative deviational variables €, W b i, W e i, representing the amount by which we are numerically under the corresponding goal. Finally, W _  W k and W k denote the penalties per unit of deviation from each goal. That is to say, they are the relative weights that establish the importance of goals to which they are attached. Constraint (2) refers to goal 1. Taking into account that Z is a lower bound on the total cost, only a € in the objective function (1). Similarly positive deviational variable g€þ is included, which is weighted by W a different target cost can be used, but in this case both positive and negative deviational variables should be included. Constraint (3) refers to goal 2. It indicates that the arrival time to customer i, ai = pi  si, can undersatisfy or oversatisfy the arrival time preference of the customer, esi . A similar statement can be made with constraint (4) and the departure time. In the objective function, weights are assigned to undersatisfy the arrival time preference and to oversatisfy the departure time preference, that is to say, penalties are assigned s to either arrivals before the preferred earliest time esi or departures Pn after Pn thek preferred latest time li . Constraint (5) allows us to formulate goal 3. The expression i¼2 qi j¼1 xP ij gives the total load allocated to vehicle k. Similarly, constraint (6) concerns goal 4. By assuming s1 = 0, ði;jÞ2A ðsi þ tij Þxkij measures the total time the driver of vehicle k is working. Since we assume that both capacity Uk and maximum daily time working tmax cannot be exceeded, only the positive deviational variables have been added up. If either Uk or tmax can be exceeded, both positive and negative deviational variables should be included and weighted accordingly in the objective function. The remaining constraints are hard constraints which allow us to express the essential features of VRP. Constraints (7) and (8) ensure that only a vehicle arrives to and leaves from customer j. Constraints (9) and (10) ensure that just used vehicles are activated. Constraints (11) and (12) impose that all vehicles leave from and return to the depot. Constraint (13) is the typical flow conservation equation that ensures the continuity of each vehicle route. Constraints (14)–(17) guarantee the feasibility of the schedule for each vehicle. Constraints (18) and (19) refer to hard time windows of customers. Finally, (20) ensure that the maximum driving time is not exceeded. The resulting model is a mixed-integer goal programming problem, which can be solved by applying usual techniques [18,21]. Nevertheless, it is worth noting that the computational time needed to exactly solve this kind of problem strongly depends on the number of 0–1 variables involved and is prohibitive from a practical point of view except for very small-sized models. When all vehicles are identical this model can be simplified, since it is not necessary to identify the vehicle which is going from the customer i to the customer j through the connection (i, j). Thus, the number of 0–1 variables decreases since it is enough to define a single 0–1 variable xij, indicating if a vehicle either travels from i to j or not. This reduction makes it possible to solve the model with a larger number of customers. A simplified version of this model was studied in [4] and, as was expected, only small-sized models with up to twelve customers were able to be solved optimally with standard optimization software in a reasonable amount of computational time. Moreover,

1726

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

the computational time needed to get an optimal solution increases very quickly with the number of customers. Clearly, this is not very appealing when using the model to handle real systems.

3. Solution methodology The problem formulation MIGP emphasizes vehicles and connections between nodes since it determines the value of the 0–1 variables xkij . Since we are interested in a set of routes we can think of reformulating the problem taken into account in terms of routes (with its associated vehicle) rather than in terms of connections between nodes. At first sight this methodology could seem not very efficient since it requires a previous computation of the set of all feasible routes. A priori, the number of feasible routes seems to be huge and it gives the impression that the enumeration of this set should be used only as a last resort for solving the problem. Nevertheless, one of the principal reasons that motivated us to propose the enumeration-followed-by optimization methodology for solving the problem is the realization that in most real medium-sized VRPTW a simple visual examination of the problem allows us to conclude that, mainly owing to distances between nodes, hard time windows, labor and driving regulations and capacity of vehicles, there is a great number of possible routes which actually are not feasible, but it is impossible to ÔreportÕ this circumstance to the model MIGP. On the other hand, the more realistic the description of a problem is, the more complicated the formal model MIGP results, but, in general, there are less feasible routes. Hence, to solve the VRPSTW considered we propose the methodology EFBO which consists of two phases. First phase enumerates all feasible routes. Since there are time windows, a route should specify the customers which are visited, the order in which they are visited, the vehicle which is used and the departure time of the vehicle from the central depot. Moreover, according to the goal programming approach, for each feasible route R, its total penalty incurred due to deviations from targets, FR, is computed. The second phase of the methodology selects, from the set of feasible routes R, the subset of best ones. For this purpose, let us define for each R 2 R the following 0–1 variable:  1 if route R is selected; XR ¼ 0 otherwise. Hence, in order to determine the routes that should be used for delivering the following set partitioning problem is solved: SP: min

X

F RX R

R2R

s.t.

X

ARi X R ¼ 1;

i ¼ 2; . . . ; n;

ð21Þ

R2R

X R 2 f0; 1g

8R 2 R;

where  ARi ¼

1

if route R visits customer i;

0

otherwise;

and constraint (21) ensures that all customers are visited exactly once. Notice that the number of variables in problem SP will usually be much larger than in problem MIGP (one variable for every feasible route). However, the number of constraints is much less (one constraint for every customer). Despite of the initial impression, we will show in Section 4 that the EFBO methodology yields good results in terms of the total

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1727

computational time involved in both phases and real medium-sized delivery problems can be well managed with currently used hardware and software. Now we are going to look at the first phase. 3.1. Phase 1: Computing feasible routes The purpose of this phase is to enumerate all possible ways of connecting the depot and the customers in order to identify which of them provide feasible routes. A route involves information on the customers which are visited, the order in which they are visited, the vehicle which is used and the departure time of the vehicle from the central depot. Let {1, i2, i3, . . . , ih, 1} be the sequence of ordered customers of a route R. Regarding working time, for the route R to be feasible, (22) has to hold true, irrespective of which vehicle serves the route. A similar statement can be made with constraint (23) and driving time. tR ¼ t1i2 þ si2 þ ti2 i3 þ si3 þ    þ sih þ tih 1 6 tmax t1i2 þ ti2 i3 þ    þ tih 1 6 ~tmax

ð22Þ ð23Þ

All possible routes not verifying both constraints are not feasible, hence can be rejected. Moreover, according to goal 4, the expression tmax  tR measures the underutilization of the route R driverÕs working time. Concerning time windows it is worth noting that either the possibility of holding a hard time window constraint or not and the total deviation from soft time windows depend only on the vehicle departure time from the central depot. Taking into account the discrete set of departure times established {d1, . . . , dm}, for each du 2 {d1, . . . , dm}, arrival time to, aij ðd u Þ, and departure time from, pij ðd u Þ, customer ij 2 R are calculated and hard time window constraints are checked. Then, if there are several departure times which verify these constraints, the departure time from the depot dR is selected as the one which minimizes deviations from soft time windows, i.e., it is obtained by solving the following problem: h h X X b ij maxf0; es  aij ðd u Þg þ e ij maxf0; pi ðd u Þ  ls g W W min ij ij j du

s.t.

j¼2

j¼2

ehij ;

j ¼ 2; . . . ; h;

pij ðd u Þ 6 lhij ;

j ¼ 2; . . . ; h.

aij ðd u Þ P

Once decided the departure time dR, the expression (24) measures, according to goal 2, the weighted deviation of route R from customer preferences h X j¼2

b ij maxf0; es  aij ðd R Þg þ W ij

h X

e ij maxf0; pi ðd R Þ  ls g. W ij j

ð24Þ

j¼2

Next, we will decide which vehicle serves the route R. In fact, we will decide which type of vehicle serves the route. It is not necessary to distinguish between vehicles of the same type because they are identical. For this purpose, we consider total load and vehicle capacity, that is to say, we examine which types of vehicle are able to serve all customers of the route (owing to characteristics of the customer location some vehicles may not be able to gain access to it) and carry the total demand of the route. Assume for the time being that there are three different types of vehicles and all of them can be used in route R. In order to simplify notations we consider that vehicle 1 is one of the first type, vehicle 2 is one of the second type and vehicle 3 is one of the third type. The total load of R is QR ¼ qi2 þ qi3 þ    þ qih . Since all the three vehicles can serve the route, QR 6 U1, QR 6 U2 and QR 6 U3. Moreover, the expression Uk  QR, k = 1, 2, 3, measures the underutilization of vehicle k capacity.

1728

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

Let cR(k) be the cost of using the vehicle k in route R, i.e. X cR ðkÞ ¼ ck þ ckij ; k ¼ 1; 2; 3.

ð25Þ

ði;jÞ2R

Since goal 1 is to minimize total operational cost and goal 3 is to avoid underutilization of vehicle capacity, we should select the vehicle kR for which the minimum in (26) is achieved, min fWcR ðkÞ þ W_ k ðU k  QR Þg;

ð26Þ

k¼1;2;3

where W stands for the penalty associated with route cost. That is to say, we associate with route R Ôthe most economicalÕ vehicle, according to cost and underutilization of vehicle capacity. Note that this approach saves us solving a VRP problem to compute the target Z that appeared in the MIGP formulation since operational cost is directly minimized. Besides, one of the advantages of the EFBO approach is that (25) can have as complicated an expression as desired since it will not affect the optimization problem to be solved in the second phase. For instance, it could be more realistic to consider that operational cost is computed by adding up three terms: an activation cost, a cost by distance traveled and a cost by load, computed as the maximum between a minimum load cost and the product of total load by cost per unit loaded. This type of cost really complicates the formulation MIGP, but can be easily managed with the EFBO approach. Hence, upon termination of the first phase, the set R of all feasible routes will be available. Each feasible route R will be characterized by an ordered sequence of customers {1, i2, i3, . . . , ih, 1}, a vehicle kR and a departure time from the central depot dR. These routes are the variables in the problem SP and, from the previous remarks, we conclude that their objective function coefficient in this problem is obtained as indicated in the following expression: h h X X b ij maxf0; es  aij ðd R Þg þ e ij maxf0; pi ðd R Þ  ls g W W F R ¼ WcR ðk R Þ þ ij ij j j¼2

 kR ðtmax  tR Þ. þ W_ kR ðU kR  QR Þ þ W

j¼2

ð27Þ

3.2. Filtering feasible routes Once the first phase has finished providing the set of feasible routes R, we have to solve the SP problem in the second phase of the procedure using Lingo or another optimization software. For this purpose, the number of variables is probably the main parameter to be taken into account. Hence, we can conclude that the efficiency of the proposed procedure lies on the number of feasible routes to be dealt with. Besides the time needed to compute them, they determine the ability of Lingo to solve the set partitioning problem in which these feasible routes are the variables. Therefore, we have improved the procedure by considering an additional step, that we have named ÔfilteringÕ. Its purpose is to reduce the number of routes provided as variables to the SP problem by removing routes which will not be in the set of optimal routes to this problem. Certainly this reduction must be done guaranteeing that we do not leave any promising feasible route behind. In the worst case, no route would be discarded. But, in general, this step will lead to a substantial reduction in the number of variables of the SP problem, as we will show in the next section. The primary idea in the filtering step is to select from the set of feasible routes involving a subset of nodes only the most economical according to (27). Let us consider the customers j1, j2 and j3. Assume for the time being that the feasible routes involving all these nodes are R1 = 1  j2  j1  j3  1, R2 = 1  j2  j3  j1  1 and R3 = 1  j3  j2  j1  1, and that F R1 6 F R2 6 F R3 . If in the optimal solution a vehicle must visit the customers j1, j2 and j3, it is clear that the route selected will be R1 since it has the less cost F. Hence, routes R2 and R3 can be discarded without affecting the optimal solution of the SP problem.

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1729

In general, given the set R obtained at the end of the first phase, the filter step analyzes all the routes involving each subset of nodes and selects the best permutation according to (27). At the end of the process e is constructed which contains all routes which are candidates for the optimal a new set of feasible routes R solution. In order to determine the routes that should be used for delivering it is enough to solve a set parf , formulated in a similar way to the SP problem, but considering the new set of feasible titioning problem SP e Notice that, in general, this new problem will have less variables, so will be easier to be solved. routes R.

4. Experiment results A computational experiment has been conducted to analyze the performance of the proposed methodology. Since it provides an optimal solution to the problem analyzed, with this experiment we aim to show our method efficiency in terms of computing time for managers who have to make regular decisions quickly. Our methodology has been tested on two different sets of problems. The first one is inspired by the delivery problem which motivated the development of the procedure. The second one is based on standard problems found in [23,1]. The experiments were performed on a PC Pentium 4 at 3.0 GHz having 2 GB of RAM under Windows XP. The set partitioning problem of the second phase was solved by using Lingo 8.0. In the first set of problems we consider a central depot serving a geographical zone within a distance radius of approximately 200 km. There are up to 69 customers and the maximum time–distance between customers is about 300 minutes. For the customer i, the service time fixed term is 15 minutes and the variable term is twice the demand in minutes. The time windows reflect town council regulations or customer requirements for delivery early in the morning, early in the afternoon or in the evening. The soft time windows are from 7 a.m. to 11 a.m., from 2 p.m. to 16.30 p.m. and from 5 p.m. to 9 p.m. The hard time windows permit a maximum deviation from the soft time windows of one hour on each side. Hence, we have clustered customers into four groups (see Table 1). Taking into account travel times and logical constraints, we have applied the incompatibilities shown in Table 1 when building feasible routes. The methodology has been tested on 60 problems grouped into six different configurations with 30, 50 and 70 nodes. The characteristics of these configurations are summarized in Table 2. They were selected with differences in the number of customers and in the percentage of customers pertaining to group (4). In previous experiences, we found that this percentage had a strong influence on the number of feasible routes at the end of the first phase, because members of this group accept being served at any moment of the day. For each configuration, ten problems were generated which differ in customer demands. They are randomly selected within the range shown in Table 2. Customer demands also have a great impact on the number of feasible routes since they influence the route length. For all problems there are two different types of vehicle, with capacity 14 and 30 load units, respectively. Fourteen possible departure times from the depot were taken into consideration, from 6 a.m. to 19 p.m., every hour. The total daily time that a driver can be working is 600 minutes and the total daily time that a driver can be driving is 480 minutes. The results of this experiment are displayed in Tables 3 and 4. The first one depicts the minimum, maximum, average and standard deviation of the number of feasible routes obtained after applying the first Table 1 Compatibilities between customers

(1) (2) (3) (4)

A customer from group

Can only be followed by a customer of groups

Service Service Service Service

(1), (2), (3), (1),

from from from from

6 1 6 6

a.m. to 12 a.m. p.m. to 17.30 p.m. p.m. to 10 p.m. a.m. to 10 p.m.

(2), (4) (3), (4) (4) (2), (3), (4)

1730

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

Table 2 Configurations for test problems

Configuration Configuration Configuration Configuration Configuration Configuration

1 2 3 4 5 6

# of customers

% of customers from group (4)

Demand range

29 29 49 49 69 69

0 20 0 20 0 20

[2, [2, [2, [2, [2, [2,

15] 15] 8] 8] 15] 15]

Table 3 First phase results from 10 instances in each configuration # of feasible routes

Configuration Configuration Configuration Configuration Configuration Configuration

1 2 3 4 5 6

First phase time (in seconds)

Min

Max

Average

Std. Dev.

Min

Max

Average

Std. Dev.

1069 2001 6314 19,285 10,024 31,124

1707 4952 9684 31,961 20,133 87,334

1384.4 3299.3 7805.8 25,475.8 15,822.0 46,642.4

186.5 942.1 1215.8 4062.0 3005.1 17,407.6

0 1 2 7 5 13

1 2 4 12 9 37

0.9 1.5 3 10.2 6.9 19.6

0.3 0.5 0.8 1.5 1.6 7.6

Table 4 Second phase results from 10 instances in each configuration using Lingo Second phase time (in seconds)

Configuration Configuration Configuration Configuration Configuration Configuration

1 2 3 4 5 6

Min

Max

Average

Std. Dev.

2 4 21 62 74 198

2 7 35 113 133 564

2 4.8 27.7 85.5 108.2 310.5

0 1.0 5.2 15.2 20.1 126.3

phase without filtering to 10 instances in each configuration. Likewise it shows the same data for the time (in seconds) required to compute all feasible routes. Table 4 shows the time taken by Lingo to solve the corresponding set partitioning problem of the second phase (minimum, maximum, average and standard deviation). As mentioned above, the number of ÔflexibleÕ customers, i.e. customers accepting being served at any moment of the day substantially increases the number of feasible routes. When there are no customers from group (4), problems with 29 customers can be completely solved, in average, within 2.9 seconds, problems with 49 customers within 30.7 seconds and problems with 69 customers within 115.1 seconds. When considering a 20% of customers from that group, these times are 6.3, 95.7 and 330.1, respectively. The second set of problems in which we have tested the procedure are based on Balakrishnan soft time window problems [1]. These problems are derived from data set R1 of Solomon test problems. For them, hard time windows are converted to soft time windows by allowing a certain percentage of time window violation Pmax, expressed as a percentage of the maximum route time allowed. For all problems in original data set R1, the customer locations, demands and time windows are randomly generated. The percentage of customers with time windows ranges from 100% to 25%. The float is homogeneous and each vehicle capac-

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1731

Table 5 Results from problems R101 to R112 with 25 customers (time in seconds) First phase

R101 R102 R103 R104 R105 R106 R107 R108 R109 R110 R111 R112

Filtering

Second phase

# of feasible routes

Time

# of feasible routes

Time

Time

1771 223,301 1,220,236 5,119,022 18,595 604,599 2,131,244 6,874,800 154,298 452,991 1,347,641 2,327,057

1 88 477 2006 8 229 834 2694 61 179 527 917

1744 33,829 92,961 205,353 12,813 74,303 146,309 260,256 52,780 93,193 125,274 217,574

1 14 70 332 2 34 132 473 9 25 85 149

2 38 90 230 11 80 160 297 50 89 121 243

ity is 200 load units. The maximum route time is 230 time units. We have applied the idea of Balakrishnan [1] to the data set R1 of SolomonÕs 25 customers benchmarks problems used by Desrochers et al. [8] and we have solved all twelve problems in this set. As indicated by Balakrishnan we took Pmax = 5%. Moreover, we selected sixteen different departure times from the depot, every fifteen time units, and there is no waiting time at any customer. Previously, we had analyzed problems with 100 and 50 customers, but the large number of feasible routes made it impossible to solve all of them, except the problem R101 with 50 customers, owing to Lingo running out of memory. It is worth noting that, as mentioned throughout the paper, our procedure is well designed to solve medium-sized problems with constraints that indeed restrict the number of feasible routes. The results of this experiment are shown in Table 5. Notice that only problems R101, R102, R105 and R109 could have been solved if filtering had not been applied because of Lingo running out of memory. In fact, this additional step achieves reductions of up to 96% at a reasonable cost in terms of computational time. As was expected, the ability to eliminate routes increased when the percentage of customers with time windows decreased. Problems with a number of feasible routes from 1700 to 150,000 have been completely solved within 4–120 seconds; problems with a number of feasible routes from 220,000 to 600,000 have been completely solved within 140–343 seconds; finally problems with a number of feasible routes from more than a million to almost seven million need between 637 and 3464 seconds.

5. Conclusions In this paper, we propose to apply a two phase approach to solve optimally real medium-sized VRPTW with hard and soft time windows and multiple objectives. It combines enumeration of feasible routes and goal programming. In the first phase feasible routes are enumerated and the total penalty incurred by each route due to deviations from targets is computed. With these data, in the second phase a set partitioning problem is solved in order to get the best set of feasible routes. The procedure provides an optimal solution to the problem analyzed. It permits computation of complicated nonlinear expressions for operational costs and the addition of a lot of constraints without affecting the complexity of the problem solved in the second phase. In fact, in the first phase, the more the number of constraints, the less the number of feasible routes, making it much easier to solve the set partitioning problem. Moreover, at the end of the first phase a filtering step can be applied which reduces the number of variables of the set partitioning problem. This process removes some routes, that will not be in the set of optimal routes, from the set of feasible routes.

1732

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

We do not recommend the use of this methodology for every problem. Its main drawback is the impossibility of solving problems with a large number of feasible routes. But rather, our aim is to propose a methodology that provides an optimal solution and works well when there are multiple objectives and constraints that prevent the use of mathematical programming techniques directly. Computational results show that, at this moment and taking into account currently used computer technology, this methodology is competitive for solving real-life medium-sized problems. It provides optimal solutions in a computational time that managers would consider more than acceptable. Acknowledgment Many thanks are due to anonymous referees for valuable suggestions to improve the presentation of the paper. References [1] N. Balakrishnan, Simple heuristics for the vehicle routing problem with soft time windows, Journal of the Operational Research Society 44 (3) (1993) 279–287. [2] M.O. Ball, T.L. Magnanti, C.L. Monma, G.L. Nemhauser (Eds.), Network Routing. Handbook in Operations Research and Management Science, vol. 8, North-Holland, Amsterdam, 1995. [3] J.F. Bard, G. Kontoravdis, G. Yu, A branch-and-cut procedure for the vehicle routing problem with time windows, Transportation Science 36 (2) (2002) 250–269. [4] H.I. Calvete, C. Gale´, M.J. Oliveros, B. Sa´nchez-Valverde, Vehicle routing problems with soft time windows: An optimization based approach, in: M.C. Lo´pez de Silanes et al. (Eds.), Proceedings of VIII Journe´es Zaragoza-Pau de Mathe´matiques Applique´es et Statistiques, vol. 31, Monografı´as del Seminario Matema´tico Garcı´a Galdeano, Univ. de Zaragoza, 2004, pp. 295– 304. [5] V. Chankong, Y.Y. Haimes, Multiobjective Decision Making, Theory and Methodology, North-Holland, New York, 1983. [6] J.F. Cordeau, M. Gendreau, G. Laporte, J.Y. Potvin, F. Semet, A guide to vehicle routing heuristics, Journal of the Operational Research Society 53 (2002) 512–522. [7] J.F. Cordeau, G. Laporte, A. Mercier, A unified tabu search heuristic for vehicle routing problems with time windows, Journal of the Operational Research Society 52 (2001) 928–936. [8] M. Desrochers, J. Desrosiers, M. Solomon, A new optimization algorithm for the vehicle routing problem with time windows, Operations Research 40 (1992) 342–354. [9] K. Fagerholt, Ship scheduling with soft time windows: An optimisation based approach, European Journal of Operational Research 131 (2001) 559–571. [10] M. Gendreau, G. Laport, J. Potvin, Vehicle routing: Modern heuristics, in: E. Aarts, J.K. Lenstra (Eds.), Local Search in Combinatorial Optimization, Wiley, London, 1997, pp. 311–336. [11] B.L. Golden, A.A. Assad (Eds.), Vehicle routing: Methods and Studies, Elsevier Science Publishers, Amsterdam, 1988. [12] G. Ioannou, M.N. Kritikos, A synthesis of assignment and heuristic solutions for vehicle routing with time windows, Journal of the Operational Research Society 55 (2004) 2–11. [13] G. Ioannou, M.N. Kritikos, G. Prastacos, A greedy look-ahead heuristic for the vehicle routing problem with time windows, Journal of the Operational Research Society 52 (2001) 523–537. [14] N. Kohl, O.B.G. Madsen, An optimization algorithm for the vehicle routing problem with time windows based upon lagrangian relaxation, Operations Research 45 (1997) 395–406. [15] A.W.J. Kolen, A.H.G. Rinnooy Kan, H.W.J.M. Trienekens, Vehicle routing with time windows, Operations Research 35 (2) (1987) 266–273. [16] Y.A. Koskosidis, W.B. Powell, M.M. Solomon, An optimization-based heuristic for vehicle routing and scheduling with soft time window constraints, Transportation Science 26 (2) (1992) 69–85. [17] H. Min, A multiobjective vehicle routing problem with soft time windows: The case of a public library distribution system, SocioEconomic Planning Science 25 (3) (1991) 179–188. [18] G. Nemhauser, L. Wolsey, Integer and Combinatorial Optimization, John Wiley and Sons, New York, 1999. [19] C. Romero, Handbook of Critical Issues in Goal Programming, Pergamon Press, Oxford, 1991. [20] R. Ruiz, C. Maroto, J. Alcaraz, A decision support system for a real vehicle routing problem, European Journal of Operational Research 153 (2004) 593–606.

H.I. Calvete et al. / European Journal of Operational Research 177 (2007) 1720–1733

1733

[21] M. Schniederjans, Goal Programming Methodology and Applications, Kluwer Publishers, Boston, 1995. [22] T.R. Sexton, Y.M. Choi, Pickup and delivery of partial loads with soft time windows, American Journal of Mathematical and Management Sciences 6 (3,4) (1986) 369–398. [23] M.M. Solomon, Algorithms for the vehicle routing and scheduling problem with time windows constraints, Operations Research 35 (1987) 254–265. [24] M.M. Solomon, J. Desrosiers, Time window constrained routing and scheduling problems, Transportation Science 22 (1) (1988) 1–13. [25] R.E. Steuer, Multiple Criteria Optimization: Theory, Computation, and Application, John Wiley and Sons, New York, 1986. [26] E. Taillard, P. Badezu, M. Gendreau, F. Guertin, J.Y. Potvin, A tabu search heuristic for the vehicle routing problem with soft time windows, Transportation Science 31 (2) (1997) 170–186. [27] M. Tamiz, D.F. Jones, C. Romero, Goal programming for decision making: An overview of the current state-of-the-art, European Journal of Operational Research 111 (1998) 569–581. [28] P. Toth, D. Vigo (Eds.), The Vehicle Routing Problem, in: SIAM Monographs on Discrete Mathematics and Applications, vol. 9, l, 2002. [29] M. Zeleny, Multiple Criteria Decision Making, McGraw-Hill, New York, 1982.