3. Tabu-Embedded Simulated Annealing - Semantic Scholar

Report 3 Downloads 118 Views
Manpower Scheduling with Time Windows Andrew Lim1, Brian Rodrigues2 and Lei Song3 1

Department IEEM, Hong Kong University of Science and Technology, Clear Water Bay, Kowloon, Hong Kong 2 School of Business, Singapore Management University, 469 Bukit Timah Road, Singapore 259756 3 Department of Computer Science, National University of Singapore, 3 Science Drive 2, Singapore 117543

Abstract In this paper, we propose a manpower allocation model with time windows which is of practical interest to serviceman scheduling operations. Specifically, this problem originates from peculiar port yard scheduling needs where demand is generated from locations in the yard for servicemen who are dispatched from a central point and where the objectives are to minimize the number of servicemen scheduled, travel distances, travel times and waiting times at each location. Although closely related to the well-known vehicle routing problem, this problem is different while its solution could provide insight to the latter. We develop solutions using metaheuristic methods, and in particular provide tabu-embedded simulated annealing and squeaky wheel optimization with local search algorithms for this problem. We apply these newly-developed metaheuristics with adaptations for solutions to the manpower allocation problem, while our analysis throws light on how these work. Computational results are reported which show the effectiveness of our approach when applied to the manpower allocation problem.

1. Introduction This work is motivated from a study of manpower allocation needs of the Port of Singapore Authority (PSA) which is a large technology corporation located in Singapore. It manages one of the busiest ports in the world handling over seventeen million TEU’s (twenty-foot equivalent units containers) annually representing about 9% of global container business. PSA is concerned with maximizing throughput at its ports in view of pressures derived from limited port size, high cargo transshipment volumes and limited physical facilities and equipment. A simplified description of the manpower allocation problem for the port is as follows: Different locations in the port yard demand varying amounts of service work to be met. A service control center dispatches servicemen to satisfy these demands. Each serviceman performs one task at one location at any one time. Each location has its own demand time window, which is defined by an “early” time and a “late” time. A serviceman has to start a task only after the “early” time and finish before the “late” time, and has to wait if he arrives before the early time. We make the assumption here that servicemen have the same travel speed and work rate and can move from location to location but will eventually be required to return to the service control center. Each demand location has a requirement for a given number of servicemen, which is in proportion to its work needs, to fulfill work within its demand time window. _____________________________________________ Permission to make digital or hard copies of all or part of this work for personal or classroom use is granted without fee provided that copies are not made or distributed for profit or commercial advantage and that copies bear this notice and the full citation on the first page. To copy otherwise, or republish, to post on servers or to redistribute to lists, requires prior specific permission and/or a fee SAC 2003, Melbourne, Florida, USA  2003 ACM 1-58113-624-2/03/03 . . . $5.00

The main purpose of this study is to develop and solve a model for this basic manpower allocation problem which, as described, takes into account practically important ground constraints. Specifically we formulate this problem as a multi-criteria problem where the primary objective is to minimize the number of servicemen used in satisfying all service demands, while the secondary objective is to minimize the total travel distance for such a schedule. Other tertiary objectives are to minimize total scheduled time and total waiting time. Our model, with the given objectives and criterion as stated here, is closely related to the well-known Vehicle Routing Problem with Time Windows (VRPTW) (see Pardalos and Resende 2002, for example). It is however, essentially different from the latter in that the demand (or customer points) require a fixed number of servicemen whereas in the VRPTW there is no such requirement. Furthermore, a location may have more than one demand which overlap in time. Because of the closeness of our manpower allocation model with the VRPTW, we cite references to work done on the latter as the most relevant to our manpower allocation model. Savelsbergh (1985) showed the VRPTW to be NP-hard. Much work has been done to derive exact algorithms to solve the problem. In general, more recent research on the VRPTW can be divided into optimization and heuristic approaches. Results obtained recently by heuristics are proposed in (Rochat, Taillard, 1995), (Chiang, Russell, 1996), (Taillard, et al., 1997), (Thangiah 1999) (Homberger and Gehring, 1999), (Cordeau et al. 2000), (Backer et al. 2000), (Rousseau et al. 2000) and (Li and Lim 2002).

