Journal of the Operational Research Society (1996) 47, 678-689

0

1996 Operational Research Society Ltd. All rights reserved. 0160· 5682/96 $12.00

Search Heuristics for Resource Constrained Project Scheduling JAE-KWAN LEE 1 and YEONG-DAE KIM 2 1 LG-EDS

Corporation, Korea and 2 Korea Advanced Institute of Science and Technology

We develop a search procedure for project scheduling problems with multiple resource constraints as well as precedence constraints. The procedure is applied to three popular search heuristics, simulated annealing, tabu search and genetic algorithms. In the heuristics, a solution is represented with a string of numbers each of which denotes priority of each activity. The priorities are used to select an activity for scheduling among competing ones. The search heuristics with this encoding method can always generate feasible neighbourhood solutions for a given solution. Moreover, this encoding method is very flexible in that problems with objective functions of a general functional form (such as a nonlinear function) and complex constraints can be considered without much difficulty. Results of computational tests on the performance of the search heuristics showed that the search heuristics, especially the simulated annealing and tabu search algorithms worked better than existing heuristics. Key words : heuristics, project management, search

INTRODUCTION The critical path method (CPM) and the project evaluation and review technique (PERT) have been widely used in government and industry project management problems. However, CPM and PERT are based on the assumption that there is an infinite amount of resource available for a project at any time. Since this restricts their use in many practical problems, a number of attempts have been made to solve the project scheduling problem with resource constraints as well as precedence constraints, i.e., the resource constrained project scheduling problem (RCPSP). There are many different versions of this problem in the literature and they can be characterized by the number of projects that are to be simultaneously scheduled (single or multiple), the objective function, characteristics of resources, and the preemption condition, among others. (A recent survey paper by Elmaghraby 1 discusses various aspects of activity networks and project scheduling.) Minimization of project duration is often used as an objective of the problem 2 - 12 , while other objectives such as minimization of the total project cost 11 , maximization of the net present value of cash flows 1 0 • 13 - 1 S, and levelling resource usage 16 are used in other research. In addition, in cases in which there are multiple projects to be scheduled, we may use minimization of the tardiness penalties 17 · 18 or other functions of individual project completion times as an objective. There are two types of resources in RCPSPs: renewable and nonrenewable resources. A renewable resource is one that is needed for processing of an activity but can be recovered when the activity is completed. On the other hand, if the total amount of a resource is limited over the project duration and the consumed resource cannot be recovered, the resource is called nonrenewable. Labour hours and money are good examples of the former and the latter, respectively. If both their per-period and total availability are limited, the resources are defined as doubly constrained. In most research, resources are assumed to be of the renewable type.11 Project scheduling is an important issue in industry, as indicated by the growth of the market for project management software, which has continued to grow at a rate of 18.6% per year 19 . However, most existing software packages (PC-based, at least) do not include sufficient details of the capability of dealing with resource constrained project scheduling problems. De Wit and Herroelen 20 report that the procedures for resource levelling and resource constrained scheduling used by many software packages are almost never explained in detail. This motivated our study to develop an efficient solution approach for RCPSPs. In this study, we consider the standard Correspondence: Y-D Kim, Department of Industrial Engineering, Korea Advanced Institute of Science and Technology, Yusong-gu, Daejon 305-701, Korea

Operational Research Society is collaborating with JSTOR to digitize, preserve, and extend access to Journal of the Operational Research Society. ® www.jstor.org

J-K. Lee and Y-D. Kim-Search Heuristics for Project Scheduling

679

version of the problem where resources are renewable, preemption is not allowed, and the objective is to minimize the project duration (makespan). There have been a number of attempts to find optimal solutions for RCPSPs, such as recently developed algorithms of Bell and Park 3 , Christofides et al. 5 , and Demeulemeester and Herroelen 7 . Since the problem is NP-complete 21 , most optimal solution algorithms employ branch and bound or enumeration approaches. There is still no efficient procedure that is computationally feasible for large and complex projects which occur in practice. Generally, it is known that optimal solution methods are appropriate only for small or moderate size projects. This lack of success in finding an efficient optimal solution algorithm motivated development of heuristic solution methods. Most existing heuristic methods use priority rules in making scheduling decisions. Priority rules resolve resource conflicts by scheduling an activity with a higher priority value earlier than other competing activities. That is, when resources become available, the methods select and schedule an activity with the highest priority value among schedulable activities, activities whose predecessors are all completed and which do not require resources more than the amounts available at the time. Generally, a priority of an activity in most priority rules is a function of time (duration), amounts of resources required and location of the activity on the network. A functional form of these attributes affects performance of a priority rule. Performance of a rule may also depend on the network topology, the numbers of activities and constraints, and the tightness of resource constraints. Although there are many priority rules, Davis and Patterson 6 argue that there is little basis, a priori, for making a choice among different heuristic procedures and that no priority rule dominates all others or performs consistently better than other rules. To cope with this feature of individual rules, we can include several priority rules in an algorithm and select the best solution among those from the rules as a final solution for the problem. Boctor4 and Khattab and Choobineh 16 suggest such composite heuristics and show that these give better solutions than heuristics that use individual rules. As another approach to RCPSPs, improvement procedures can be applied after a feasible schedule is obtained. Bell and Han 2 , Li and Willis 9 and Yang et alY develop multi-pass heuristic methods that employ improvement procedures, and show that single-pass heuristics are outperformed by these multi-pass heuristic methods. Also, local search methods can be applied after an initial schedule is obtained. Sampson and Weiss 22 suggest a local search method in which neighbourhood solutions are obtained from the current feasible solution by increasing or decreasing start and finish times of an activity and a move is made to obtain an improved solution. On the other hand, Oguz and Bala 23 use a more direct local search method, called MIXED heuristic. This method, being a general purpose algorithm designed for integer programming problems, uses subgradient and n-dimensional search techniques. However, these local search methods terminate when they find a local minimum without escaping from it to pursue the global minimum. To overcome the weakness of such greedy type search methods,, i.e., being trapped in a local minimum, three search heuristics are often used. These are simulated annealing, tabu search, and genetic algorithms, which have specific mechanisms of their own for escaping from a local minimum. The basic ideas behind the mechanisms of these algorithms are fairly general so that they can be used in various combinatorial optimization problems. Recently, Icmeli and Erenguc 14 present a tabu search algorithm for RCPSPs in which a move is made to a neighbourhood solution obtained by starting an activity one time unit early or late. However, infeasible solutions may be generated with this method and this must be taken into account when selecting a move. In many search algorithms for constrained optimization problems including the algorithm of Icmeli and Erenguc, constraints are not specifically considered when generating neighbourhood solutions but incorporated into the objective function as some penalty terms. In this paper, we present a new method that can always give a feasible neighbourhood solution using a simple encoding scheme for representing solutions of the problem, and employ this method in the three search algorithms. We develop a scheduling procedure that can be applied to simulated annealing, tabu search and genetic algorithms for RCPSPs and compare these algorithms with existing heuristics after obtaining appropriate values for various parameters needed in these search algorithms. It is assumed that: the duration of each activity is known and deterministic; each activity requires a

