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

Authors: Jian Xiong; Ying-wu Chen; Ke-wei Yang; Qing-song Zhao; Li-ning Xing
Journal: Mathematical Problems in Engineering
 (2012)

---

## Abstract

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.


---

## Body

## 1. 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 asrobust 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 toquality 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.


## 2. Problem Description

### 2.1. 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 activitiesN=N0,N1,…,Nn,Nn+1, where each activity has to be executed without interruption. N0 and Nn+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=(R1,R2,…,RK), where Rk(k=1,…,K) represents the availability of each resource type k in each time period. For each activity Ni(i=1,…,n), the duration is denoted as di, and the requested resources vector over the entire duration is represented as ri=(ri1,…,riK). Note that for dummy activities N0 and Nn+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 Ni are denoted as Succi and Predi, 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 sti and fti, respectively. A schedule S consists of the vector of start time of the nondummy activities, denoted as S=(st1,…,stn). The makespan of the schedule S is denoted as M(S).


### 2.2. 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 stiu and ftiu, respectively. The actual duration is denoted as diu and is calculated as follows:(2.1)diu={Round((1+St)di)+1,di,ifstiu∈[OT,ET]orstiu+di∈[OT,ET],otherwise,
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:(2.2)Robμ(S)=1exp⁡(Mμ(S)-M(S)/M(S)).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:(2.3)Rob(S)=E[Robμ(S)]=1Nμ∑m=1NμRobμm(S),
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 asStd[Robμ(S)] and calculated as follows:(2.4)Std[Robμ(S)]=1Nμ∑m=1Nμ(Robμm(S)-E[Robμ(S)])2.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:(2.5)obj.:(1)min⁡f1=M(S),(2)max⁡f2=E[Robμ(S)],(3)min⁡f3=Std[Robμ(S)]s.t.:(1)max⁡{stj+dj∣Nj∈Predi}≤sti,∀i,(2)max⁡{stju+dju∣Nj∈Predi}≤stiu,∀i,(3)∑i=1nrik(t)≤Rk,∀tandk=1,2,…,K.
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 Ni 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 individualI1 dominates individual I2 if I1 does at least as well as I2 on all objectives, and strictly better on at least one. Thus, in this paper, we formally define the schedule dominance as follows.Definition 2.1.
A scheduleS1dominates scheduleS2ifM(S1)≤M(S2)∧E[Robμ(S1)≥E[Robμ(S2)]∧Std[Robμ(S1)]≤Std[Robμ(S2)], with strict inequality holding for at least one objective measure.The nondominance concept can be defined as follows.Definition 2.2.
A scheduleS*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.Rk:
the availability ofkth type of resource, k=1,…,K.rik:
the amount of thekth type of resource requested by activity Ni, i=1,…,n.di:
the duration of activityNi in the deterministic situation.sti:
the start time of activityNi in the deterministic situation.fti:
the finish time of activityNi in the deterministic situation.diu:
the duration of activityNi in the presence of perturbation.stiu:
the start time of activityNi in the presence of perturbation.ftiu:
the finish time of activityNi in the presence of perturbation.S:
the schedule,S=(st1,…,stn).M(S):
the makespan of scheduleS in the deterministic situation.Mμ(S):
the makespan of scheduleS under a specific realization of uncertainty μ.St:
the perturbation strength.OT:
the occurrence time of the perturbation.ET:
the end time of the perturbation.


## 2.1. 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 activitiesN=N0,N1,…,Nn,Nn+1, where each activity has to be executed without interruption. N0 and Nn+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=(R1,R2,…,RK), where Rk(k=1,…,K) represents the availability of each resource type k in each time period. For each activity Ni(i=1,…,n), the duration is denoted as di, and the requested resources vector over the entire duration is represented as ri=(ri1,…,riK). Note that for dummy activities N0 and Nn+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 Ni are denoted as Succi and Predi, 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 sti and fti, respectively. A schedule S consists of the vector of start time of the nondummy activities, denoted as S=(st1,…,stn). The makespan of the schedule S is denoted as M(S).


## 2.2. 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 stiu and ftiu, respectively. The actual duration is denoted as diu and is calculated as follows:(2.1)diu={Round((1+St)di)+1,di,ifstiu∈[OT,ET]orstiu+di∈[OT,ET],otherwise,
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:(2.2)Robμ(S)=1exp⁡(Mμ(S)-M(S)/M(S)).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:(2.3)Rob(S)=E[Robμ(S)]=1Nμ∑m=1NμRobμm(S),
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 asStd[Robμ(S)] and calculated as follows:(2.4)Std[Robμ(S)]=1Nμ∑m=1Nμ(Robμm(S)-E[Robμ(S)])2.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:(2.5)obj.:(1)min⁡f1=M(S),(2)max⁡f2=E[Robμ(S)],(3)min⁡f3=Std[Robμ(S)]s.t.:(1)max⁡{stj+dj∣Nj∈Predi}≤sti,∀i,(2)max⁡{stju+dju∣Nj∈Predi}≤stiu,∀i,(3)∑i=1nrik(t)≤Rk,∀tandk=1,2,…,K.
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 Ni 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 individualI1 dominates individual I2 if I1 does at least as well as I2 on all objectives, and strictly better on at least one. Thus, in this paper, we formally define the schedule dominance as follows.Definition 2.1.
A scheduleS1dominates scheduleS2ifM(S1)≤M(S2)∧E[Robμ(S1)≥E[Robμ(S2)]∧Std[Robμ(S1)]≤Std[Robμ(S2)], with strict inequality holding for at least one objective measure.The nondominance concept can be defined as follows.Definition 2.2.
A scheduleS*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.Rk:
the availability ofkth type of resource, k=1,…,K.rik:
the amount of thekth type of resource requested by activity Ni, i=1,…,n.di:
the duration of activityNi in the deterministic situation.sti:
the start time of activityNi in the deterministic situation.fti:
the finish time of activityNi in the deterministic situation.diu:
the duration of activityNi in the presence of perturbation.stiu:
the start time of activityNi in the presence of perturbation.ftiu:
the finish time of activityNi in the presence of perturbation.S:
the schedule,S=(st1,…,stn).M(S):
the makespan of scheduleS in the deterministic situation.Mμ(S):
the makespan of scheduleS under a specific realization of uncertainty μ.St:
the perturbation strength.OT:
the occurrence time of the perturbation.ET:
the end time of the perturbation.


## 3. 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.


### 3.1. Chromosome Representation and Population Initialization

#### 3.1.1. 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:(3.1)I=(j0,j1,…,jn,jn+1).
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 ofn 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.Algorithm 1:Pseudocode of the Serial Schedule Generation Scheme (SSGS).
i=1whilei≤ndosti=0,sti′=0,ftmax=0Obtain the maximum finish time of the predecessors:ftmax=max{ftj∣Nj∈Predi}Determine the earliest time point,sti′, which satisfying the resource con-straints andsti′≥ftmax.Setsti=sti′.end while


#### 3.1.2. Population Initialization

In the operator of population initialization, the individuals are generated serially via a two-step 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 be eligible if all of its predecessors have been scheduled. Formally, let P(i) be the position of activity Ni in the activity list, if we look at the position q(q=0,1,…,n+1), an activity Nq* is said to be an eligible activity at the position q if P(j)<q, for all j,Nj∈Predq*. The set which consists of all eligible activities at position q is denoted as EAq.


### 3.2. 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, r1 and r2, r1<r2, 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 P1 and P2 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 (r1 and r2). For example, CH1 is constructed as (j1CH1,…,jr1CH1,jr1+1CH1,…,jr2CH1,jr2+1CH1,…,jnCH1), where the first part (j1CH1,…,jr1CH1) and the third part (jr2+1CH1,…,jnCH1) are inherited from the parent P1 and the second part (jr1+1CH1,…,jr2CH1) is taken from the parent P2. 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 P1 and P2 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 r1=4 and r2=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).Figure 1
An example of crossover operation.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 jb′ be the last predecessor of activity ja′ and let activity jc′ be the first successor of ja′ 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 (j0′,…,jd-1′,ja′,jd′,…,ja-1′,ja+1′,…,jn+1′). If d>a, it is replaced by (j0′,…,ja-1′,ja+1′,…,jd-1′,ja′,jd+1′,…,jn+1′). Figure 2 shows an example of mutation operator.Figure 2
An example of mutation operation.


### 3.3. Hybrid Multiobjective Evolutionary Algorithm for RCPSP

#### 3.3.1. 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. To date, there 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 and Murata [25, 26] were among the first to implement such a hybridization. [27] improved the performance of multiobjective 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 single-objective combinatorial optimization problems, experience shows that this heuristic information is useful to the convergence performance of the algorithms [29–32]. Different from single-objective 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.Figure 3
The conceptual framework of the hybrid multiobjective evolutionary algorithm for RCPSP.


#### 3.3.2. Heuristic Information

In this paper, we focus on the extraction and utilization ofPosition 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 Ni 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 I1, 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).Table 1
Example of 6 schedules with 6 nondummy activities.
I1I2I3I4I5I6Schedule 1123456Schedule 2132465Schedule 3213456Schedule 4231546Schedule 5124365Schedule 6312645Table 2
Example of thePosition Priority Information matrix.
I1I2I3I4I5I6Activity 1321000Activity 2222000Activity 3122100Activity 4001320Activity 5000123Activity 6000123


#### 3.3.3. 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 EA2=(N3,N4), 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 N3 and N4 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–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 Indi be an individual in the population, Indinew the individual after the local search. If Indinew dominates Indi, then replace the Indi with Indinew in the population. Similar with the crossover and the mutation operators, the local search is conducted with a probability rate, denoted as pls. 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.Algorithm 2:Pseudocode of the local search procedure.
i=1,p=0,m=0whilei≤population size  dop= randomly generated value between the interval [0,1]ifp≤pls  thenm= randomly generated value between the interval [1,n].Copy the[1,m] part of the individual Indi to the individual Indinew.Generate the[m+1,n] part of the individual Indinew by the operatorMIG.Evaluate the individualIndinew.if  Indinew dominates Indi thenReplaceIndi with Indinew.end ifend ifi=i+1end while


## 3.1. Chromosome Representation and Population Initialization

### 3.1.1. 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:(3.1)I=(j0,j1,…,jn,jn+1).
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 ofn 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.Algorithm 1:Pseudocode of the Serial Schedule Generation Scheme (SSGS).
i=1whilei≤ndosti=0,sti′=0,ftmax=0Obtain the maximum finish time of the predecessors:ftmax=max{ftj∣Nj∈Predi}Determine the earliest time point,sti′, which satisfying the resource con-straints andsti′≥ftmax.Setsti=sti′.end while


### 3.1.2. Population Initialization

In the operator of population initialization, the individuals are generated serially via a two-step 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 be eligible if all of its predecessors have been scheduled. Formally, let P(i) be the position of activity Ni in the activity list, if we look at the position q(q=0,1,…,n+1), an activity Nq* is said to be an eligible activity at the position q if P(j)<q, for all j,Nj∈Predq*. The set which consists of all eligible activities at position q is denoted as EAq.


## 3.1.1. 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:(3.1)I=(j0,j1,…,jn,jn+1).
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 ofn 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.Algorithm 1:Pseudocode of the Serial Schedule Generation Scheme (SSGS).
i=1whilei≤ndosti=0,sti′=0,ftmax=0Obtain the maximum finish time of the predecessors:ftmax=max{ftj∣Nj∈Predi}Determine the earliest time point,sti′, which satisfying the resource con-straints andsti′≥ftmax.Setsti=sti′.end while


## 3.1.2. Population Initialization

In the operator of population initialization, the individuals are generated serially via a two-step 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 be eligible if all of its predecessors have been scheduled. Formally, let P(i) be the position of activity Ni in the activity list, if we look at the position q(q=0,1,…,n+1), an activity Nq* is said to be an eligible activity at the position q if P(j)<q, for all j,Nj∈Predq*. The set which consists of all eligible activities at position q is denoted as EAq.


## 3.2. 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, r1 and r2, r1<r2, 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 P1 and P2 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 (r1 and r2). For example, CH1 is constructed as (j1CH1,…,jr1CH1,jr1+1CH1,…,jr2CH1,jr2+1CH1,…,jnCH1), where the first part (j1CH1,…,jr1CH1) and the third part (jr2+1CH1,…,jnCH1) are inherited from the parent P1 and the second part (jr1+1CH1,…,jr2CH1) is taken from the parent P2. 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 P1 and P2 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 r1=4 and r2=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).Figure 1
An example of crossover operation.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 jb′ be the last predecessor of activity ja′ and let activity jc′ be the first successor of ja′ 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 (j0′,…,jd-1′,ja′,jd′,…,ja-1′,ja+1′,…,jn+1′). If d>a, it is replaced by (j0′,…,ja-1′,ja+1′,…,jd-1′,ja′,jd+1′,…,jn+1′). Figure 2 shows an example of mutation operator.Figure 2
An example of mutation operation.


## 3.3. Hybrid Multiobjective Evolutionary Algorithm for RCPSP

### 3.3.1. 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. To date, there 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 and Murata [25, 26] were among the first to implement such a hybridization. [27] improved the performance of multiobjective 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 single-objective combinatorial optimization problems, experience shows that this heuristic information is useful to the convergence performance of the algorithms [29–32]. Different from single-objective 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.Figure 3
The conceptual framework of the hybrid multiobjective evolutionary algorithm for RCPSP.


### 3.3.2. Heuristic Information

In this paper, we focus on the extraction and utilization ofPosition 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 Ni 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 I1, 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).Table 1
Example of 6 schedules with 6 nondummy activities.
I1I2I3I4I5I6Schedule 1123456Schedule 2132465Schedule 3213456Schedule 4231546Schedule 5124365Schedule 6312645Table 2
Example of thePosition Priority Information matrix.
I1I2I3I4I5I6Activity 1321000Activity 2222000Activity 3122100Activity 4001320Activity 5000123Activity 6000123