1.2 Computational Complexity and Objectives The manpower allocation model is clearly NP-hard. It suffices to consider the minimum travel distance objective with the demands at locations set to unity without time windows to see that this problem has larger complexity than that for the traveling salesman problem. In general, methods applied to similar problems, such as the VRPTW include heuristics and metaheuristics. The former include construction heuristics, improvement heuristics and composite heuristics. In particular, route construction heuristics build feasible solutions by inserting an unrouted location into a current partial route until all locations are routed. The sequential insertion algorithm proposed by Solomon (Solomon, 1987) belongs to the class of construction heuristics. Route improvement heuristics modify a solution by performing local searches for better neighboring solutions within its neighborhoods generated by node/edge-swapping operators. A composite heuristic usually mixes route construction and route improvement procedures. Metaheuristics are a new generation of heuristics including tabu search, simulated annealing, evolutionary algorithms and hybrids. The objective of this study is to find an efficient way to solve the manpower allocation problem with time windows with the use of heuristics and metaheuristics. We will use a relatively new technique, “Squeaky Wheel Optimization” (SWO), for this problem. Because this is a new problem, we generate problem

instances independently. To compare results from SWO, we implement another algorithm which is Tabu-Embedded Simulated Annealing. Our paper is organized as follows. In Section 2 we give a problem statement and model formulation. In Section 3, we introduce tabu-embedded simulated annealing, an algorithm that differs from traditional implementation, and in Section 4 we implement Squeaky Wheel Optimization with Local Search to the manpower allocation model. In Section 5, we provide experimental results and in Section 6, we have our conclusion.

Decision variables: 0 if there is no arc form i to j and 1 otherwise, i ≠ j, i, j = 0,1,…,N; k =1,…, M arrival time at node i wait time at node i total number of servicemen used

Parameters: N ci c0 cij tij di ei li si mk

total number of demand nodes customers customer i, i = 1,…,N service center distance cost between node i and j travel time cost between nodes i and j demand in numbers of servicemen at node i earliest arrival time at node i latest arrival time at node i service time at node i maximum shift time allowable for serviceman k

minimize

N

(2.1) N

M

∑∑∑ c i = 0 j = 0 k =1

ij

xijk

M

minimize

N

N

ij

(2.3)

xijk

M

∑∑∑ w x i = 0 j = 0 k =1

i

(2.4)

ijk

Subject to: M

N

∑∑ x k =1 j = 0

M

ijk

N

∑∑ x

ijk

= di

for i = 0,1,…,N

(2.5)

(t i + t ij + si + wi ) ≤ t j

for j = 1,…,N

(2.6)

for i = 0,1,…,N

(2.8)

for k =1,…, M

(2.9)

for k = 1,…, M

(2.10)

(2.7)

ei ≤ (t i + wi ) ≤ li N

N

∑∑ x i =0 j =0

N

∑x i =1

i 0k

ijk

(t ij + s i + wi ) ≤ mk

= 1,

N

∑x j =1

0 jk

=1

The objectives given in (2.1) to (2.4) are in deceasing priority. Constraint (2.5) represents the demand from node i for di servicemen. Constraints (2.6) to (2.9) represent time- window constraints and (2.10) ensure that each serviceman begins from and returns to the service center. We note that, in comparison to the VRPTW, constraint (2.5) replaces the usual vehicle capacity constraints, or variants, (see, for example, Pardalos and Resende 2002) found in the latter. Here, servicemen are capacitated only by their service rate which are subsumed in the demand requirements. The following notations will be used throughout the paper. G (V, A) G is a digraph generated from problem data. V = VN ∪ {c0 } is the node set of G, where VN = (c1,…,cN) represents the demand locations and node c0 denotes the service center. A is the set of arcs weighted by travel and distance costs Cost (s) objective function of solution s Dist (s) total travel distance of solution s SchedTime(s) total schedule time of solution s Wait (s) total waiting time of solutions s s* a local optimal solution.

