A Hybrid Multiobjective Genetic Algorithm for Robust Resource-Constrained Project Scheduling with Stochastic Durations

We study resource-constrained project scheduling problems with perturbation on activity durations. With the consideration of robustness and stability of a schedule, we model the problem as a multiobjective optimization problem. Three objectives—makespan minimization, robustness maximization, and stability maximization—are simultaneously considered. We propose a hybrid multiobjective evolutionary algorithm H-MOEA to solve this problem. In the process of the H-MOEA, the heuristic information is extracted periodically from the obtained nondominated solutions, and a local search procedure based on the accumulated information is incorporated. The results obtained from the computational study show that the proposed approach is feasible and effective for the resource-constrained project scheduling problems with stochastic durations.


Introduction
As for its practical importance in a wide range of real-world application areas, project scheduling problems PSPs have received considerable attention in the operation research field.PSPs usually address two matters: resource and time 1 .Given scarce resources, the problem is modelled as a resource-constrained project scheduling problem RCPSP , which is a general type of PSPs that contains job-shop, flow-shop, and open-shop problems as special cases 2 .As RCPSP is an NP-hard problem 1 , heuristics or metaheuristics are preferred to solve it.For details about the models and methods of the RCPSP, readers are referred to several comprehensive surveys 3-7 .
It is very common that the execution of a schedule will be in the presence of uncertain factors.With the consideration of uncertainty, the problem can be modelled as stochastic RCPSP.One of the main avenues to solve the stochastic project scheduling problems is to develop robust schedules that can deal with uncertain circumstances.Robustness maximization is often taken as an objective in stochastic project scheduling problems.This methodology is referred to as robust scheduling.Van de Vonder et al. 8 presented concepts of quality robustness and solution robustness and addressed the issue of achieving a tradeoff between these two robustness measures.
For stochastic RCPSP, the conventional objective, makespan minimization, and the new introduced objective, robustness maximization, are in conflict with each other to a certain degree.In other words, stochastic RCPSP is inherently a multiobjective optimization problem.It is interesting to note that although stochastic project scheduling has considerable importance in practical applications, the works on solving this problem by multiobjective evolutionary approaches are rather scant 9 .Al-Fawzan and Haouari 10 proposed to measure robustness of a given schedule via the total amount of free slack of all activities.The free slack was also referred to as floating time by Abbasi et al. 11 .Chtourou and Haouari 12 followed the approach of Al-Fawzan and Haouari but assigned weights to the free slack of an activity according to the number of its successors and/or the sum of its required resources.Lambrechts et al. 13 studied the stochastic RCPSP with uncertainty resource availability.The authors developed proactive strategies to build a robust schedule that meets the project deadline and minimizes the schedule instability cost, indicated as the deviations between the planned and the actually realized activity starting times during the execution.Ballestín and Leus 14 investigated various objective functions related to timely project completion of RCPSP with stochastic activity durations.The authors developed a Greedy Randomized Adaptive Search Procedure GRASP to produce high-quality solutions.Ashtiani et al. 15 introduced a new class of scheduling policies for solving RCPSP with stochastic activity durations.The authors underlined the value of preprocessing in stochastic scheduling.In the process, a subset of sequencing choices at the beginning of the planning horizon was made and the rest of the scheduling decisions to future points were relegated in time.
In the present paper, we focus on the RCPSP with stochastic activity durations.The term "robustness" used in the remainder of this paper refers to quality robustness, which represents the schedule's ability to cope with the perturbation on activity durations.We measure the robustness as the difference between the planned makespan and the expectation of the makespan in the real execution under stochastic environment.By considering a schedule's stability, the problem is modelled as a three-objective optimization problem, where makespan minimization, robustness maximization, and stability maximization are simultaneously taken into account.A hybrid multiobjective genetic algorithm incorporating a local search procedure is proposed to solve the problem.