### 3.3.3. 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 EA2=(N3,N4), 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 N3 and N4 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–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 Indi be an individual in the population, Indinew the individual after the local search. If Indinew dominates Indi, then replace the Indi with Indinew in the population. Similar with the crossover and the mutation operators, the local search is conducted with a probability rate, denoted as pls. 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.Algorithm 2:Pseudocode of the local search procedure.
i=1,p=0,m=0whilei≤population size  dop= randomly generated value between the interval [0,1]ifp≤pls  thenm= randomly generated value between the interval [1,n].Copy the[1,m] part of the individual Indi to the individual Indinew.Generate the[m+1,n] part of the individual Indinew by the operatorMIG.Evaluate the individualIndinew.if  Indinew dominates Indi thenReplaceIndi with Indinew.end ifend ifi=i+1end while


## 3.3.1. 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. To date, there 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 and Murata [25, 26] were among the first to implement such a hybridization. [27] improved the performance of multiobjective 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 single-objective combinatorial optimization problems, experience shows that this heuristic information is useful to the convergence performance of the algorithms [29–32]. Different from single-objective 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.Figure 3
The conceptual framework of the hybrid multiobjective evolutionary algorithm for RCPSP.


## 3.3.2. Heuristic Information

In this paper, we focus on the extraction and utilization ofPosition 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 Ni 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 I1, 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).Table 1
Example of 6 schedules with 6 nondummy activities.
I1I2I3I4I5I6Schedule 1123456Schedule 2132465Schedule 3213456Schedule 4231546Schedule 5124365Schedule 6312645Table 2
Example of thePosition Priority Information matrix.
I1I2I3I4I5I6Activity 1321000Activity 2222000Activity 3122100Activity 4001320Activity 5000123Activity 6000123