3. Tabu-Embedded Simulated Annealing

The multi-criteria model for the MAPTW with the above definitions is then: Objectives: minimize M

N

t o = w0 = s o = 0

This section describes the Manpower Allocation Problem with Time Windows (MAPTW). This problem is represented by a set of servicemen and a set of demand nodes. It is assumed that there are a total of N + 1 locations, including the service center and that there are enough servicemen to meet demand. Each arc in the network represents a connection between nodes and indicates the direction of travel and each route starts at the service center visits a customer and can proceed from there to visit other customers before returning to the service center. A travel distance cost cij and a travel time cost tij are associated with each arc. Each customer in the route has a demand di of numbers of servicemen required who must provide service within the time window specified for the customer location given by an earliest and a latest time. Servicemen who arrive early wait. The decision variables are xijk for i,j = 0,1,…,N and k = 1,…, M, which is 1 if serviceman k travels from customer i to customer j, and 0 otherwise. The decision variable ti denotes the time of arrival of the serviceman at location i, and wi is the wating time at i. The model is as follows:

ti wi M

N

∑∑∑ t i = 0 j = 0 k =1

k =1 i = 0

2. Problem statement and Formulation

xijk

minimize

(2.2)

This section describes preliminary algorithms devised to solve the MAPTW which adapts simulated annealing with tabuembedding. Since it is difficult to choose suitable configuration parameters for simulated annealing, we will introduce an algorithm which differs from traditional implementation. Instead of performing the simulated annealing procedure on the probabilistically accepted solution repeatedly until termination, we force the simulated annealing procedure to restart from the current best solution after several simulated annealing iterations. The algorithm terminates after K restarts without improvement (Li and Lim 2001, Li and Lim 2002). In addition to the K-restart strategy, the algorithm is strengthened by recording visited solutions into a tabu set so as to avoid cycling. In effect, our algorithms attempts to make better use of the intensification strategy on the best solution, while

diversifying exploration using a tabu-embedded simulated annealing strategy. The algorithm is divided into two steps: Step1: Obtaining an initial solution for the problem using a heuristic. Step2: Starting a local search from the initial solution, using tabu-embedded simulated annealing with a K-restart strategy. In the Step 1, we use the insertion heuristic proposed by Solomon to produce initial solutions. The initial solutions are then used to be the input of the local search procedure in the second step.

3.1 Neighborhood Generating Mechanisms In our method, three edge-swapping operators are used to generate local search neighborhoods. These operators are: the shift operator, the exchange operator and the rearrange operator (Li and Lim et al. 2002). Savelsbergh and Tailard et al. also used similar edge-exchange operators. Among these operators, the shift operator and the exchange operator are used for generating neighborhoods for local searches in the tabu-embedded simulated annealing procedure. On the other hand, the rearrange operator is used only for post-processing on local optima obtained by local searches.

solutions in a tabu list. (This differs from some researchers who use edge-moves in their tabu structure.)

3.2 The Algorithms The pseudo code of the tabu-embedded simulated annealing procedure is described as follows: INPUT: a solution

OUTPUT: a local optimal solution

2. 3. 4.

s* = s . Set the variable NoImprovment = 0 . While NoImprovment < L DO

Set the initial best solution

a.

Use TABU_EMBEDDED_METROPOLIS_PROC (s) to obtain a non-tabu feasible neighboring solution

s' . b.

Using the shift operator, non-null route segments, which contain at least one node, can be moved from one route to another route under the condition of not violating time window constraints. For each pair of selected routes, say, Route i and Route j, the shift operator is used in two directions, so that two sub-neighborhoods can be obtained by shifting route segments between the two routes.

Perform a DESCENT_LOCAL_SEARCH within neighborhoods obtained by basic operations and set

c.

IF

ii. Set ELSE set NoImprovment