Deterministic Resource-Constrained Project Scheduling
RCPSP is the optimization problem to schedule the interrelated activities tasks, operations of a project with consideration of minimizing the makespan while given precedence constraints among activities and resource constraints are satisfied.
RCPSP can be formally described as follows: a project consists of a set of activities N N 0 , N 1 , . . ., N n , N n 1 , where each activity has to be executed without interruption.N 0 and N n 1 are dummy activities which, respectively, represent the start and end activity.The number of available resource types is K and the vector of resources availability is denoted as R R 1 , R 2 , . . ., R K , where R k k 1, . . ., K represents the availability of each resource type k in each time period.For each activity N i i 1, . . ., n , the duration is denoted as d i , and the requested resources vector over the entire duration is represented as r i r i1 , . . ., r iK .Note that for dummy activities N 0 and N n 1 , the durations and the required resources are both 0. An activity-on-node direct acyclic graph, denoted as G N, A , is employed to represent the precedence relations amongst the activities.The direct successors and predecessors of an activity N i are denoted as Succ i and Pred i , respectively.The activities are interrelated by two types of constraints: precedence constraint and resource constraint.The former constraint indicates that each activity should be scheduled after all of its predecessor activities are accomplished.The latter constraint ensures that for each type of resource, the amount of being occupied cannot exceed the total availability at any time point.Each activity has start and finish times, denoted as st i and ft i , respectively.A schedule S consists of the vector of start time of the nondummy activities, denoted as S st 1 , . . ., st n .The makespan of the schedule S is denoted as M S .

Resource-Constrained Project Scheduling with Stochastic Duration
However, in reality, it is quite common that execution of a schedule may be confronted with uncertain factors.Although there might be multiple uncertain factors affecting a schedule's execution, we just consider the situation that the activity durations are subject to perturbation, which is also a most common uncertainty addressed in the related literature.For example, in an airlift task scheduling, the extreme bad weather will possibly have a negative impact on the execution of the tasks.
Under the effect of perturbation, the durations of some activities will be delayed and, finally, the perturbation might result in delay to the whole schedule.The delay amount is corresponding to the magnitude of perturbation, which is depicted by the following two elements: perturbation strength and perturbation period.The former element represents the strength of a perturbation, which is denoted as St and is normalized in the interval 0, 1 .The closer St is to 1, the more strength the perturbation is and vice versa.The perturbation period is used to indicate the duration of perturbation and is determined by the occurrence time and the end time of the perturbation, denoted as OT and ET , respectively.In the presence of perturbation, actual start time and finish time of an activity are represented as st u i and ft u i , respectively.The actual duration is denoted as d u i and is calculated as follows: where Round • is the function of obtaining the integer part of a real value.We assume that if any part of the execution of an activity is affected by the perturbation, then the whole activity is considered to be affected.
In this paper, we assume that perturbation strength, perturbation occurrence time and end time follow normal probability distribution, which is a most often used distribution in real-world application.In a specific realization of the perturbation, denoted as µ, the actual makespan of a schedule S is denoted as M µ S .In the presence of perturbation, robust predictive schedules are preferred to cope with the impact on the durations.Here we indicate the robustness of a schedule under a specific realization of uncertainty, denoted as Rob µ S , by the difference between the deterministic makespan before disruption, and the actual makespan after execution 13, 16 .The robustness under a specific realization of the uncertainty is formally defined as follows: In the measurement of robustness under uncertainty with a probability distribution, an analytical closed form of the effective fitness function is usually not available 17 .In order to tackle this issue, we employ simulation method, Monte Carlo in concrete, to approximate the effective function of the robustness, denoted as the expectation value with a sample size N µ .Then, the robustness value of a schedule in the presence of uncertainty, Rob S , is calculated as follows: where µ m is the mth realization of the uncertainty in the sample space.
In the robust design techniques, another important issue in need to be considered is the stability of the solutions.We employ the standard deviation to indicate the stability measure of the solutions, denoted as Std Rob µ S and calculated as follows: With the consideration of both magnitudes of robustness measure and stability of the solution, the problem can be modelled as an uncertainty-based multiobjective design problem 18 , instead of combining these two measures by a linear weighted sum method, such as in 16 .What should be noted is that the difference between predesigned and actual makespan by itself has little relevance since a lower level of difference can be easily achieved by predesigning a schedule with a longer makespan.In 16 , the authors defined the robustness of a schedule as the linear combination of the expected delay and the expected makespan via weighted sum method.Lambrechts et al. 13 just measured the robustness as a weighted deviation of the starting times.However, the weight coefficient or weight vector varies in different decision makers.A final solution is much dependent on the decision makers' preferences.Thus, we model the RCPSP with uncertain durations as a multiobjective optimization problem where the concurrently considered three objectives are makespan minimization, robustness maximization, and stability maximization standard deviation minimization .The optimization model is formally given as follows: obj.: 1 min f 1 M S , r ik t ≤ R k , ∀t and k 1, 2, . . ., K.