## 3.3.3. 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 EA2=(N3,N4), 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 N3 and N4 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–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 Indi be an individual in the population, Indinew the individual after the local search. If Indinew dominates Indi, then replace the Indi with Indinew in the population. Similar with the crossover and the mutation operators, the local search is conducted with a probability rate, denoted as pls. 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.Algorithm 2:Pseudocode of the local search procedure.
i=1,p=0,m=0whilei≤population size  dop= randomly generated value between the interval [0,1]ifp≤pls  thenm= randomly generated value between the interval [1,n].Copy the[1,m] part of the individual Indi to the individual Indinew.Generate the[m+1,n] part of the individual Indinew by the operatorMIG.Evaluate the individualIndinew.if  Indinew dominates Indi thenReplaceIndi with Indinew.end ifend ifi=i+1end while


## 4. Case Study

### 4.1. 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 availabilityR=(48,53,42,50). The precedence relation amongst activities is depicted as in Figure 4. In the network, the nodes N0 and N21 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.Table 3
Durations and required resources of the nondummy activities of the hypothetic problem.
NidiriN117(19, 20, 10, 19)N230(13, 7, 9, 13)N323(14, 7, 7, 20)N421(19, 20, 20, 20)N529(13, 15, 5, 20)N626(14, 19, 6, 16)N729(18, 12, 18, 6)N822(11, 18, 10, 5)N921(12, 16, 7, 18)N1024(6, 7, 19, 16)N1126(13, 9, 6, 15)N1222(18, 8, 9, 11)N1325(8, 13, 16, 12)N1419(16, 13, 6, 20)N1524(10, 15, 6, 18)N1625(17, 20, 19, 7)N1720(14, 22, 13, 12)N1818(13, 17, 20, 10)N1927(21, 22, 19, 17)N2020(16, 24, 17, 14)Figure 4
The precedence relation of tasks.


### 4.2. 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.


### 4.3. 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 Figure5, where the distribution of the solutions under perturbation strengths 0.15, 0.25, and 0.35 is shown in Figures 5(a), 5(b), and 5(c), respectively. The well-distributed solutions shown in these figures indicate that the proposed three objectives, named makespan, robustness, and stability, are in conflict with each other.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.
(a)(b)(c)In order to visualize the relation between the makespan and the robustness, we projected the nondominated solutions in the robustness-makespan space (see Figure6). 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.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.
(a)(b)(c)Another important issue of RCPSP is genotype-phenotype analysis. We investigated two nondominated solutions obtained under the perturbation 0.25. SolutionA, IA=(0,1,6,3,13,2,7,4,16,5,9,10,12,19,8,14,15,20,17,11,18,21), with objective values makespan= 210, robustness= 0.8951, and standard deviation= 0.0859; solution B, IB=(0,3,1,6,4,16,13,2,5,7,9,12,8,14,11,10,19,17,15,18,20,21), with objective values makespan = 210, robustness = 0.8444, and standard deviation = 0.0342. One can see that although these two solutions had the same makespans, solution A had a better measure on robustness than solution B. Figures 7 and 8, respectively, show the genotype-phenotype mapping of solutions A and B. It can be seen from the figures that for both schedules, since there were many available resources, more than one precedence-feasible activities were allowed to be executed concurrently.Figure 7
Snapshots of obtained scheduleA, with objective values makespan= 210, robustness= 0.8951, and standard deviation= 0.0859.Figure 8
Snapshots of obtained scheduleA, with objective values makespan = 210, robustness= 0.8444, and standard deviation= 0.0342.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.


### 4.4. 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 ofset 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).


#### 4.4.1. Test Problem Instances

Besides the hypothetical problem described in the above section (denoted as P0), we tested our approach on the benchmark problems generated withProgen, 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.Table 4
Parameter settings of the 20 generated instances of benchmark.
IndexNNCKRFIndexNNCKRFP1201.840.75P11302.520.75P2201.841P12302.521P3202.140.75P13601.840.75P4202.141P14601.841P5202.520.75P15602.140.75P6202.521P16602.141P7301.840.75P17602.520.75P8301.841P18602.521P9302.140.75P191201.840.75P10302.141P201201.841