3.1.4 Objectives and Cost Function in Local Search

e. 5.

Return

Set

s

s =s *

'

= NoImprovment + 1

.

*

The procedure TABU_EMBEDDED_METROPOLIS_PROC is as follows: ALGORITHM TABU_EMBEDDED_METROPOLIS_PROC (Solution s) 1. WHILE NOT accept a neighboring solution s DO

The priority-order of the objectives for solving the MAPTW as given (2.1) to (2.4) above. To reflect the order of priority, we use the following cost function:

a.

Set ∆ = Cost ( s ' ) − Cost ( s)

b.

IF

∆≤0

prob = 1

i. THEN set variable

Cost ( s) = α .M + β .Dist ( s) + γ .SchedTime( s) + λ.WaitTime( s)

Here, the penalty weight factors are set to be: α > β > γ > λ . Using this setting, even if initial solutions obtained by the route construction heuristic contain a large number of servicemen, this number will be reduced by local searches. In addition, we can change the setting of these factors to get different priority orders of the four objectives. To introduce tabu search in local searches, we record the four objective values in a data structure with four fields. Since the probability that two different solutions have the same record values is very small, it is reasonable to regard two solutions as the same if their record values are not different. We store the records of the

'

THEN

s =s NoImproviment = 0 *

i. Set

3.1.2 The Exchange Operator

The rearrange operator can be used to reposition route segments within the same route while maintaining the time-window constraints.

s' Cost ( s ' ) < Cost ( s )

the result to

d.

3.1.3 The Rearrange Operator

s*

ALGORITHM TABU_EMBEDDED_SIMULATED_ANNEALING(Solution s) 1. Generate the initial solution s.

3.1.1 The Shift Operator

Using the exchange operator, non-null route segments can be exchanged between two routes under the condition of not violating time-window constraints.

s

prob = e random(0,1) ≤ prob ' i. Accept s ii. ELSE

c.

2.

IF

RETURN

s

−∆ / T

THEN

ii. Update

annealing

iii.

into tabu list

T = δT ' Record s

temperature:

'

The procedure DESCENT_LOCAL_SEARCH is as follows:

ALGORITHM DESCENT_LOCAL_SEARCH (Solution s)

s* = s

1.

Set

2.

Select the best solution s with minimum objective cost within a defined neighborhood of s

3.