2.5
In the above optimization model, objective 1 represents the minimization of makespan, objectives 2 and 3 , respectively, indicate the maximization of robustness and stability.Constraints 1 and 2 are precedence constraints which ensure the precedence feasibility in both deterministic and uncertain environment, that is, the activity N i cannot be started until all of its predecessors are completed.Constraint 3 is a resource constraint, indicating that the total amount of the resources occupied at any time point cannot exceed the resource availability.
Within the framework of multiobjective optimization, the concept of dominance is defined as follows: an individual I 1 dominates individual I 2 if I 1 does at least as well as I 2 on all objectives, and strictly better on at least one.Thus, in this paper, we formally define the schedule dominance as follows.
with strict inequality holding for at least one objective measure.
The nondominance concept can be defined as follows.

Definition 2.2. A schedule S * is an individual in the population. It is said to be a nondominated schedule if it is not dominated by any other individual in the population.
For the sake of clarity, we list the notations which are used to represent the RCPSP as follows.
R k : the availability of kth type of resource, k 1, . . ., K. r ik : the amount of the kth type of resource requested by activity N i , i 1, . . ., n.    Mathematical Problems in Engineering st u i : the start time of activity N i in the presence of perturbation.ft u i : the finish time of activity N i in the presence of perturbation.S: the schedule, S st 1 , . . ., st n .
M S : the makespan of schedule S in the deterministic situation.
M µ S : the makespan of schedule S under a specific realization of uncertainty µ.
St: the perturbation strength.
OT : the occurrence time of the perturbation.
ET : the end time of the perturbation.

Methodology
As the three objectives are in conflict with each other, a multiobjective optimization approach is required to solve the problem.To do this, we employ one of the classic multiobjective evolutionary algorithms, named NSGA-II 19 , as the basic optimizer for the problem.In the procedure of NSGA-II, the parent population and the offspring are combined and sorted in order to generate the next-generation population.The combined population is classified into different ranks of nondomination via a nondominated sorting mechanism.In order to maintain the diversity among nondominated solutions, a crowding-distance assignment mechanism is used in the selection of solutions in the same nondomination rank.However, the chromosome and genetic operators are redesigned for the RCPSP.

Chromosome Representation
In the RCPSP, a solution must contain the information about the sequence of activity execution.We employ the activity sequence proposed by Hartmann 20 to represent an individual, I, denoted as follows: The activity sequence is a precedence feasible permutation of the set of activities and the element at each index is the ID of the activity to be executed.
In order to generate a schedule, the Serial Schedule Generation Scheme SSGS is employed to perform this transformation from chromosome representation.SSGS consists of n stages for a network of n nondummy nodes, and in each stage one activity from the activity sequence is selected and executed at its earliest possible time, considering precedence and resource constraints 21 .Algorithm 1 shows the pseudocode of the procedure of SSGS.