680

Journal of the Operational Research Society Vol. 47, No. 5

known and deterministic level of resource of each type throughout its processing; and the level of resource availability of each resource type remains constant.

THE OVERALL PROCEDURE In the search heuristics suggested in this study, a priority based scheduling method is used for schedule generation. The priority based scheduling method uses priority rules to resolve conflicts among activities competing for common resources. Priority values for the activities are often determined using information of the project network and activities. The method consists of two stages. First, a list of activities is prepared according to a selected criterion or heuristic rule (priority rule) so that higher priorities are given to the activities which are ranked higher or considered urgent by the selected criterion. The list determines the order in which the activities are to be scheduled. Then, the activities are considered one by one in that order for the assignment of their start times and finish times subject to precedence and resource constraints. A procedure for the priority scheduling method is given below. Procedure (Priority scheduling) Step 1. Determine priorities of activities. Set the current time tnow to be 0. Initialize the resource availability vector of which the ith component (r;) denotes the amount of type-i resource available at tnow. Step 2. If all the activities are scheduled, stop. Otherwise, reset the set E of eligible activities, i.e. activities whose predecessors have all been completed by tnow and which require no more than ri of type-i resource for all i. If E is empty, go to Step 3. Otherwise, go to Step 4. Step 3. Reset tnow to the earliest time at which an activity that is currently being processed becomes completed. Update the resource availability vector. Go back to Step 2. Step 4. Select an activity with the highest priority among those in E and schedule it by assigning the start and the finish times and update the resource availability vector. Go back to Step

2. With the above priority scheduling method, a complete schedule can be obtained once the priority values of the activities are given (and if there is a tie-breaking rule for when two or more activities competing for resources have the same priority value). To find a good schedule in a priority scheduling method, a good set of priority values is required. We use search methods to find such priority values. As mentioned earlier, we test three popular search heuristics in this study, simulated annealing, tabu search, and genetic algorithms. Our search heuristics use a solution encoding scheme in which a solution is represented with a string of numbers. Positions of the numbers in the string correspond to indexes of the activities and values of the numbers denote priorities of the corresponding activities. In the search heuristics, we have to generate neighbourhood solutions for a given solution. With the priority scheduling method, a new solution (in the neighbourhood of the current solution) can be obtained by changing the priority values. In our search heuristics, which use the above encoding scheme, neighbourhood solutions are obtained by perturbing the string of numbers with some method, that is, we obtain a string of numbers for a new solution by changing the string associated with the current solution. Since we change only priority values in this method and we can always generate a feasible schedule for all sets of priority values, new solutions are always feasible. An initial schedule or a set of initial schedules are needed to start the search algorithms suggested in this paper. In other words, an initial set or initial sets of priority values are required. In our algorithms, priority values in the initial set(s) are generated randomly from a uniform distribution between 0 and 1. With a random generator, there is very little probability that two or more activities have the same priority value, and therefore a specific tie-breaking rule may not be needed. Advantages of our approach can be summarized as follows. First, search algorithms following this approach do not generate infeasible solutions and hence do not need a feasibility check.

J-K. Lee and Y-D. Kim-Search Heuristics for Project Scheduling

681