#### 4.4.2. 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 proposedH-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 SC1 represents the set coverSC(H-MOEA, MOEA) and SC2 denotes the set coverSC(MOEA, H-MOEA).Table 5
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same number of generations.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3158±0.06310.1158±0.01510.4375±0.06750.0983±0.00930.4550±0.09670.0883±0.0076P10.2492±0.04870.1508±0.01760.3900±0.09310.1908±0.03670.2942±0.06590.1542±0.0270P20.2933±0.04010.1792±0.02450.3167±0.06470.1683±0.03270.3483±0.05640.2058±0.0427P30.4383±0.04870.2000±0.02000.4667±0.09010.0825±0.01240.2992±0.03540.1200±0.0150P40.2942±0.05300.1350±0.02450.3400±0.06000.1542±0.02210.4517±0.09740.1400±0.0264P50.2892±0.05740.1250±0.01450.4283±0.08170.1608±0.02740.4267±0.08200.1758±0.0416P60.4700±0.04060.1008±0.01450.5475±0.07590.1458±0.02560.6250±0.08700.0558±0.0068P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3517±0.05520.2583±0.06530.3542±0.06450.2800±0.05340.3608±0.06940.2667±0.0302P90.3858±0.07470.2317±0.03540.2600±0.02900.2592±0.02590.3864±0.07430.2218±0.0189P100.4625±0.06530.1442±0.03960.4458±0.07770.1733±0.03450.4408±0.09790.1442±0.0299P110.5025±0.09810.1650±0.04320.3608±0.06290.1500±0.02810.4127±0.05030.2967±0.0363P120.3258±0.04940.2258±0.04460.3025±0.05830.2650±0.04550.3325±0.06470.3033±0.0510P130.3033±0.04270.2367±0.04150.2508±0.04410.3242±0.07800.3783±0.06730.2433±0.0610P140.3408±0.06360.2567±0.04560.3483±0.07160.2992±0.05590.2908±0.05620.2974±0.0758P150.3983±0.08030.2492±0.04810.3358±0.08140.2325±0.06630.3175±0.06120.2608±0.0549P160.4392±0.07810.2383±0.05350.3292±0.06270.2708±0.03980.3667±0.09270.3017±0.0765P170.3017±0.06930.3992±0.07270.3500±0.06800.3225±0.07270.4075±0.09740.2642±0.0780P180.3400±0.04750.2833±0.04470.3508±0.08340.3433±0.07000.3792±0.05980.2683±0.0440P190.3842±0.05280.2817±0.07200.3983±0.09120.2375±0.06940.3375±0.08040.2533±0.0740P200.3858±0.06490.1375±0.02710.3508±0.10010.2392±0.07200.3075±0.05160.2358±0.0592In Table5, 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 coverSC2 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.25, H-MOEA and MOEA approximately had the same 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.


#### 4.4.3. 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. 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 coverRSC and relative computation time RCT, to indicate the performance of the proposed H-MOEA. RSC and RCT are calculated as follows:(4.1)RSC=SC(H-MOEA,MOEA)SC(MOEA,H-MOEA),(4.2)RCT=CTH-MOEACTMOEA,
where CTH-MOEA denotes the computation time of the proposed H-MOEA and CTMOEA represents the computation time of the normal MOEA.We tested the effect of local search ratepls 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 pls. The values of relative set cover RSC with variation of pls 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.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
The obtained relative set cover (RSC) with different local search rates for P0, P8, and P15 with a perturbation strength 0.15.The definition of relative set coverRSC 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.