Population Initialization
In the operator of population initialization, the individuals are generated serially via a twostep procedure 20, 21 : the first step is identifying the eligible activities, then one of the them is randomly selected in the second step.At a position in the activity list, an activity is said to i 1 while i ≤ n do st i 0, st i 0, ft max 0 Obtain the maximum finish time of the predecessors: ft max max{ft j | N j ∈ Pred i } Determine the earliest time point, st i , which satisfying the resource constraints and st i ≥ ft max .Set st i st i .end while Algorithm 1: Pseudocode of the Serial Schedule Generation Scheme SSGS .
be eligible if all of its predecessors have been scheduled.Formally, let P i be the position of activity N i in the activity list, if we look at the position q q 0, 1, . . ., n 1 , an activity N q * is said to be an eligible activity at the position q if P j < q, for all j, N j ∈ Pred q * .The set which consists of all eligible activities at position q is denoted as EA q .

Genetic Operators
The crossover and the mutation operators are crucial for the performance of the genetic algorithms.The crossover operator combines two parents to generate two children.Similar with the situation in nature, crossover ensures that the children maintain some good characteristics of their parents.We employ a two-point position-based crossover operator 22, 23 in this research.For each crossover operator, two integer random numbers, r 1 and r 2 , r 1 < r 2 , are selected from the interval 1, n .Usually, the two-point position-based crossover for scheduling problem can be divided into two versions 22 : in the first version, the activities outside the two selected points are inherited from one parent to the child, while in the other version, the set of activities between two randomly selected points is inherited from one parent to the child.In this research, we employ the first version of crossover operator.Let P 1 and P 2 be the two selected parents, CH1 and CH2 be the two new generated children.It is clear that the children consist of three parts separated by the two crossover points r 1 and r 2 .For example, CH1 is constructed as j , where the and the third part j CH1 r 2 1 , . . ., j CH1 n are inherited from the parent P 1 and the second part j CH1 r 1 1 , . . ., j CH1 r 2 is taken from the parent P 2. What should be noted is that in order to maintain the precedence feasibility of the new generated children, in the inheritance of the second and the third parts, all elements that already exist in the previous parts of children should be eliminated.The construction of CH2 is the same with CH1.In 21 , the authors proved this crossover operator can create precedence feasible offspring.Figure 1 shows an example of the crossover operation.Suppose a project has 10 nondummy activities, the parents P 1 and P 2 are 1, 2, 3, 4, 5, 6, 7, 8, 9, 10 and 1, 4, 5, 2, 3, 7, 6, 9, 8, 10 , respectively.The crossover points are r 1 4 and r 2 7.According to the two-point crossover operation, the new generated children are CH1 1, 2, 3, 4, 5, 7, 6, 8, 9, 10 and CH1 1, 4, 5, 2, 3, 6, 7, 9, 8, 10 .
Mutation is another important operator which is used to strengthen the algorithm's ability of exploring the unexplored area of the search space.We use the mutation operator presented by 21 .The mutation operator can be described as follows.An random integer  number, a, is generated from interval 1, n .Let activity j b be the last predecessor of activity j a and let activity j c be the first successor of j a in the activity list.Another random integer number, d, is generated from the interval b 1, c − 1 .If d < a, the activity list is replaced by j 0 , . . ., j d−1 , j a , j d , . . ., j a−1 , j a 1 , . . ., j n 1 .If d > a, it is replaced by j 0 , . . ., j a−1 , j a 1 , . . ., j d−1 , j a , j d 1 , . . ., j n 1 .Figure 2 shows an example of mutation operator.