Secondly, these algorithms do not need any preliminary information about the problem, which makes them more flexible and more easily adapted to a wide range of problems. Thirdly, the approach can be applied to RCPSPs with objective functions of various forms such as resource leveling and maximization of the net present value as well as duration minimization. Also, it can deal with various types of resource constraints and other constraints if they exist. Finally, its notion is very simple and easy to implement. Project managers would be able to use this approach conveniently in complicated real project management problems. In the project network, activities are represented with nodes, which are numbered (indexed) in the following way. First, nodes are layered so that nodes which have no predecessors are included in layer 1, and immediate successors of layer-i nodes which have all their predecessors already included in layers 1, 2, ... , i are included in layer i + 1. Then, serial numbers are given to the nodes according to the layer, i.e., nodes in a lower layer before those in a higher layer. The reason why the nodes are indexed in this way is explained in the next section. SIMULATED ANNEALING Simulated annealing (SA) algorithms have been used successfully for various combinatorial problems. (Applications of SA to various optimization problems are surveyed in Koulamas et al. 24 .) SA algorithms can be viewed as an enhanced version of local optimization or iterative improvement algorithms, in which an initial solution is repeatedly improved by making small local alterations until no such alteration yields a better solution. SA attempts to avoid entrapment in a local optimum by sometimes accepting a move that deteriorates the value of the objective function g( · ). Starting from an initial solution, SA generates a new solution S' in the neighbourhood of the original solution S. Then, the change in the objective function value, 3. = g(S') - g(S), is calculated. For a minimization problem, if 3. < 0, then the transition to the new solution is accepted. If 3. ~ 0, then the transition to the new solution is accepted with a specified probability usually denoted by the function, exp(- 3./t), where t is a control parameter called the temperature. At a high temperature, most of moves are accepted and the system moves almost freely in its solution space, while it behaves like a greedy search algorithm at a very low temperature. SA algorithms generally start with a high temperature and then the temperature is gradually lowered. At each temperature, the search is carried out for a certain number of iterations, called the epoch length, to attain equilibrium. As described in the previous section, a solution in the SA algorithm suggested in this paper is encoded as a string of numbers representing priorities of activities. The neighbourhood of a solution can be obtained through a certain alteration of the string associated with the current solution. We use the interchange method in which two activities are selected and priorities of the selected activities are exchanged, i.e., two numbers in the string are selected and exchanged with each other. In the method, a first activity is selected randomly and then another (second) activity is selected randomly among ac.;tivities whose indexes are between i1 - rxm and i1 + rxm, where i 1 is the index of the first activity, m is the maximum number of predecessors or successors of an activity over all the activities, and rx is a parameter to be determined. This is to limit the number of alternatives for a move and to prevent the algorithm from searching unnecessary alternatives. Note that with the convention of node numbering described in the previous section, nodes with a large difference in indexes (node numbers) can hardly compete for resources. We set rx = 4 after a series of preliminary experiments. In general, SA algorithms are specified by several parameters and/or methods, namely the initial temperature t 0 , the epoch length L, the rule specifying how the temperature is reduced, and the termination condition. The first three define a so-called cooling schedule, which may affect performance or convergence rate, while the last checks whether the annealing process is frozen. The initial temperature is often chosen in such a way that the fraction of accepted uphill moves in an abbrevated trial annealing run is approximately F 0 , a parameter. That is, after the mean increase of objective function values (3.) is computed with uphill moves only, the value of t 0 is obtained from exp(- 3.jt 0 ) = F 0 . In our algorithm, the number of transitions made in the trial run is set to be equal to the total number of activities.

682

Journal of the Operational Research Society Vol. 47, No. 5