#### 4.4.4. 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.Table6 shows the computational results obtained under the same examined solutions, while the results under the same computation time are reported in Table 7. After careful examination of the data, it could be found that the proposed H-MOEA outperformed normal MOEA on small size instances, while the performance of H-MOEA decreased on relative large scale cases. By looking at the instances P0–P12, the data suggests H-MOEA performed better than normal MOEA in most instances, under both the same examined solutions (see Table 6) and the same computation time (see Table 7). For relative large scale instances (P13–P20), shown as in Table 6, H-MOEA still outperformed MOEA in most cases under the same examined solutions. While under the same computation time, the performances of H-MOEA and MOEA were comparable on instances P13–P20, suggested as in Table 7. The main reason of performance decrease of H-MOEA on large scale instances could be considered as limitation on global search ability of the algorithm. As indicated in [37], the number of generations is decreased when the available computation time is limited. As a result, the exploration in search space is not fully utilized.Table 6
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same examined solutions.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3692±0.09390.0983±0.00830.4425±0.09000.1108±0.01470.4117±0.07560.1017±0.0136P10.3275±0.05590.1267±0.02390.3992±0.09660.1425±0.02420.2800±0.05360.1842±0.0313P20.2792±0.03800.1458±0.01700.3025±0.04910.1633±0.02190.3233±0.04760.2242±0.0446P30.3392±0.04140.2383±0.02870.4067±0.08830.1392±0.02530.3392±0.07260.1175±0.0207P40.2833±0.04690.1883±0.04230.2775±0.03350.1925±0.03170.4492±0.10740.2408±0.0659P50.3100±0.06590.1492±0.02050.4158±0.07480.1908±0.03070.4150±0.05420.1583±0.0273P60.4433±0.05820.1292±0.02930.5600±0.12260.2292±0.06350.4800±0.07320.1175±0.0188P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3075±0.05260.2942±0.04870.3208±0.07770.3333±0.07740.3075±0.05790.3117±0.0492P90.2808±0.04280.2617±0.03220.3433±0.05760.2408±0.03350.3627±0.03480.2769±0.0198P100.3117±0.07410.2233±0.02990.3758±0.05010.1775±0.02280.3483±0.05850.2092±0.0564P110.5037±0.10870.1650±0.04620.3642±0.05240.1558±0.02970.3128±0.02650.2364±0.0127P120.2883±0.05970.2983±0.06090.3167±0.06190.2700±0.04470.3542±0.05400.3500±0.0780P130.2458±0.05230.3183±0.07460.2475±0.03480.3458±0.08450.3733±0.08950.2883±0.0681P140.3642±0.07910.2525±0.05340.3175±0.08260.3242±0.06320.3726±0.06350.3214±0.0468P150.3558±0.09460.3033±0.06190.3162±0.04860.3210±0.03310.4058±0.05960.3374±0.0288P160.3225±0.05030.2783±0.04550.3433±0.08800.2792±0.06290.2700±0.04800.3750±0.0840P170.2717±0.05740.3667±0.06610.3292±0.05750.3033±0.04880.2833±0.06540.3767±0.0676P180.3483±0.07070.2417±0.05220.3150±0.06740.3925±0.06870.3842±0.06940.2258±0.0473P190.3458±0.07270.3575±0.07900.3650±0.07790.3142±0.05720.3150±0.06090.3125±0.0863P200.3254±0.06490.2384±0.03410.2704±0.08010.2137±0.05200.2674±0.06130.2059±0.0462Table 7
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same computational time.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.2800±0.05510.1367±0.00950.3800±0.04280.1083±0.00820.4775±0.07330.0892±0.0117P10.3025±0.06840.1633±0.02880.3858±0.08450.0917±0.01150.2250±0.04070.2517±0.0365P20.2075±0.01280.2200±0.01630.2850±0.04620.1833±0.02130.3233±0.04240.2217±0.0280P30.4133±0.07860.2025±0.02850.3850±0.06820.1067±0.01300.3108±0.04270.1233±0.0069P40.2617±0.05010.1917±0.02900.3183±0.04710.1842±0.03050.3575±0.06160.2092±0.0606P50.2125±0.01610.1525±0.01860.3650±0.06910.1942±0.04060.3833±0.08670.1975±0.0425P60.4475±0.07460.1283±0.03620.5933±0.07460.1475±0.03230.4800±0.07440.1417±0.0282P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.2983±0.06450.2808±0.05990.4275±0.06940.2208±0.04250.3133±0.04450.3250±0.0488P90.2883±0.04860.2967±0.04680.2833±0.04150.2308±0.01790.3015±0.03120.2453±0.0203P100.3524±0.04520.29432±0.02980.3453±0.06740.30343±0.04460.4018±0.09060.3472±0.0439P110.4725±0.10440.2175±0.06180.4450±0.08400.1717±0.03080.4531±0.03680.2741±0.0153P120.3008±0.04890.2617±0.05550.3000±0.03720.3025±0.05750.3383±0.06740.3433±0.0684P130.2158±0.04720.4250±0.05940.2267±0.04000.3825±0.06510.3117±0.05790.3158±0.0816P140.3967±0.07420.2583±0.04660.3458±0.06990.3513±0.09920.3369±0.02640.2867±0.0224P150.3350±0.04880.2750±0.04590.3452±0.04310.3501±0.05230.3940±0.05110.4012±0.0634P160.3417±0.08650.3400±0.09710.3233±0.05250.2692±0.06230.3267±0.08580.3350±0.0696P170.2875±0.06560.2904±0.05570.2575±0.06500.2467±0.06550.3183±0.05880.3258±0.0854P180.3342±0.06100.3450±0.07750.2575±0.06380.2592±0.06080.2908±0.05740.3000±0.0497P190.2883±0.04510.3267±0.07450.2683±0.07120.3458±0.07860.2567±0.05570.4400±0.0887P200.2978±0.04470.3004±0.06310.2503±0.06010.2294±0.05230.3772±0.07190.3685±0.0793


## 4.1. 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 availabilityR=(48,53,42,50). The precedence relation amongst activities is depicted as in Figure 4. In the network, the nodes N0 and N21 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.Table 3
Durations and required resources of the nondummy activities of the hypothetic problem.
NidiriN117(19, 20, 10, 19)N230(13, 7, 9, 13)N323(14, 7, 7, 20)N421(19, 20, 20, 20)N529(13, 15, 5, 20)N626(14, 19, 6, 16)N729(18, 12, 18, 6)N822(11, 18, 10, 5)N921(12, 16, 7, 18)N1024(6, 7, 19, 16)N1126(13, 9, 6, 15)N1222(18, 8, 9, 11)N1325(8, 13, 16, 12)N1419(16, 13, 6, 20)N1524(10, 15, 6, 18)N1625(17, 20, 19, 7)N1720(14, 22, 13, 12)N1818(13, 17, 20, 10)N1927(21, 22, 19, 17)N2020(16, 24, 17, 14)Figure 4
The precedence relation of tasks.


## 4.2. 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.


## 4.3. 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 Figure5, where the distribution of the solutions under perturbation strengths 0.15, 0.25, and 0.35 is shown in Figures 5(a), 5(b), and 5(c), respectively. The well-distributed solutions shown in these figures indicate that the proposed three objectives, named makespan, robustness, and stability, are in conflict with each other.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.
(a)(b)(c)In order to visualize the relation between the makespan and the robustness, we projected the nondominated solutions in the robustness-makespan space (see Figure6). 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.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.
(a)(b)(c)Another important issue of RCPSP is genotype-phenotype analysis. We investigated two nondominated solutions obtained under the perturbation 0.25. SolutionA, IA=(0,1,6,3,13,2,7,4,16,5,9,10,12,19,8,14,15,20,17,11,18,21), with objective values makespan= 210, robustness= 0.8951, and standard deviation= 0.0859; solution B, IB=(0,3,1,6,4,16,13,2,5,7,9,12,8,14,11,10,19,17,15,18,20,21), with objective values makespan = 210, robustness = 0.8444, and standard deviation = 0.0342. One can see that although these two solutions had the same makespans, solution A had a better measure on robustness than solution B. Figures 7 and 8, respectively, show the genotype-phenotype mapping of solutions A and B. It can be seen from the figures that for both schedules, since there were many available resources, more than one precedence-feasible activities were allowed to be executed concurrently.Figure 7
Snapshots of obtained scheduleA, with objective values makespan= 210, robustness= 0.8951, and standard deviation= 0.0859.Figure 8
Snapshots of obtained scheduleA, with objective values makespan = 210, robustness= 0.8444, and standard deviation= 0.0342.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.


## 4.4. 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 ofset 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).


### 4.4.1. Test Problem Instances

Besides the hypothetical problem described in the above section (denoted as P0), we tested our approach on the benchmark problems generated withProgen, 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.Table 4
Parameter settings of the 20 generated instances of benchmark.
IndexNNCKRFIndexNNCKRFP1201.840.75P11302.520.75P2201.841P12302.521P3202.140.75P13601.840.75P4202.141P14601.841P5202.520.75P15602.140.75P6202.521P16602.141P7301.840.75P17602.520.75P8301.841P18602.521P9302.140.75P191201.840.75P10302.141P201201.841