Framework of the H-MOEA
As a nature-inspired metaheuristic, multiobjective evolutionary algorithms MOEAs have been used successfully in the past to solve the multiobjective problems, especially the complex optimization problem.genetic local search, and a comparative experiment study was implemented on the 0/1 knapsack problem 28 .
In order to deal with the stochastic RCPSPs more effectively, a hybrid multiobjective evolutionary algorithm H-MOEA is proposed.In the natural-inspired metaheuristics, the individuals in the population contain a rich information about the solution.In several singleobjective combinatorial optimization problems, experience shows that this heuristic information is useful to the convergence performance of the algorithms 29-32 .Different from singleobjective optimization, the useful information can be obtained from the nondominated set, instead of from a single solution.In the procedure of H-MOEA, the information contained in the nondominated set is extracted periodically and utilized with a Pareto dominance relation based local search 33 .Figure 3 shows the conceptual framework of the proposed H-MOEA.From the figure, it can be seen that H-MOEA consists of two main functional modules: a search module MOEA and a knowledge module.The first module is used to search through the vast solution space and identify the nondominated solutions, while the knowledge module takes charge of obtaining useful information throughout the ongoing optimization process and applying this knowledge to guide the subsequent search process.

Heuristic Information
In this paper, we focus on the extraction and utilization of Position Priority Information which indicates the useful knowledge about the position of an activity in the schedule.The information used to establish an appropriate schedule priority for different activities is stored in a matrix PPI with size n × n, where n is the number of nondummy activities.The element in the matrix PPI i, j denotes the appearance times of an activity N i at the position j.Tables 1  and 2 show an example of the construction of the priority information matrix.Table 1 shows 6 different schedules of a project with 6 nondummy activities.If we look at the first position I 1 , for the 6 different schedules, the scheduled activities at this position are 1, 1, 2, 2, 1, 3, respectively.In other words, for the 6 different activities, the appearance times at the first position are 3, 2, 1, 0, 0 and 0, respectively.Then, the first column of the matrix PPI is accordingly constructed see Table 2 .

Local Search with the Heuristic Information
In H-MOEA, the feasible schedule or part of it can be constructed based on the obtained heuristic information.Then the two-step individual generation procedure 20, 21 can be modified as follows: the first step is identifying the eligible activities, then one of them is selected according to the Position Priority Information.In the second step of the modified individual generation procedure, the heuristic information stored in the matrix is transformed into the selection probability and Roulette Selection is employed to identify the next activity to be scheduled.For example, we suppose the Eligible Activity Set in the second position EA 2 N 3 , N 4 , from the matrix PPI shown in Table 2, and we can see that PPI 3, 3 and PPI 4, 3 are 2 and 1, respectively.Then the selection probability for activities N 3 and N 4 are 2/ 2 1 and 1/ 2 1 , respectively.We denote this Modified Individual Generation procedure as MIG.
In single-objective project scheduling problems, experience showed that good candidate schedules are usually found "fairly close" to other good schedules 34 the local search is conducted with a probability rate, denoted as p ls .The procedure of the local search is illustrated in Algorithm 2. In each generation, the local search is implemented on each individual with the predefined local search probability.For each local search procedure, just one neighbor is generated for the individual.However, the local search procedure employed in the current research is the simplest form.We are aware of that several issues are important in the hybridization of MOEAs and local search, such as balance between genetic search and local search 37, 38 size of examined neighborhood and choice of solutions to which local search is applied 39 .The main focus of this paper is extracting and utilizing the knowledge in RCPSP, and the above-mentioned issues and comparative study with respect to RCPSP will be addressed in our future study to improve the effectiveness of knowledge utilization and the performance of H-MOEA.

Hypothetic Problem Description
In this section, we first studied a hypothetical problem to illustrate the proposed approach.Suppose that a project consists of 20 nondummy activities and there are 4 types of resources with availability R 48, 53, 42, 50 .The precedence relation amongst activities is depicted as in Figure 4.In the network, the nodes N 0 and N 21 are dummy nodes which, respectively, indicate the start node and the end node.The parameters about durations and required resources for each activity are listed in Table 3.