In SA algorithms, the temperature should be decreased in such a way that the cooling process would not take too long. The method which appears to be the most common in the current literature specifies the temperature with tk = ctk- 1> during the kth epoch (k = 1, 2, ...), where cis a constant called the cooling ratio with a value less than 1. Higher cooling ratios correspond to a slower cooling process, and therefore more moves are required before the process is terminated. The epoch length is often set to be equal to the size of neighbourhood (the number of neighbourhood solutions) times a parameter. In this paper, the epoch length is set to be the number of activities times a parameter q. (In this case, a large value may be needed for q since the size of neighbourhood is proportional to the square of the number of activities.) For different methods for decreasing the temperature, see Osman and Potts 25 . As a criterion to terminate the algorithm we use one given by Johnson et al. 26 , in which a counter is incremented by one when an epoch is completed with the fraction (or percentage) of accepted moves less than a predetermined limit (denoted as F min). This counter is reset to 0 when a new incumbent solution is found. The algorithm can be terminated when the counter reaches a given maximum count. In this study, we select values of these two parameters so that the resulting SA algorithm gives good solutions within a reasonable amount of computation time. After a series of preliminary experiments, F min and the maximum count were set to 0.005 and 15, respectively. For more details of SA algorithms, see Johnson et al. 26 • 27 , Osborne and Gillett 28 , Kouvelis and Chiang 29 and Stern 30 . TABU SEARCH Tabu search (TS) is another search procedure designed to escape from a local optimum in pursuit of a global optimal solution. This search procedure incorporates flexible memory functions to forbid transitions in solutions (moves) that reinstate certain attributes of past solutions for the purpose 31 . It starts with an initial feasible solution, and for each solution (S) a new solution (S') is obtained with a function, called a move, that transforms S into S'. A simple TS makes the move to the best neighbouring solution even though it is worse than the given solution, while taking steps to ensure that the search procedure does not re-visit solutions previously visited. This is accomplished by introducing restrictions on possible moves which discourage reversal and/or repetition of selected moves. That is, we define a set of moves that are tabu (forbidden) now, and these moves are sorted in a set called the tabu list. Elements in the list define all tabu moves that cannot be applied to the current solution. The size of the list is bounded by a parameter sL, called the tabu list size. If the number of elements in the list is equal to sL, one must remove an element in it, generally the oldest one, before adding a move to the tabu list. In other words, we use the short term memory when selecting a move at each iteration. Intermediate and long term memory functions may also be incorporated to intensify and diversify the search, which will be discussed later. Glover 31 · 32 and Laguna and Glover 33 •34 give more detailed descriptions and discussions of various aspcts of TS algorithms. The TS algorithm suggested in this paper uses the same methods as those used in the SA algorithm when generating an initial solution and neighbourhood solutions: that is, an initial solution is generated randomly and the interchange method is used for neighbourhood generation. To select two activities in the interchange method, we select a first activity randomly and then select another activity randomly among activities whose indexes are between i 1 - 4m and i 1 + 4m, where i 1 is the index of the first activity, m is the maximum number of predecessors or sucessors of an activity over all the activities. Note that we use the same method and the same value for the parameter IX as those used in the SA algorithm. Differently from SA algorithms, classical TS algorithms generate all neighbourhood solutions for a given solution to select a move. In our TS algorithm, on the other hand, a subset of neighbourhood solutions is generated. Only nN neighbourhood solutions are generated randomly with the above neighbourhood generation method (without replacement) and a move that results in the best solution among non-tabu moves is selected. There may be several ways of defining tabu moves. In our TS algorithm, a move, interchanging priorities of activities i and j, is defined as a tabu move if interchange of i and j has been done in a recent iteration. This short term memory is used to avoid moves that will cancel the effect of

J-K. Lee and Y-D. Kim-Search Heuristics for Project Scheduling

683

moves made recently. In general, short term memory is used to prevent duplication of conditions encountered in the recent past. In many situations, however, short term memory alone is not sufficient to drive the search to regions where improved solutions can be found, and a complementary use of longer term memory may be required. The suggested algorithm uses long term memory with information of how often each move has been made in the past. The use of frequency information is shown to be effective in other settings when applied in cases in which no non-tabu move improves the solution 33 . We use frequencybased long term memory in our algorithm as follows. When a move of interchanging priorities of activities i and j is made, we update fii, which identifies the number of times priority values of activities i and j have been exchanged. The values of fii are used in selecting a move. Among non-tabu moves, a move is selected that results in a new solution S' with the smallest value of g(S') + Pfii, where g(S') is the solution value (project duration) of S', and p = p if g(S') is worse (greater) than the current solution value and p = 0 otherwise. Here, p is a parameter called the penalty factor that must be determined for each application of TS algorithms. Despite all these, our TS algorithm allows a move even though it is a tabu move if it yields a solution that is better than the incumbent solution, the best solution found so far. As a criterion to terminate the algorithm, two rules are commonly used. One is to terminate if the iteration count (the number of moves made) reaches a predetermined number, and the other is to terminate if no improvement has been made of a certain number of iterations. The latter method is used in our algorithm. In the algorithm, a counter is incremented by one when an iteration is completed with no improvement and this counter is reset to 0 when a new incumbent solution is found. The search procedure is stopped when the counter reaches I max, a parameter to be selected. GENETIC ALGORITHM Genetic algorithms (GAs) are search algorithms based on the mechanics of natural selection and natural genetics. They combine the concept of survival of the fittest among solutions with a structured yet randomized information exchange and offspring creation 35 . They have been used for many combinatorial optimization problems such as assembly line balancing problems 36 , scheduling problems 37 , and assignment problems 38 . A GA does not evaluate and improve a single solution but analyses and modifies a set of solutions simultaneously. The ability of a GA, operating on many solutions simultaneously and gathering information from all current solutions to direct the search, reduces the possibility of being trapped in a local optimum. In a simple GA, a solution (called chromosome) is represented by a string of numbers (genes). A chromosome's potential as a solution is measured by a fitness function which evaluates a chromosome with the objective function value. A judiciously selected set of chromosomes is called a population and the population at a given time is called a generation. Usually, the number of chromosomes remains fixed from generation to generation. The fundamental mechanism of GAs consists of three main operations: reproduction, crossover and mutation. Chromosomes resulting from these operations applied to those in the current population, often known as offspring or children, form a population of the next generation. The procedure of generating a next generation is repeated for a certain number of times, which can be determined in various ways. An application of a GA is characterized by an encoding scheme, an initial population, genetic operators (reproduction, crossover, and mutation), a fitness function, and a termination rule. We use the same encoding scheme as in the previous two search algorithms, that is, a solution is encoded as a string of numbers representing priorities of activities. To initiate the algorithm, an initial population can be obtained either with a random method or a well-adapted method. Since a random method is used in the other two search algorithms tested in this study, we generate an initial population randomly in our GA. (In fact, it is reported that use of well-adapted populations provides little advantage in spite of a little faster convergence 39 .) As mentioned earlier, three operations are performed to produce a next generation with better chromosomes. Reproduction is a process in which individual chromosomes are copied according to their fitness, i.e., chromosomes with higher fitness values tend to have more of their copies at

684

Journal of the Operational Research Society Vol. 47, No. 5

the next generation. This is done by randomly selecting and duplicating a chromosome with a probability which is proportional to the fitness value of the chromosome. Since our objective is to minimize the project duration while fitness is usually maximized in GAs, we cannot use the objective function values directly for the fitness values. Various methods can be used to transform an original objective function value to a fitness value, such as subtracting the objective function value from an arbitrary large constant or taking the reciprocal of the objective function value. In our implementation of the algorithm, we compute the fitness value of solution i with exp(- hv;), where V; is the objective function value of the solution and h is a parameter to be chosen to make the fitness value lie in a particular range (for example between 0 and 1). In our algorithm, it was set at 0.004. After the above reproduction operation is performed, the operation of crossover is applied to the resulting set of chromosomes. Crossover introduces new chromosomes by recombining current ones. Two chromosomes are randomly selected from members of newly reproduced chromosomes in the mating pool, and a crossover point is selected randomly in between two consecutive numbers in the string. Then, with a probability called crossover rate, two new chromosomes can be created by swapping substrings after the crossover point. This process is repeated until all chromosomes in the current population are mated. Although reproduction and crossover effectively search and recombine desirable characteristics in solutions, occasionally they may lose some potentially useful genetic material. In genetic algorithms, mutation is performed to compensate for such an irrecoverable loss. It introduces a random change to a chromosome by altering the value of a gene. Mutation is often applied to each of the genes of chromosomes in the current population with a probability called the mutation rate. It is noted that a low mutation rate gives good results in empirical studies 35 • In our implementation, for each chromosome two genes are selected randomly and their values are changed with each other with a probability equal to the mutation rate. This is applied to all chromosomes in the population. As a termination condition, we can limit the total number of generations. Alternatively, we can terminate the algorithm when an incumbent solution cannot be changed (a solution cannot be improved) for a predetermined number of generations, which is applied in this study. In our algorithm, the termination condition is controlled so that the CPU time for GA is approximately the same as that for SA or TS for a same size problem. Among various parameters, the number of generations allowed (without an improvement in solutions) before termination and the number of chromosomes in a population affect the computation time, and values for these parameters will be determined through a test on a number of problems. Values of other parameters, crossover rate and mutation rate, which significantly affect performance of a GA, are to be selected to obtain the best solution within the computation time set by the former two parameters. COMPUTATIONAL EXPERIMENTS In this section, we test the search algorithms on randomly generated test problems and 110 problems assembled by Patterson40 . We generated the random test problems as follows. Since it is known that performance (solution quality and CPU time) of an algorithm is affected by various factors that characterize RCPSPs such as the number of activities, the number of resources types and the tightness of resource constraints, several levels for these factors may have to be considered. We select three levels for the number of activities (50, 100 and 150), three levels for the number of resource types (1, 4 and 7), and two levels for the tightness of resource constraints (tight and not-tight). We generated ten problems for each of all combinations of the levels for the three factors. In each test problem, time duration for an activity was generated from a discrete uniform distribution with range [1, 10]. For each activity, the numbers of predecessors and successors were selected randomly from a discrete uniform distribution with range [1, 5] except for the source and sink nodes, and predecessors of each activity were selected randomly from activities with indexes smaller than the activity. After this process was completed, the activities were re-indexed according to the node-numbering method explained earlier. In the test problems with tight resonance

J-K. lee and Y-D. Kim-Search Heuristics for Project Scheduling

685

constraints, the resource requirement of each activity was generated from a uniform distribution with range (0, 1), while the range of (0, 2/3) was used in those with not-tight constraints. Amounts of resource available and resource required are assumed to be normalized so that the total resource availability is 1 unit for each resource. For better implementation of the search algorithms described in the previous sections, we have to select appropriate values for the parameters. Therefore, we first select a set of parameter values that give good results for each algorithm through a series of computational experiments on a subset of (small-sized) test problems generated with the method described above. Then, the three algorithms with the most appropriate parameter values are compared with each other and with existing heuristics. The existing heuristics included in the comparison on the randomly generated problems are the minimum slack method, the iterative method, and the SEARCH method. The minimum slack method uses slack, i.e., the difference between late start time and early start time determined by the critical path method, for priority determination. An activity with the least slack time is selected for scheduling when there is resource conflict. The iterative method, which is developed by Li and Willis 9 , employs a new procedure called backward scheduling which schedules each activity one by one starting from the last activity to the first activity. Forward and backward schedulings are carried out iteratively until there is no further improvement in the project completion time. The SEARCH method developed by Khattab and Choobineh 16 is a composite heuristic which includes eight priority rules. This method schedules activities with each of the rules and selects the best of the schedules obtained. (For more detail of the algorithms, see the corresponding references.) They are selected because they are simple to use (minimum slack method) or claimed to be better than other algorithms. In the test, we use the project duration (makespan) as the performance measure with which the heuristics are compared. All the heuristics were coded in C language and run on a personal computer with an 80486 processor. Selection of parameter values

First, we select values of the parameters needed in the SA algorithm. These are F 0 , a parameter for setting initial temperature, c, the cooling ratio, and q, a parameter for setting the epoch length. See the third section for detailed description of how these parameters are used in the SA algorithm. We tested four levels for F 0 (0.3, 0.5, 0.7 and 0.9), three levels for c (0.95, 0.97 and 0.99) and three levels for q (10, 15 and 20). After a series of computational experiments, it was found that the SA algorithm with F 0 = 0.7, c = 0.99, and q = 20 gave better solutions than the others without requiring much longer computation time, and hence it was selected for the comparison with other algorithms. Next, a few values are tested for parameters to be specified in the TS algorithm. These parameters are the size of the tabu list (sL), the number of neighbourhood solutions (nN) to be generated before 'Selecting a move, the penalty factor (p) to take into account the long term memory when selecting a move, and the maximum iteration count for termination (I maJ First, a value for I max was selected so that the computation time for the TS algorithm became approximately the same as that for the SA algorithm for the same sized problems. The selected value for /max was 100. Given this termination condition, we tested five levels for sL (3, 4, 6, 8 and 10), three levels for p (5, 10 and 15) and three levels for nN (n, 2n and 4n), where n is the number of activities. We selected the combination (sL, p, nN) = (6, 10, 2n) for the TS algorithm to be compared with other algorithms. Parameter values required in the GA were determined as follows. First, the parameter for the termination condition and the population size were determined in such a way that the resulting GA required approximately the same computation time as the SA algorithm and gave reasonably good solutions. As a result, the population size was set at 20 and the algorithm was terminated when the incumbent solution could not be improved for 30 generations. To find the most appropriate values for the crossover rate and the mutation rate in this algorithm, we tested four levels for the crossover rate (0.6, 0.7, 0.8 and 0.9) and five levels for the mutation rate (0.1, 0.2, 0.3, 0.4 and 0.5) and selected 0.7 and 0.4 for the two rates.

686

Journal of the Operational Research Society Vol. 47, No. 5

Comparison of the algorithms

We compared the three search algorithms and three existing algorithms on 180 test problems generated with the method described earlier. The algorithms were compared with each other by a measure called the (relative) performance ratio, which is the ratio of the solution of each algorithm to that of the minimum slack method. Solutions of the minimum slack method are used as benchmark solutions since they are clearly the worst among solutions of all algorithms. Also used as measures were the average CPU time and the number of problems each algorithm yielded the shortest duration, i.e., frequency with which each algorithm gave the best solution. Table 1 shows performance of the heuristics for different numbers of activities. Relative solution quality was not affected by the number of activities. As expected, the search heuristics, SA, TS and GA consistently outperformed the existing heuristic methods in all problem sizes. Among the search heuristics, GA gave slightly worse solutions than the other two, especially for large-sized problems. It may be due to premature termination of the genetic algorithm in large sized problems as can be seen from computation time required for the algorithm, which is shorter than those for the others. Different results may be obtained if a different termination criterion is used in the GA. Computation time increases almost linearly in most algorithms, as the number of activities increases. TABLE 1. Performance of the methods for different numbers of activities

Number of activities

Iterative

SEARCH

SA

TS

GA

50

NBS APR SDPR CPUT

13 0.970 0.038 0.21

16 0.970 0.035 0.23

50 0.945 0.043 7.01

45 0.946 0.041 8.88

38 0.950 0.040 8.01

100

NBS APR SDPR CPUT

6 0.969 0.031 0.22

7 0.975 0.029 0.78

43 0.944 0.040 12.28

43 0.944 0.040 17.58

23 0.952 0.033 12.05

150

NBS APR SDPR CPUT

2 0.972 0.025 0.29

4 0.973 0.024 1.89

34 0.943 0.037 20.58

38 0.943 0.036 26.66

12 0.958 0.028 15.97

overall

NBS APR SDPR CPUT

21 0.970 0.031 0.24

27 0.972 0.030 0.97

127 0.944 0.040 13.29

126 0.945 0.039 17.71

73 0.953 0.034 12.01

NBS: number of problems (frequency) for which each algorithm found the best solution. APR: average performance ratio (ratio of the solution of each algorithm to that of the Minimum Slack algorithm). SDPR: standard deviation of performance ratio. CPUT: average CPU times given in seconds.

Similar results can be seen from Table 2, which shows the effects of the number of constraints on the performance of the algorithms. CPU time increases almost linearly as the number of constraints increases and the search heuristics consistently outperformed the existing heuristics. However, the relative superiority of the search heuristics diminishes as the number of constraints increases. A reason may be that there is not much possibility of improving solutions when there are many constraints (feasible region for solutions is small). Note that in the priority scheduling method, which is employed in the search algorithms, there are not many alternatives each time when selection has to be made, if there are many constraints. Therefore, improvement methods including the iterative method and the SEARCH method may not be able to improve solutions much. Similar arguments can be made for the results obtained from using different levels for the tightness of the constraints, which are shown in Table 3. That is, there are not many alternatives when selection has to be made if the constraints are tight, so there is less opportunity for improvements. This table shows more clearly that good improvement methods can have higher relative superiority when constraints are not tight. It shows not only the superiority of existing improve-

J-K. Lee and Y-D. Kim-Search Heuristics for Project Scheduling

687

TABLE 2. Performance of the methods for different numbers of constraints Number of constraints

Iterative

SEARCH

SA

TS

GA

NBS APR SDPR CPUT

7 0.964 0.036 0.13

7 0.966 0.033 0.91

38 0.936 0.036 10.78

41 0.938 0.036 10.58

27 0.945 0.034 9.20

4

NBS APR SDPR CPUT

2 0.970 0.028 0.18

4 0.974 0.027 0.98

42 0.944 0.035 13.32

46 0.943 0.036 17.95

17 0.956 0.027 12.07

7

NBS APR SDPR CPUT

12 0.977 O.Q28 0.41

16 0.977 0.028 1.02

47 0.951 0.046 15.77

39 0.952 0.044 24.58

29 0.960 0.039 14.75

ment methods over the minimum slack method but also the superiority of better improvement methods over less effective methods in cases of loose constraints. Longer CPU time was required when constraints were not tight since more improvements were made. TABLE 3. Performance of the methods for different levels of tightness of constraints Tightness of constraints

Iterative

SEARCH

SA

TS

GA

tight

NBS APR SDPR CPUT

16 0.982 0.028 0.15

23 0.982 0.027 0.98

68 0.965 0.032 10.72

69 0.965 0.032 15.54

50 0.969 0.030 11.93

not tight

NBS APR SDPR CPUT

5 0.959 0.031 0.33

4 0.963 0.029 0.96

59 0.923 O.Q35 15.86

57 0.924 0.035 19.87

23 0.938 0.031 12.09

Finally, we compared the two best performing search heuristics on the 110 test problems assembled by Patterson 41 . The result is given in Table 4. Since the most appropriate values for the parameters needed in the search algorithm may vary according to the characteristics of problem data and the data of the 110 problems are different from those randomly generated in this study, we selected different values for the parameters. In the table, SA2 and TS2 denote the algorithms for which the parameter values were newly selected for the 110 problems while SAl and TSl are those with the same set of parameter values used in the previous tests. Also included in the table are the results of the heuristics developed by Bell and Han 2 and Simpson and Weiss 22 . (These are recently developed heuristics tested on Patterson's 110 problems.) As can be seen from the table, our search algorithms without additional tuning process (SAl and TSl) outperformed both of the existing heuristics. After being appropriately tuned (SA2 and TS2) for the test problems, the two search algorithms worked much better than the existing heuristics. CPU times of the search algorithms required for a problem did not exceed 15 seconds (average was approximately 6 seconds). Overall, search heuristics tested in this study gave better solutions than the existing heuristics with allowable additional computation time. With the current implementation of the search algoTABLE 4. Comparison of heuristics on Patterson's 110 problems

ADP SDDP NOS

SAl

SA2

TSl

TS2

Bell & Han 2

Sampson & Weiss 22

0.95% 1.96% 81

0.57% 1.56% 91

1.05% 1.85% 76

0.75% 1.71% 81

2.6% 3.1% 49

1.98% 2.62% 61

ADP: average deviation percentage from optimal project durations. SDDP: standard deviation of the deviation percentages from optimal project durations. NOS: number of problems (frequency) for which each algorithm found the optimal solution.

688

Journal of the Operational Research Society Vol. 47, No. 5

rithms, the SA algorithm was slightly better than the other two search algorithms. The SA algorithm gave almost the same solution quality as the TS algorithm with a little shorter computation time on average and gave better solutions than GA with almost the same computation time. However, one should note that the search algorithms as tested in this study may not be appropriately tuned for comparison among each other. Especially, termination conditions for the TS and GA were devised so that the resulting computation times for the three algorithms were approximately the same in a number of test problems included in the preliminary experiments. Moreover, they were set so that the computation time would be short enough for the algorithms to be able to deal with large size problems, and therefore they may have been stopped too early in many test problems. If more appropriate parameter values are used, the search algorithms can be improved and different results may be obtained. Also, other methods can be used to compare the search heuristics. For example, we can compare the search methods with solutions obtained within exactly the same computation time, with the speed of improvement in solutions, or with some function of the solution quality and computation time. CONCLUDING REMARKS In this paper, we considered the resource constrained project scheduling problem (RCPSP) and suggested three search heuristics for the problem, simulated annealing, tabu search, and genetic algorithms. A new encoding scheme was used in the search algorithms. In the scheme, a solution is represented with a string of numbers each of which denotes priority of its corresponding activity. From these priorities, a complete schedule can be obtained with the priority scheduling method. Unlike other improvement methods for RCPSPs, our algorithms can always generate feasible neighbourhood solutions. The search algorithms can be considered flexible in that they can deal with problems with objective functions of various forms and complex resource constraints. The three search algorithms were compared with existing algorithms for RCPSPs after a series of tests for finding the most appropriate values for various parameters needed in the algorithms. Results of the comparison on randomly generated test problems showed that our search algorithms outperformed existing heuristics, the minimum slack mehod, the iterative technique and the SEARCH algorithm. The results also revealed that search heuristics could more easily improve solutions in problems with less tight resource and/or precedence constraints. Also, results of the test on Patterson's 110 problems showed that solutions obtained from SA and TS were close to the optimal solutions. Our research can be extended in several ways. Various problems can be considered such as cash flow maximization problems and problems with additional constraints other than the resource constraints. While it may be very difficult to solve these problems with an analytical approach, the search heuristics can be adapted easily to such situations. They can be applied to multi-project scheduling problems as well since they do not depend on the network topology. Also, other methods can be used for neighbourhood generation in the search heuristics. Our search heuristics have limits in that only priorities of activities are changed to improve solutions. That is, we can further improve solutions by delaying the start of an activity on purpose although it can be started, since the set of solutions from the priority scheduling method does not always include an optimal solution as the non-delay schedules are not dominant in job shop scheduling. (See Sprecher et a/. 41 for a definition of the non-delay schedule in RCPSPs.) Acknowledgements-We would like to thank Colin Bell, Eric Demeulemeester, and Willy Herroelen for providing us with a data file for Patterson's 110 test problems and their branch and bound code for the resource constrained project scheduling problems.

REFERENCES 1. S. E. ELMAGHRABY (1995) Activity nets: a guided tour through some recent developments. Eur. J. Opl Res. 82, 383-408. 2. C. E. BELL and J. HAN (1991) A new heuristic solution method in resource-constrained project scheduling. Naval Res. Logist. 38, 315-331. 3. C. E. BELL and K. PARK (1990) Solving resource constrained project scheduling problems by A• search. Naval Res. Logist. 37, 61-84. 4. F. F. BOCTOR (1990) Some efficient multi-heuristic procedures for resource-constrained project scheduling. Eur. J. Opl Res. 49, 3-13.

J-K. Lee and Y-D. Kim-Search Heuristics for Project Scheduling

689

5. N. CHRISTOFIDES, R. ALVAREZ-VALDES and J. M. TAMARIT (1987) Project scheduling with resource constraints: a branch and bound approach. Bur. J. Opl Res. 29, 262-273. 6. E. W. DAVIS and J. H. PATTERSON (1975) A comparison of heuristic and optimum solutions in resource-constrained project scheduling. Mgmt Sci. 21, 944-955. 7. E. DEMEULEMEESTER and W. HERROELEN (1992) A branch-and-bound procedure for the multiple resource-constrained project scheduling problem. Mgmt Sci. 38, 1803-1818. 8. E. A. ELSAYED and N. Z. NASR (1986) Heuristics for resource-constrained scheduling. Int. J. Prod. Res. 24, 299-310. 9. K. Y. LI and R. J. WILLIS (1992) An iterative scheduling technique for resource-constrained project scheduling. Bur. J. Opl Res. 56, 370-379. 10. J. H. PATTERSON, B. F. TALBOT, R. SLOWINSKI and J. WEGLARZ (1990) Computational experience with a backtracking algorithm for solving a general class of precedence and resource-constrained scheduling problems. Bur. J. Opl Res. 49, 68-79. 11. B. TALBOT (1982) Resource-constrained project scheduling with time-resource tradeoffs: the nonpreemptive case. Mgmt Sci. 28, 1197-1210. 12. T. YANG, J. P. IGNIZIO and J-S. SoNG (1989) An exchange heuristic algorithm for project scheduling with limited resources. Eng. Optimization 14, 189-205. 13. R. H. DoERSCH and J. H. PATTERSON (1977) Scheduling a project to maximize its present value: a zero-one programming approach. Mg111t Sci. 23, 882-889. 14. 0. IcMELI and S. S. ERENGUC (1994) A tabu search procedure for the resource constrained project scheduling problem with discounted cash flows. Comps and Opns Res. 21, 841-853. 15. A. H. RussELL (1970) Cash flows in networks. Mgmt Sci. 16, 357-373. 16. M. M. KHATTAB and F. CHOOBINEH (1991) A new approach for project scheduling with a limited resource. Int. J. Prod. Res. 30, 185-198. 17. I. KuRTULUS and E. W. DAVIS (1982) Multi-project scheduling: categorization of heuristic rules performance. Mgmt Sci. 28, 161-172. 18. J. H. PATTERSON (1973) Alternate methods of project scheduling with limited resources. Naval Res. Logist. Q. 20, 767-784. 19. R. WALLACE and W. HALVERSON (1992) Project management: a critical success factor or a management fad. Ind. Eng. 24 (4), 8-50. 20. J. DE WIT and W. HERROELEN (1990) An evahiation of microcomputer-based software packages for project management. Bur. J. Opl Res. 49, 102-139. 21. J. BLAZEWICZ, J. K. LENSTRA and A. H. G. RINNOOY KAN (1983) Scheduling projects to resource constraints: classification and complexity. Discrete Applied Maths 5, 11-24. 22. S. E. SAMPSON and E. N. WEISS (1993) Local search techniques for the generalized resource constrained project scheduling problem. Naval Res. Logist. 40, 665-675. 23. 0. 0Guz and H. BALA (1994) A comparative study of computational procedures for the resource constrained project scheduling problem. Bur. J. Opl Res. 72, 406-416. 24. C. KouLAMAS, S. R. ANTONY and R. JEAN (1994) A survey of simulated annealing applications to operations research problems. Omega 22, 41-56. 25. I. H. OsMAN and C. N. PoTTS (1989) Simulated annealing for permutation flow-shop scheduling. Omega 17, 551-557. 26. D. S. JOHNSON, C. R. ARGON, L. A. McGoECH and C. SCHEVON (1989) Optimization by simulated annealing: an experimental evaluation; Part I, graph partitioning. Opns Res. 37, 865-892. 27. D. S. JOHNSON, C. R. ARAGON, L. A. McGoECH and C. SCHEVON (1991) Optimization by simulated annealing: an experimental evaluation; Part II, graph coloring and number partitioning. Opns Res. 39, 378-406. 28. L. J. OsBORNE and B. E. GILLETT (199) A comparison of two simulated annealing algorithms applied to the directed Steiner problem on networks. ORSA J. Computing 3, 213-225. 29. P. KouVELIS and W. CHIANG (1992) A simulated annealing procedure for single row layout problems in flexible manufacturing systems. Int. J. Prod. Res. 30, 717-732. 30. J. M. STERN (1992) Simulated annealing with a temperature dependent penalty function. ORSA J. Computing 4, 311319. 31. F. GLOVER (1989) Tabu search-part I. ORSA J. Computing 1, 190-206. 32. F. GLOVER (1990) Tabu search-part II. ORSA J. Computing 2, 4-32. 33. M. LAGUNA and F. GLOVER (1993) Bandwidth packing: a tabu search approach. Mgmt Sci. 39,492-500. 34. M. LAGUNA and F. GLOVER (1993) Integrating target analysis and tabu search for improved scheduling systems. Expert Systems with Applications 6, 287-297. 35. D. E. GOLDBERG (1989) Genetic Algorithms in Search, Optimization and Machine Learning. Addison-Wesley, Reading, MA. 36. E. J. ANDERSON and M. C. FERRIS (1994) Genetic algorithms for combinatorial optimization: the assembly line balancing problem. ORSA J. Computing 6, 161-173. 37. C. G. GUPTA, Y. P. GUPTA and A. KuMAR (1993) Minimizing flow time variance in a single machine system using genetic algorithms. Bur. J. Opl Res. 70, 289-303. 38. G. LEVITIN and J. RUBINOVITZ (1993) Genetic algorithm for linear and cyclic assignment problem. Computers and Opns Res. 20, 575-586. 39. G. E. LIEPINS, M. R. HILLIARD, M. PALMER and M. MoRROW (1987) Greedy genetics. In Proc. 2nd Int. Conf. on Genetic Algorithms and their Applications (J..J. GREFENSTETTE, Ed.), pp 90-99. Lawrence Erlbaum Associates, Hillsdale, NJ. 40. J. H. PATTERSON (1984) A comparison of exact approaches for solving the multiple constrained resources, project scheduling problem. Mgmt Sci. 30 (7), 854-867. 41. A. SPRECHER, R. KouscH and D. ANDREAS (1995) Semi-active, active and non-delay schedules for the resourceconstrained project scheduling problem. Bur. J. Opl Res. 80, 94-102.

Received March 1995; accepted October 1995 after one revision