### 4.4.2. 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 proposedH-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 SC1 represents the set coverSC(H-MOEA, MOEA) and SC2 denotes the set coverSC(MOEA, H-MOEA).Table 5
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same number of generations.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3158±0.06310.1158±0.01510.4375±0.06750.0983±0.00930.4550±0.09670.0883±0.0076P10.2492±0.04870.1508±0.01760.3900±0.09310.1908±0.03670.2942±0.06590.1542±0.0270P20.2933±0.04010.1792±0.02450.3167±0.06470.1683±0.03270.3483±0.05640.2058±0.0427P30.4383±0.04870.2000±0.02000.4667±0.09010.0825±0.01240.2992±0.03540.1200±0.0150P40.2942±0.05300.1350±0.02450.3400±0.06000.1542±0.02210.4517±0.09740.1400±0.0264P50.2892±0.05740.1250±0.01450.4283±0.08170.1608±0.02740.4267±0.08200.1758±0.0416P60.4700±0.04060.1008±0.01450.5475±0.07590.1458±0.02560.6250±0.08700.0558±0.0068P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3517±0.05520.2583±0.06530.3542±0.06450.2800±0.05340.3608±0.06940.2667±0.0302P90.3858±0.07470.2317±0.03540.2600±0.02900.2592±0.02590.3864±0.07430.2218±0.0189P100.4625±0.06530.1442±0.03960.4458±0.07770.1733±0.03450.4408±0.09790.1442±0.0299P110.5025±0.09810.1650±0.04320.3608±0.06290.1500±0.02810.4127±0.05030.2967±0.0363P120.3258±0.04940.2258±0.04460.3025±0.05830.2650±0.04550.3325±0.06470.3033±0.0510P130.3033±0.04270.2367±0.04150.2508±0.04410.3242±0.07800.3783±0.06730.2433±0.0610P140.3408±0.06360.2567±0.04560.3483±0.07160.2992±0.05590.2908±0.05620.2974±0.0758P150.3983±0.08030.2492±0.04810.3358±0.08140.2325±0.06630.3175±0.06120.2608±0.0549P160.4392±0.07810.2383±0.05350.3292±0.06270.2708±0.03980.3667±0.09270.3017±0.0765P170.3017±0.06930.3992±0.07270.3500±0.06800.3225±0.07270.4075±0.09740.2642±0.0780P180.3400±0.04750.2833±0.04470.3508±0.08340.3433±0.07000.3792±0.05980.2683±0.0440P190.3842±0.05280.2817±0.07200.3983±0.09120.2375±0.06940.3375±0.08040.2533±0.0740P200.3858±0.06490.1375±0.02710.3508±0.10010.2392±0.07200.3075±0.05160.2358±0.0592In Table5, 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 coverSC2 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.25, H-MOEA and MOEA approximately had the same 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.


### 4.4.3. 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. 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 coverRSC and relative computation time RCT, to indicate the performance of the proposed H-MOEA. RSC and RCT are calculated as follows:(4.1)RSC=SC(H-MOEA,MOEA)SC(MOEA,H-MOEA),(4.2)RCT=CTH-MOEACTMOEA,
where CTH-MOEA denotes the computation time of the proposed H-MOEA and CTMOEA represents the computation time of the normal MOEA.We tested the effect of local search ratepls 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 pls. The values of relative set cover RSC with variation of pls 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.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
The obtained relative set cover (RSC) with different local search rates for P0, P8, and P15 with a perturbation strength 0.15.The definition of relative set coverRSC 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.


### 4.4.4. 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.Table6 shows the computational results obtained under the same examined solutions, while the results under the same computation time are reported in Table 7. After careful examination of the data, it could be found that the proposed H-MOEA outperformed normal MOEA on small size instances, while the performance of H-MOEA decreased on relative large scale cases. By looking at the instances P0–P12, the data suggests H-MOEA performed better than normal MOEA in most instances, under both the same examined solutions (see Table 6) and the same computation time (see Table 7). For relative large scale instances (P13–P20), shown as in Table 6, H-MOEA still outperformed MOEA in most cases under the same examined solutions. While under the same computation time, the performances of H-MOEA and MOEA were comparable on instances P13–P20, suggested as in Table 7. The main reason of performance decrease of H-MOEA on large scale instances could be considered as limitation on global search ability of the algorithm. As indicated in [37], the number of generations is decreased when the available computation time is limited. As a result, the exploration in search space is not fully utilized.Table 6
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same examined solutions.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3692±0.09390.0983±0.00830.4425±0.09000.1108±0.01470.4117±0.07560.1017±0.0136P10.3275±0.05590.1267±0.02390.3992±0.09660.1425±0.02420.2800±0.05360.1842±0.0313P20.2792±0.03800.1458±0.01700.3025±0.04910.1633±0.02190.3233±0.04760.2242±0.0446P30.3392±0.04140.2383±0.02870.4067±0.08830.1392±0.02530.3392±0.07260.1175±0.0207P40.2833±0.04690.1883±0.04230.2775±0.03350.1925±0.03170.4492±0.10740.2408±0.0659P50.3100±0.06590.1492±0.02050.4158±0.07480.1908±0.03070.4150±0.05420.1583±0.0273P60.4433±0.05820.1292±0.02930.5600±0.12260.2292±0.06350.4800±0.07320.1175±0.0188P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3075±0.05260.2942±0.04870.3208±0.07770.3333±0.07740.3075±0.05790.3117±0.0492P90.2808±0.04280.2617±0.03220.3433±0.05760.2408±0.03350.3627±0.03480.2769±0.0198P100.3117±0.07410.2233±0.02990.3758±0.05010.1775±0.02280.3483±0.05850.2092±0.0564P110.5037±0.10870.1650±0.04620.3642±0.05240.1558±0.02970.3128±0.02650.2364±0.0127P120.2883±0.05970.2983±0.06090.3167±0.06190.2700±0.04470.3542±0.05400.3500±0.0780P130.2458±0.05230.3183±0.07460.2475±0.03480.3458±0.08450.3733±0.08950.2883±0.0681P140.3642±0.07910.2525±0.05340.3175±0.08260.3242±0.06320.3726±0.06350.3214±0.0468P150.3558±0.09460.3033±0.06190.3162±0.04860.3210±0.03310.4058±0.05960.3374±0.0288P160.3225±0.05030.2783±0.04550.3433±0.08800.2792±0.06290.2700±0.04800.3750±0.0840P170.2717±0.05740.3667±0.06610.3292±0.05750.3033±0.04880.2833±0.06540.3767±0.0676P180.3483±0.07070.2417±0.05220.3150±0.06740.3925±0.06870.3842±0.06940.2258±0.0473P190.3458±0.07270.3575±0.07900.3650±0.07790.3142±0.05720.3150±0.06090.3125±0.0863P200.3254±0.06490.2384±0.03410.2704±0.08010.2137±0.05200.2674±0.06130.2059±0.0462Table 7
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same computational time.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.2800±0.05510.1367±0.00950.3800±0.04280.1083±0.00820.4775±0.07330.0892±0.0117P10.3025±0.06840.1633±0.02880.3858±0.08450.0917±0.01150.2250±0.04070.2517±0.0365P20.2075±0.01280.2200±0.01630.2850±0.04620.1833±0.02130.3233±0.04240.2217±0.0280P30.4133±0.07860.2025±0.02850.3850±0.06820.1067±0.01300.3108±0.04270.1233±0.0069P40.2617±0.05010.1917±0.02900.3183±0.04710.1842±0.03050.3575±0.06160.2092±0.0606P50.2125±0.01610.1525±0.01860.3650±0.06910.1942±0.04060.3833±0.08670.1975±0.0425P60.4475±0.07460.1283±0.03620.5933±0.07460.1475±0.03230.4800±0.07440.1417±0.0282P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.2983±0.06450.2808±0.05990.4275±0.06940.2208±0.04250.3133±0.04450.3250±0.0488P90.2883±0.04860.2967±0.04680.2833±0.04150.2308±0.01790.3015±0.03120.2453±0.0203P100.3524±0.04520.29432±0.02980.3453±0.06740.30343±0.04460.4018±0.09060.3472±0.0439P110.4725±0.10440.2175±0.06180.4450±0.08400.1717±0.03080.4531±0.03680.2741±0.0153P120.3008±0.04890.2617±0.05550.3000±0.03720.3025±0.05750.3383±0.06740.3433±0.0684P130.2158±0.04720.4250±0.05940.2267±0.04000.3825±0.06510.3117±0.05790.3158±0.0816P140.3967±0.07420.2583±0.04660.3458±0.06990.3513±0.09920.3369±0.02640.2867±0.0224P150.3350±0.04880.2750±0.04590.3452±0.04310.3501±0.05230.3940±0.05110.4012±0.0634P160.3417±0.08650.3400±0.09710.3233±0.05250.2692±0.06230.3267±0.08580.3350±0.0696P170.2875±0.06560.2904±0.05570.2575±0.06500.2467±0.06550.3183±0.05880.3258±0.0854P180.3342±0.06100.3450±0.07750.2575±0.06380.2592±0.06080.2908±0.05740.3000±0.0497P190.2883±0.04510.3267±0.07450.2683±0.07120.3458±0.07860.2567±0.05570.4400±0.0887P200.2978±0.04470.3004±0.06310.2503±0.06010.2294±0.05230.3772±0.07190.3685±0.0793