Parameter Settings
We used a population size of 40.Crossover and mutation probabilities were 0.98 and 0.1, respectively.The evolution was terminated after 200 generations.Please note that the crossover probability is defined for each individual.In the calculation of robustness, the sample size was set as 50.For the proposed H-MOEA, the rate of local search was set as 0.45.In all  of the experiments, we assumed that the perturbation interval covered the whole execution process of the schedules, that is, all of the activities were affected by the perturbation.For the hope of eliminating the effect of random generator, each algorithm was repeated 30 times with different random seeds.

Experiment Results
We studied the scenarios under uncertainties with different perturbation strengths 0.15, 0.25, and 0.35.The obtained nondominated solutions in three-dimensional space are depicted as In order to visualize the relation between the makespan and the robustness, we projected the nondominated solutions in the robustness-makespan space see Figure 6 .Figures 6 a , 6 b , and 6 c , respectively, show the projection of the obtained results under perturbation strengths 0.15, 0.25, and 0.35.It can be seen from these figures that some dominated solution were included in the robustness-makespan space.This is because these solution provided a better performance on stability.By looking at the measure of robustness, we can see that, by the increase of the perturbation strength, the robustness value decreased from Figure 6 a to Figure 6 c as expected.
Let us recall that we indicated the robustness of a schedule as the difference between the deterministic makespan before disruption, and the actual makespan after execution.In some literature, such as 10 , a surrogate measure-total amount of free slacks-was employed to indicate the robustness of a schedule.It seems that a schedule with higher measure of free slack sum would cope with the perturbation better.Actually, this situation cannot be held in some cases.For example, by looking at the obtained solutions A and B, it can be seen from Figure 7 that the total amount of the free slacks of schedule A is 4 4 6 14 denoted as blue blocks in the figure , while the value for schedule B is 3 2 6 10 21 see Figure 8 .However, in real execution in the presence of perturbation, solution A, compared to solution  B, had a better performance of absorbing the uncertainty in the environment, thus, had a lower measure of difference with the planned execution in deterministic environment.

Performance Analysis
In this section, we evaluated the performance of the proposed H-MOEA.Here we use MOEA to represent the normal MOEA.We employ the measure of set cover 40 , denoted as SC, to implement the comparison among different approaches.Let A and B be the two obtained nondominated sets, SC A, B is defined as the rate of solutions in B which are dominated by the solutions in A. It is clear that a higher value of SC A, B indicates a better performance of the set A. Note, however, SC A, B is not necessarily equal to SC B, A .

Test Problem Instances
Besides the hypothetical problem described in the above section denoted as P0 , we tested our approach on the benchmark problems generated with Progen, which is a problem generator developed by Kolisch et al. 41 .The network complexity NC and resource factor RF are two important parameters to Progen.The former parameter indicates the average number of immediate successors of an activity.The latter parameter, ranged between 0 and 1, is used to control the percentage of resource types required by an activity.In the generation of   the problem instances, we assumed that the durations, the resource availability, and the required resources for each activity were generated between the predefined intervals, satisfying the resource constraint that the maximum of the required resource cannot exceed the resource availability.Table 4 shows the parameter settings of the 20 generated instances.

Performance Comparison under the Same Number of Generations
The algorithm parameters remained the same as in the above section.We ran the algorithm with all of the generated instances under three different perturbation strengths 0.15, 0.25, and 0.35 .First, we compared the performances of the proposed H-MOEA with normal MOEA under the same number of generations, that is, both algorithms were terminated after 200 generations.The results are reported in the Table 5, where SC 1 represents the set cover SC H-MOEA, MOEA and SC 2 denotes the set cover SC MOEA, H-MOEA .In Table 5, the mean value of the set cover which is better than its reverse figure is presented in bold face.From the figures, one can clearly see that the H-MOEA considerably outperformed the normal one in most cases.For the instances with 20 P0-P6 and 30 P7-P12 nondummy activities, H-MOEA had a significant better performance than MOEA in all scenarios with three different perturbations.For example, if we look at the experimental results of instance P1, the mean values of SC1 under the perturbations with strengths 0.15, 0.25, and 0.35 are 0.2492, 0.3900, and 0.2942, respectively, while the reverse figures SC2 under the perturbations with strengths 0.15, 0.25, and 0.35 are only 0.1508, 0.1908, and 0.1542, respectively.The same set of figures for P6 are 0.4700, 0.5475, 0.6250 and 0.1008, 0.1458, 0.0558, respectively.For the instances with 60 nondummy activities, H-MOEA performed better than normal MOEA in most cases.However, one can see that for the instances P13 under perturbation strength 0.25 and P17 under perturbation strength 0.15, MOEA outperformed H-MOEA, with the mean values of set cover SC 2 0.3242 and 0.3992, respectively, and the reverse figures are 0.2508 and 0.3017.Also, for the instances P14 under perturbation strength 0.35 and P18 under perturbation strength 0.  performance.Finally, we investigated two instances with 120 nondummy activities.The experimental results also suggest that H-MOEA performed better than MOEA under deferent perturbation strength levels.