IF Cost( s ' ) < Cost ( s ) THEN

'

a. b. 4.

s* = s ' GOTO Step 2

ELSE Return

s*

In these algorithms, the cost function

Cost (s )

is defined as:

Cost ( s ) = Dist ( s ) + ϕ × WaitTime( s ), whereϕ = 0.01 × Dist ( s ) . We note that in the procedure TABU_EMBEDDED_SIMULATED_ANNEALING, from the viewpoint of tabu search, the TABU_EMBEDDED_METROPOLIS_PROC in Step 4 (a) above provides an annealing-probabilistic diversification strategy for escape from local optima. In addition, Step 4 (b) is an intensification strategy.

3.2.1 General Framework of the Algorithm The general framework of the algorithm is a two-phase heuristic. In the first phase, it reads data from a problem data file and then obtains an initial solution, s, using Solomon’s Insertion Heuristic. Then, in second phase, the K-restart search procedure starts from this initial solution s, which is by default the current best solution. If a better solution is found at the kth. iteration, where k ≤ K , then the current best solution is overwritten by the new better solution and the search is restarted and the iteration counter is reset to zero. Otherwise, the program terminates after k > K . In fact, the best-based K-restart strategy acts as an intensification strategy in favor of the latest best solution.

4. Squeaky Wheel Optimization Squeaky Wheel Optimization (SWO) is a general approach to optimization (Joslin and Clements, 1998.). The core of SWO is a Construct-Analyze-Prioritize cycle. A solution is constructed by a greedy algorithm, making decisions in an order determined by priorities assigned to the elements of the problem. That solution is then analyzed to find the elements of the problem that are “trouble makers”. The priorities of the “trouble makers” are then increased, causing the greedy constructor to deal with them sooner on the next iteration.

4.1 Squeaky Wheel Optimization with Local Search The SWO algorithm has been proven to be effective to solve some job scheduling and graph coloring problems and outperform Tabu Search and Integer Programming (Joslin and Clements, 1998.). However, it has some drawbacks: (1) Although SWO tends to avoid getting trapped in local optima by changing the sequence in prioritizer, the same blame values will still make SWO trapped in small cycles. The two existing ways to solve it is to introduce some randomness and restart periodically with a new initial sequence; (2) Problems sometimes have tasks that must be handled badly in order to achieve a good overall solution. However, the analyzer always assigns high blames to those “sacrificial” tasks to

make the constructor handles them well in the next iteration. Thus the good overall quality may not be achieved. To reduce these drawbacks of SWO, we propose an enhanced framework of SWO. Instead of having one unique initial solution, we use a generator to randomly generate a problem sequence after certain iterations. This is used to reduce drawback (1). Another enhancement is that after we get a greedy solution from the constructor, we perform a local search on the solution and try to find a better solution among its neighbors. This improvement can reduce the drawback (2) but will increase the calculation time.

4.2 Implementation of SWOL We describe here the implementation of the five main components of SWOL: Generator The generator randomly picks one problem element each time and uses a greedy algorithm to select the best insertion position into the solution. A task is added by selecting a position relative to the tasks already in the routes. A task may be inserted between any two tasks already in the routes or at the beginning or end of any route. For each of the possible insertion positions in the schedule, relative to the tasks already in each route, the constructor calculates the effect on the objective function, and the selected task is placed at the best-scoring location. Ties are broken randomly. Constructor The constructor generates a solution by inserting a task one at a time, in the order they occur in the priority sequence using the same algorithm as the Generator does. After all tasks have been assigned, the constructor applies local search to the solution, attempting to improve the score by reordering the tasks that were assigned to it. The following is the pseudo code of Constructor: Local Search The Local Search in SWOL just uses the LOCAL_DESCENT_SEARCH algorithm used in Tabu-Embedded Simulated Annealing which will try to search for better solutions in the neighborhoods generated by three basic operations. During the implementation, because the full local search takes too much time, it is not feasible to do a full local search for each construction by Constructor. What we do is to perform partial local search, which searches the neighborhoods only for a certain time, for each construction and performs a full local search after certain iterations. Analyzer The analyzer will assign blame to each demand task. How large the blame value depends on how the current assignment of the demand affects the solution. The analyzer first assigns highest blame value to those demands in the shortest routes. This approach is aimed at finding “fault” with servicemen who may be underutilized. If very few locations are served by a serviceman, it may be possible to schedule them in other routes, hence reducing the number of servicemen needed. For each demand in the same route, the analyzer calculates a lower bound of the traveling distance that each demand could contribute to the solution. The blame assigned to the demand is its “excess cost,” the difference between its actual costs and its minimum possible cost.

Prioritizer Once the blame has been assigned, the prioritizer modifies the previous sequence of tasks by moving demands with non-zero blame factors forward in the sequence. Demands are moved forward a distance that increases with the magnitude of the blame. The prioritizer will reorder the customers in the priority sequence by applying a multi-level stable sort, ordering them first by vehicle number blame, then traveling distance blame. General Framework If the Constructor-Local Search-Analyzer-Prioritizer cycle runs a number of iterations without any improvement, a full local search is performed to the neighborhoods of the best solution obtained and for next iteration the generator generates another initial solution. This whole process runs times a number of times.

5. Experimental Results & Comparisons Problem instances are generated here for experiments with the algorithms developed for the MAPTW. There are some important factors about the instances that affect the behavior of routing and scheduling algorithms: geographical data; the number of tasks serviced by a serviceman; the percentage of time-constrained customers and tightness and positioning of the time windows. To reflect these factors in the problem instances, we generate different sets of testing data, which are similar as Solomon’s VRPTW problem instances (Solomon, 1987). Two important factors are considered. They are geographical data and the time needed to satisfy one demand. The geographical data can be randomly distributed or clustered or mixed and the time required to finish one task can be long or short. There are six sets of instances in total: 1.MP_C1 Type: the geographical data are clustered in these problem sets and the service time for each task is long (from 20 to 50) so that they have short scheduling horizon and allow only a few demand fulfilments per serviceman (about 7 to 12). 2.MP_C2 Type: the geographical data are clustered and the service time for each task is short (from 2 to 10) so they have a long scheduling horizon permitting many tasks (about 30) to be done by the same serviceman. 3. MP_R1 Type: the geographical data are randomly generated and they have short scheduling horizon. 4. MP_R2 Type: the geographical data are randomly generated and they have long scheduling horizon. 5. MP_RC1 Type: the geographical data are mixed. That is, some locations can randomly generated and some are clustered. The service time is long for RC1 type. 6. MP_RC2 Type: the geographical data are mixed and the service time is short. Each type has four files. The location coordinates are identical for these four files. The problems differ with respect to the width of the time windows. Some have very tight time windows, while others have time windows which are hardly constraining. In terms of time window density, that is, the percentage of customers with time windows, problems have 25, 50, 75 and 100 % time-window density.

5.1 Experiments The hardware used to run our experiments is the PC Pentium III 900MHz. The algorithm was developed using Java programming language. The OS used is Win2000 Professional. Because the problem instances are generated by random, we can’t know the optimal solutions. What can we do is to compare the results among different algorithms. Here three algorithms are compared, which are pure Greedy algorithm (GA), TabuEmbedded Simulated Annealing (TSA) and Squeaky Wheel Optimization with Local Search (SWOL). We find that pure greedy algorithm performed the worst while TSA can always improve the travel distance greatly. SWO can generate good enough results in most of the cases and for most instances, while SWOL can produce a solution which has fewer servicemen when compare to TSA.

5.2 5.2.1

Analysis Tabu-Embedded Simulated Annealing

There are two important aspects for each algorithm: diversification and intensification. TSA is in fact a depth-like search which starts from a neighbor in

N (s) , and stretches the search away from N (s) by a TABUEMBEDDED METROPOLIS PROCEDURE (TMP) and DESCENT LOCAL SEARCH (DLS). Actually, only with this depth-like method can solutions outside of N (s) be explored in each restart of TSA on s. Thus, the local searches around s are diversified. On the other hand, the K-restart strategy guides TSA to start local searches within N (s) in breadth-like way. It intensifies local searches in favor of N (s) . However, the diversification of TSA is not enough, because it only depends on TMP and DLS to explore other solution spaces. If the optimal solution is surrounded by some infeasible solutions in the solution space, TSA will fail to find this solution. 5.2.2 Squeaky Wheel Optimization with Local Search In an analysis, it is useful to think of SWO as searching two coupled spaces. One search space is the familiar solution space, and the other is priority space. Moves in the solutions spaces are made indirectly, via the re-prioritization that results from analyzing the prior solution. Similarly, successive prioritizations are generated by constructing and analyzing a solution, and then using the blame that results from that analysis to modify the previous prioritization. One consequence of the coupled search spaces is that a small change in the sequence of elements generated by the prioritizer may correspond to a large change in thecorresponding solution generated by the constructor, compared to the solution from the previous iteration. Moving an element forward in the sequence may change its state in the resulting solution. In addition, any elements that now occur after it in the sequence must accommodate that element’s state. The result is a large move that is “coherent” in the sense that it is similar to what we might expect from moving the higher priority task and then propagating the effects of that change by moving lower priority tasks as needed. This single move may correspond to a large number of moves for a search algorithm that

only looks at local changes to the solution. Throughout this method, SWO achieves diversification. However, intensification of SWO is not good since for each solution in the solution space, SWO will move to another solution without searching the neighbors. Local Search is then used to append SWO to gain intensification. In the experiments, we find that SWOL can always find out solutions which contain the minimum number of servicemen. However, the travel distance is not optimized. This is because Local Search only searches the neighbors generated by three operators which may not be broad enough.

5. Conclusion In this paper, we proposed a new problem that has practice value to the industry. We then proposed, implemented a number of good methods to solve the problem. The experiments showed that TSA and SWOL applied to the MAPTW perform well when compared with a greedy approach. Because this is a new problem there have been no benchmarks for comparisons. However, results found here were better that those already available for PSA, and close to optimal for small data sets for which optima can be easily found and verified. We have also applied newly-developed metaheuristics to this problem with adaptations and analyzed these applications. Our experience with applications of metaheuristics has been that it is fairly straightforward to implement SWO because there are usually fairly obvious ways to construct greedy solutions, and to analyze a solution to assign blame to some of the elements. Naïve implementations of SWO tend to perform reasonably well. We have found the view of SWO as performing a “coupled search” over two different search spaces to be helpful in characterizing the kinds of moves that SWO makes in each of the search spaces, and the effect this has on avoiding local optima. This we hope gives us a deeper understanding of what makes SWO work and the ability to more effectively design SWO algorithms. Although the ability to make large, coherent moves is a strength of the approach, it is also a weakness. SWO is poor at making small “tuning” moves in the solution space, so that we need to combine a local search with SWO.

References Backer, B., Furnon, V. Shaw, P., Kilby, P. Prosser, P., (2000), Solving vehicle routing problem using constraint programming and metaheuristics, Journal of Heuristics 6, pp. 501-523. Chiang, W., Russell, R., (1997), A reactive tabu search metaheuristics for the vehicle routing problem with time window, INFORMS Journal on Computing 9, pp. 417-430. Cordeau, J., Laporte, G., Mercier, A., (2000), A unified tabu search heuristic for vehicle routing problem with time windows, Technical Report. CRT-00-03, Center for Research on Transportation, Montreal, Canada. Eilon, S., Watson-Gandy C. and Christofides, N., (1971), Distribution Management: Mathematical Modeling and practical Analysis. Hafner, New York. Homberger, J., Gehring, H., (1999), Two evolutionary metaheuristics for the vehicle routing problem with time window, INFOR 37 pp. 165-172.

Joslin, D. and Clements, David P., (1998) “Squeaky Wheel” Optimization, AAAI/IAAI 1998, pp. 340-346. Li, H.B. and Lim, A. (2002) Local Search with Annealing-like Restarts to solve the VRPTW, European Journal of Operational Research, accepted for publication. Li, H.B. and Lim, A. (2001), A Metaheuristic for solving the pickup and delivery problem with time windows, The 13th IEEE International Conference on Tools with Artificial Intelligence (ICTAI-2001, Dallas, USA), 151-158. Pardalos, P. and Resende, M., Editors, (2002), Handbook of Applied Optimization, Oxford University Press. Rochat, Y., Taillard, E., (1995), Probabilistic diversification and intensification in local search for vehicle routing, Journal of Heuristics 1 pp. 147-167. Rousseau, J., Gendreau, M., Pesant, G. (2002 ) Using constraintbased operators to solve the vehicle routing problem with time windows, Journal of Heuristics 8(1), pp. 43-58. Savelsbergh, M., (1985), Local search for routing problem with time windows, Annals of Operations Research 4, pp. 285-305. Solomon, M., (1987) Algorithms for the vehicle routing and scheduling problems with time window constrains, Operations Research 35, pp. 254-264. Thangiah, S.R., (1999), A hybrid genetic algorithm, simulated annealing and tabu search heuristics for the vehicle routing problem with time windows, In Chambers, Editor, Practical Handbook of Genetic Algorithms: Complex Coding Systems, Vol. III, CRC Press, Boca Raton, pp. 253-277. Taillard, E., Badeau, P. Gendreau, M. Geurtin, F., Potvin, J., (1997) A tabu search heuristic for the vehicle routing problem with time windows, Transportation Science 31 170-186.