## 4.4.1. Test Problem Instances

Besides the hypothetical problem described in the above section (denoted as P0), we tested our approach on the benchmark problems generated withProgen, 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.Table 4
Parameter settings of the 20 generated instances of benchmark.
IndexNNCKRFIndexNNCKRFP1201.840.75P11302.520.75P2201.841P12302.521P3202.140.75P13601.840.75P4202.141P14601.841P5202.520.75P15602.140.75P6202.521P16602.141P7301.840.75P17602.520.75P8301.841P18602.521P9302.140.75P191201.840.75P10302.141P201201.841


## 4.4.2. 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 proposedH-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 SC1 represents the set coverSC(H-MOEA, MOEA) and SC2 denotes the set coverSC(MOEA, H-MOEA).Table 5
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same number of generations.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3158±0.06310.1158±0.01510.4375±0.06750.0983±0.00930.4550±0.09670.0883±0.0076P10.2492±0.04870.1508±0.01760.3900±0.09310.1908±0.03670.2942±0.06590.1542±0.0270P20.2933±0.04010.1792±0.02450.3167±0.06470.1683±0.03270.3483±0.05640.2058±0.0427P30.4383±0.04870.2000±0.02000.4667±0.09010.0825±0.01240.2992±0.03540.1200±0.0150P40.2942±0.05300.1350±0.02450.3400±0.06000.1542±0.02210.4517±0.09740.1400±0.0264P50.2892±0.05740.1250±0.01450.4283±0.08170.1608±0.02740.4267±0.08200.1758±0.0416P60.4700±0.04060.1008±0.01450.5475±0.07590.1458±0.02560.6250±0.08700.0558±0.0068P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3517±0.05520.2583±0.06530.3542±0.06450.2800±0.05340.3608±0.06940.2667±0.0302P90.3858±0.07470.2317±0.03540.2600±0.02900.2592±0.02590.3864±0.07430.2218±0.0189P100.4625±0.06530.1442±0.03960.4458±0.07770.1733±0.03450.4408±0.09790.1442±0.0299P110.5025±0.09810.1650±0.04320.3608±0.06290.1500±0.02810.4127±0.05030.2967±0.0363P120.3258±0.04940.2258±0.04460.3025±0.05830.2650±0.04550.3325±0.06470.3033±0.0510P130.3033±0.04270.2367±0.04150.2508±0.04410.3242±0.07800.3783±0.06730.2433±0.0610P140.3408±0.06360.2567±0.04560.3483±0.07160.2992±0.05590.2908±0.05620.2974±0.0758P150.3983±0.08030.2492±0.04810.3358±0.08140.2325±0.06630.3175±0.06120.2608±0.0549P160.4392±0.07810.2383±0.05350.3292±0.06270.2708±0.03980.3667±0.09270.3017±0.0765P170.3017±0.06930.3992±0.07270.3500±0.06800.3225±0.07270.4075±0.09740.2642±0.0780P180.3400±0.04750.2833±0.04470.3508±0.08340.3433±0.07000.3792±0.05980.2683±0.0440P190.3842±0.05280.2817±0.07200.3983±0.09120.2375±0.06940.3375±0.08040.2533±0.0740P200.3858±0.06490.1375±0.02710.3508±0.10010.2392±0.07200.3075±0.05160.2358±0.0592In Table5, 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 coverSC2 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.25, H-MOEA and MOEA approximately had the same 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.


## 4.4.3. 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. 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 coverRSC and relative computation time RCT, to indicate the performance of the proposed H-MOEA. RSC and RCT are calculated as follows:(4.1)RSC=SC(H-MOEA,MOEA)SC(MOEA,H-MOEA),(4.2)RCT=CTH-MOEACTMOEA,
where CTH-MOEA denotes the computation time of the proposed H-MOEA and CTMOEA represents the computation time of the normal MOEA.We tested the effect of local search ratepls 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 pls. The values of relative set cover RSC with variation of pls 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.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
The obtained relative set cover (RSC) with different local search rates for P0, P8, and P15 with a perturbation strength 0.15.The definition of relative set coverRSC 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.