Remarks on H-MOEA
In the proposed H-MOEA, the information contained in the obtained approximating nondominated solutions were extracted and utilized periodically.The experimental results discussed above show that H-MOEA has a considerable better performance than normal MOEA.This indicates that the heuristic information accumulated in the previous runs is very useful to guide the subsequent search process.Since we employed a local search procedure in the process of MOEA, it was unavoidable that the computation time of the proposed H-MOEA increased to a certain degree.In other words, the improvement of convergence performance was with the sacrifice of computation time.where CT H-MOEA denotes the computation time of the proposed H-MOEA and CT MOEA represents the computation time of the normal MOEA.We tested the effect of local search rate p ls on the proposed H-MOEA by varying the local search rate from 0.1 to 0.75 on the problems P0 20 nondummy activities , P8 30 nondummy activities , and P15 60 nondummy activities under perturbation strength 0.15. Figure 9 shows the values of relative computation time RCT with variation of local search rate p ls .The values of relative set cover RSC with variation of p ls are reported in Figure 10.All the values in the figures are the mean values of the obtained results after 30 repeated runs.It can be clearly seen, from Figure 9, that, by the increase of local search rate, more computation time was needed to perform the H-MOEA.
The definition of relative set cover RSC in 4.1 shows that a value greater than 1 indicates a better performance for H-MOEA.By looking at Figure 10, it can be seen that even a lower local search rate 0.1 could yield a better performance for the instances with 20 nondummy activities P0 and 30 nondummy activities P8 , with RSC values about 1.30 and 1.13, respectively.It is interesting to notice that the values of relative set cover RSC did not improve linearly with the increase of local search rate see Figure 10 .In some cases, for example, for instance with 20 nondummy activities P0 , the highest value of RSC did not correspond to the highest local search rate.One of the reasons for this is that a higher local search rate might limit the exploration ability of the algorithm in the search space.Thus, with another consideration, computation time, the local search rate should not be too high.In this research, we commonly consider that a local search rate between 0.3 and 0.5 is suitable to the proposed algorithm and the problem domain.In real application on other complex problems, a tradeoff between convergence performance and computational time is possibly needed.

Performance Comparison under the Same Computational Effort
As indicated above, due to local search at each generation, much more solutions are examined in proposed H-MOEA using more computation time.In order to further investigate the performance of proposed H-MOEA, we, respectively, compared the performances of H-MOEA and MOEA under the same number of examined solutions and computational time.
To do this, we executed the algorithms as follows: at first, we ran the MOEA with 200 generations and recorded the number of examined solutions in MOEA and its computational time.Then, the proposed H-MOEA was executed with the following two termination conditions: 1 the number of examined solutions was achieved; 2 the computational time was achieved.