## 4.4.4. 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.Table6 shows the computational results obtained under the same examined solutions, while the results under the same computation time are reported in Table 7. After careful examination of the data, it could be found that the proposed H-MOEA outperformed normal MOEA on small size instances, while the performance of H-MOEA decreased on relative large scale cases. By looking at the instances P0–P12, the data suggests H-MOEA performed better than normal MOEA in most instances, under both the same examined solutions (see Table 6) and the same computation time (see Table 7). For relative large scale instances (P13–P20), shown as in Table 6, H-MOEA still outperformed MOEA in most cases under the same examined solutions. While under the same computation time, the performances of H-MOEA and MOEA were comparable on instances P13–P20, suggested as in Table 7. The main reason of performance decrease of H-MOEA on large scale instances could be considered as limitation on global search ability of the algorithm. As indicated in [37], the number of generations is decreased when the available computation time is limited. As a result, the exploration in search space is not fully utilized.Table 6
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same examined solutions.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.3692±0.09390.0983±0.00830.4425±0.09000.1108±0.01470.4117±0.07560.1017±0.0136P10.3275±0.05590.1267±0.02390.3992±0.09660.1425±0.02420.2800±0.05360.1842±0.0313P20.2792±0.03800.1458±0.01700.3025±0.04910.1633±0.02190.3233±0.04760.2242±0.0446P30.3392±0.04140.2383±0.02870.4067±0.08830.1392±0.02530.3392±0.07260.1175±0.0207P40.2833±0.04690.1883±0.04230.2775±0.03350.1925±0.03170.4492±0.10740.2408±0.0659P50.3100±0.06590.1492±0.02050.4158±0.07480.1908±0.03070.4150±0.05420.1583±0.0273P60.4433±0.05820.1292±0.02930.5600±0.12260.2292±0.06350.4800±0.07320.1175±0.0188P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.3075±0.05260.2942±0.04870.3208±0.07770.3333±0.07740.3075±0.05790.3117±0.0492P90.2808±0.04280.2617±0.03220.3433±0.05760.2408±0.03350.3627±0.03480.2769±0.0198P100.3117±0.07410.2233±0.02990.3758±0.05010.1775±0.02280.3483±0.05850.2092±0.0564P110.5037±0.10870.1650±0.04620.3642±0.05240.1558±0.02970.3128±0.02650.2364±0.0127P120.2883±0.05970.2983±0.06090.3167±0.06190.2700±0.04470.3542±0.05400.3500±0.0780P130.2458±0.05230.3183±0.07460.2475±0.03480.3458±0.08450.3733±0.08950.2883±0.0681P140.3642±0.07910.2525±0.05340.3175±0.08260.3242±0.06320.3726±0.06350.3214±0.0468P150.3558±0.09460.3033±0.06190.3162±0.04860.3210±0.03310.4058±0.05960.3374±0.0288P160.3225±0.05030.2783±0.04550.3433±0.08800.2792±0.06290.2700±0.04800.3750±0.0840P170.2717±0.05740.3667±0.06610.3292±0.05750.3033±0.04880.2833±0.06540.3767±0.0676P180.3483±0.07070.2417±0.05220.3150±0.06740.3925±0.06870.3842±0.06940.2258±0.0473P190.3458±0.07270.3575±0.07900.3650±0.07790.3142±0.05720.3150±0.06090.3125±0.0863P200.3254±0.06490.2384±0.03410.2704±0.08010.2137±0.05200.2674±0.06130.2059±0.0462Table 7
The mean value and the variance of theset coversSC(H-MOEA,MOEA)(SC1) and SC(MOEA,H-MOEA)(SC2) within the scenarios of perturbation strengths 0.15, 0.25, and 0.35, under the same computational time.
St = 0.15St = 0.25St = 0.35SC1SC2SC1SC2SC1SC2P00.2800±0.05510.1367±0.00950.3800±0.04280.1083±0.00820.4775±0.07330.0892±0.0117P10.3025±0.06840.1633±0.02880.3858±0.08450.0917±0.01150.2250±0.04070.2517±0.0365P20.2075±0.01280.2200±0.01630.2850±0.04620.1833±0.02130.3233±0.04240.2217±0.0280P30.4133±0.07860.2025±0.02850.3850±0.06820.1067±0.01300.3108±0.04270.1233±0.0069P40.2617±0.05010.1917±0.02900.3183±0.04710.1842±0.03050.3575±0.06160.2092±0.0606P50.2125±0.01610.1525±0.01860.3650±0.06910.1942±0.04060.3833±0.08670.1975±0.0425P60.4475±0.07460.1283±0.03620.5933±0.07460.1475±0.03230.4800±0.07440.1417±0.0282P70.3417±0.06900.1992±0.04610.4250±0.07620.1958±0.05790.5083±0.09170.1550±0.0479P80.2983±0.06450.2808±0.05990.4275±0.06940.2208±0.04250.3133±0.04450.3250±0.0488P90.2883±0.04860.2967±0.04680.2833±0.04150.2308±0.01790.3015±0.03120.2453±0.0203P100.3524±0.04520.29432±0.02980.3453±0.06740.30343±0.04460.4018±0.09060.3472±0.0439P110.4725±0.10440.2175±0.06180.4450±0.08400.1717±0.03080.4531±0.03680.2741±0.0153P120.3008±0.04890.2617±0.05550.3000±0.03720.3025±0.05750.3383±0.06740.3433±0.0684P130.2158±0.04720.4250±0.05940.2267±0.04000.3825±0.06510.3117±0.05790.3158±0.0816P140.3967±0.07420.2583±0.04660.3458±0.06990.3513±0.09920.3369±0.02640.2867±0.0224P150.3350±0.04880.2750±0.04590.3452±0.04310.3501±0.05230.3940±0.05110.4012±0.0634P160.3417±0.08650.3400±0.09710.3233±0.05250.2692±0.06230.3267±0.08580.3350±0.0696P170.2875±0.06560.2904±0.05570.2575±0.06500.2467±0.06550.3183±0.05880.3258±0.0854P180.3342±0.06100.3450±0.07750.2575±0.06380.2592±0.06080.2908±0.05740.3000±0.0497P190.2883±0.04510.3267±0.07450.2683±0.07120.3458±0.07860.2567±0.05570.4400±0.0887P200.2978±0.04470.3004±0.06310.2503±0.06010.2294±0.05230.3772±0.07190.3685±0.0793


## 5. 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 knowledge-based 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.



---

*Source: 786923-2012-03-14.xml*