Conclusion
In this paper, we study a resource-constrained project scheduling problem in the presence of perturbation on activity durations.Robust scheduling is employed as the methodology to solve this problem.We measure the robustness of a schedule as the expectation of the difference between the actual makespan after execution and the planned makespan.The stability of the schedule under the perturbation is indicated by the statistical measure of standard deviation.Then, the problem is modelled as a multiobjective optimization problem, where makespan minimization, robustness maximization, and stability maximization standard deviation minimization are considered simultaneously.
We employ multiobjective evolutionary algorithm MOEA to obtain the Pareto optimal solutions of the problem.The normal MOEA is improved by incorporating a knowledgebased local search procedure.We call the improved MOEA as Hybrid MOEA H-MOEA .In the process of the proposed H-MOEA, two additional procedures are integrated: information extraction and information utilization.In the former procedure, heuristic information is obtained periodically from the obtained approximating nondominated individuals.In the latter procedure, the obtained information is utilized in the local search to improve the individuals in the population.Experimental study shows that the proposed approach is feasible and effective for the resource-constrained project scheduling problem with stochastic durations and that the proposed H-MOEA performs better than normal MOEA.In future work, we plan to employ the hybrid multiobjective evolutionary algorithm to investigate the project scheduling problem under dynamic environment.It is also worthwhile to focus on the design of more efficient hybrid MOEAs for solving large-scale problems in future work.
d i : the duration of activity N i in the deterministic situation.
st i : the start time of activity N i in the deterministic situation.
ft i : the finish time of activity N i in the deterministic situation.
d ui : the duration of activity N i in the presence of perturbation.

Figure 1 :
Figure 1: An example of crossover operation.

Figure 2 :
Figure 2: An example of mutation operation.

Figure 3 :
Figure 3: The conceptual framework of the hybrid multiobjective evolutionary algorithm for RCPSP.

Figure 4 :
Figure 4: The precedence relation of tasks.

Figure 5 :
Figure 5: The obtained nondominated solutions in the three-dimensional space with different perturbation strengths: a perturbation strength 0.15; b perturbation strength 0.25; c perturbation strength 0.35.

Figure 6 :
Figure 6: The projection of the nondominated solutions in the robustness-makespan space with different perturbation strengths: a perturbation strength 0.15; b perturbation strength 0.25; c perturbation strength 0.35.

Figure 9 :
Figure 9: The obtained relative computation time RCT with different local search rates for P0, P8, and P15 with a perturbation strength 0.15.

Figure 10 :
Figure 10: The obtained relative set cover RSC with different local search rates for P0, P8, and P15 with a perturbation strength 0.15.
To date, t26re are several representative MOEAs.For the brief view history of the MOEAs and the mechanism of each MOEA, readers are referred to the survey 24 .In order to speed up the convergency of the optimization algorithm, local search procedure is usually hybridized with the MOEAs.Ishibuchi andMurata 25,26were among the first to implement such a hybridization.27 improved the performance of multiobjective

Table 2 :
Example of the Position Priority Information matrix.
-36 .Encouraged by this experience, we employ a local search procedure around the neighborhood of the schedule in the multiobjective version of the PSPs.Let Ind i be an individual in the population, Ind new i the individual after the local search.If Ind new i dominates Ind i , then replace the Ind i with Ind new i in the population.Similar with the crossover and the mutation operators, i 1, p 0, m 0 while i ≤ population size do p randomly generated value between the interval 0, 1 if p ≤ p ls then m randomly generated value between the interval 1,n .Copy the 1, m part of the individual Ind i to the individual Ind new Pseudocode of the local search procedure.

Table 3 :
Durations and required resources of the nondummy activities of the hypothetic problem.

Table 4 :
Parameter settings of the 20 generated instances of benchmark.

Table 5 :
The mean value and the variance of the set covers SC H-MOEA, MOEA SC 1 and SC MOEA, H-MOEA SC 2 within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same number of generations.
25, H-MOEA and MOEA approximately had the same Thus, it is an important issue to investigate the relation between the local search rate and the performance of H-MOEA.Here, we employ two measures, relative set cover RSC and relative computation time RCT , to indicate the performance of the proposed H-MOEA.RSC and RCT are calculated as follows: