A resource flow-based branch-and-bound algorithm to solve fuzzy stochastic resource-constrained project scheduling problem

In this paper, a resource flow-based branch-and-bound procedure is designed to solve the well-known resource-constrained project scheduling problem under the mixed uncertainty of fuzziness and randomness (FS-RCPSP). The objective is to minimize the expected makespan of the project subject to precedence and resource constraints. The proposed branch-and-bound can be employed to obtain optimal solutions and also can be truncated in order to find promising near optimal solutions. The depth-first strategy is utilized for constructing the search tree, and earliest start time (EST) concept is adopted for selecting a node for further branching while traversing the tree down to the leaves. The performance of developed branch-and-bound is benchmarked against CPLEX and SADESP across an extensive set of 960 problems. The results returned by the proposed algorithm show experimentally its effectiveness to solve the FS-RCPSP.


Introduction
One of the important issues when planning a project is the scheduling of some identified activities. Our work is concerned with a well-known problem in the area of scheduling, named resource-constrained project scheduling problem (RCPSP), under the uncertain environment. The objective of RCPSP is to minimize the project makespan subject to precedence and resource constraints. Each activity needs a fixed amount of renewable resources in its processing time. The availability of resources and the activity durations are assumed to be known and fixed, and the preemption is not allowed for the activities. A real-life example to this problem is a construction project where there are precedence relations between the activities and the availability of some resources, for example, human resources, in some time periods is limited. Many attempts have been made to propose exact, heuristic, and metaheuristic procedures to solve the RCPS problem with deterministic parameters. Examples of such research efforts are Brucker et al. (1998), Damay et al. (2007), Pantouvakis and Manoliadis (2006), Ying et al. (2009), Mahdi Mobini et al. (2009), Zamani (2011), Sebt et al. (2013), Kumar and Vidyarthi (2016), and Sheng et al. (2019).
In most of the real-world situations, since the scheduling parameters of a project cannot be exactly estimated, it would be more appropriate if in addition to precedence and resource constraints, uncertainties are also taken into account in the scheduling process. The majority of research efforts to cope with the uncertain information in the project scheduling problems introduce reactive, proactive (robust), stochastic, or fuzzy scheduling approaches (Herroelen and Leus 2005). In the reactive scheduling, the uncertainty is not considered directly at the time of generating a schedule, but the baseline schedule is repaired or rescheduled when disruptions occur. However, in proactive scheduling, a baseline schedule protected as much as possible against disruptions is generated. Some recently published papers in reactive and/or proactive project scheduling are Sadeh et al. (1993), Van de Vonder et al. (2007), Lambrechts et al. (2008), Leus and Herroelen (2004), Artigues et al. (2013), and Bruni et al. (2017). In contrast with the proactive scheduling, in the stochastic scheduling no baseline schedule is produced and all the scheduling decisions are usually made based on a dynamic decision process, called policy, at the time of execution. In the stochastic project scheduling, it is assumed that activities have been performed several times in the previous projects, and therefore, there is enough historical data to provide a probability distribution for each activity. However, in fuzzy project scheduling, due to the lack of historical data, activity durations are estimated by some experts. Since these estimations are vague and imprecise, activity durations are represented by fuzzy numbers. Some papers being published in the area of Stochastic-RCPSP (S-RCPSP) or Fuzzy-RCPSP (F-RCPSP) are as follows: Ashtiani et al. (2011), Fang et al. (2015, Rostami et al. (2018), Atli and Kahraman (2012), and finally Atli and Kahraman (2013), Zoulfaghari et al. (2016). In real-world projects, especially in construction projects, it is possible to have both kinds of uncertainties (i.e., randomness and fuzziness) in durations of activities. In these projects, the duration of some activities (activities with enough historical data about the duration) can be modeled as a probability distribution and that of others (activities without enough historical data about the duration) can only be estimated imprecisely by experts. The corresponding scheduling problem for these projects is called fuzzy stochastic project scheduling problem, and its resource-constrained version (FS-RCPSP) is addressed in the current paper.
Some researchers have coped with the situations where randomness and fuzziness appear simultaneously and have proposed mathematical tools to model this kind of uncertainty. Some of these tools are as follows: probability of a fuzzy event (Zadeh 1968), probabilistic set (Hirota 1981), random fuzzy variable (Liu 2002a(Liu , 2002b, and fuzzy random variable (Puri and Ralescu 1986;Kruse and Meyer 1987;Buckley 2005). Note that the mathematical tools 'random fuzzy variable' and 'fuzzy random variable' are different from each other. Roughly speaking, random fuzzy variables are fuzzy variables taking 'random variable' values and fuzzy random variables are random variables taking 'fuzzy variable' values (Zhu and Liu 2005). An early application of the concept of a fuzzy random variable was in scheduling. Itoh and Ishii (2005) assumed that uncertainty only exists in activity due dates and employed fuzzy random variables for modeling these uncertainties in an n-job machine. A fuzzy stochastic project scheduling problem for the first time was addressed by Ke and Liu (2007). Using random fuzzy variables for modeling of uncertainties, they proposed three models for solving the problem under study: expected cost minimization model, (a,b)-cost minimization model, and chance maximization model. They also proposed a hybrid intelligent algorithm for solving these models. Huang et al. (2009) proposed an expected cost model with random fuzzy variables to handle the software project scheduling problem and solved this model using an intelligent genetic algorithm. Nematian et al. (2010) were the first researchers to tackle the Fuzzy Stochastic-RCPSP (FS-RCPSP) by considering the ready time, duration, and deadline of activities to be fuzzy random variables and converting the resulted model with fuzzy random variables into a mixed-integer linear programming model. Xu and Zhang (2012) studied a more complicated problem, namely resource-constrained multiple project scheduling problem when randomness and fuzziness co-exist. The problem was formulated as a multiobjective model with fuzzy random variables for minimizing the total project time and the total tardiness penalty of all projects. This model was solved by a hybrid genetic algorithm with a fuzzy logic controller. Chen and Zhang (2016) presented a proactive strategy for Preemption-RCPSP with fuzzy random duration and resource availabilities. They proposed a mathematical formulation model for this problem with the objective function of maximizing the sum of free slack and then transformed it into an equivalent crisp model. Finally, they demonstrated the efficiency of the model and proactive strategy just by a numerical example. Alipouri et al. (2017) addressed the RCPSP under the fuzzy random environment using the fuzzy probability theory of Buckley (2005). They proposed a mathematical programming model with fuzzy random variables exploiting the concept of the resource flow network. Then, this model was transformed into an MILP model with crisp parameters and variables. They implemented this MILP model in AIMMS (2014) modeling software and solved 960 benchmark problems using the CPLEX 12.6.0.1 solver of this software. Their results showed that the CPLEX solver in AIMMS could work well on most of the benchmark problems; however, this solver could not find any integer solution for about 10% of the problems with 30 activities and 31% of the problems with 60 activities. In addition, the computational efforts required for solving some of the benchmark problems were high. All these drawbacks can prevent practitioners to use CPLEX solver in real world applications for solving FS-RCPSP. More recently, Alipouri et al. (2019) proposed a selfadaptive differential evolution hyper-heuristic algorithm, named SADESP, to solve the same problem and model tackled by Alipouri et al. (2017). They also solved the same 960 benchmark problems as Alipouri et al. (2017) and reported the results for 500, 5000, 50,000, and 100,000 computed schedules by the algorithm. Unlike the CPLEX, SADESP was able to find integer solutions for all 960 instances, such that with 100,000 computed schedules, SADESP worked well (i.e., found better, equal or close objective values) than CPLEX on about 96.88% (i.e., 465 out of 480) of the J30 problems and 95.21% (i.e., 457 out of 480) of the J60 problems.
All the above-mentioned drawbacks of CPLEX in solving the MILP model of Alipouri et al. (2017) motivated us to introduce a new exact search methodology based on branch-and-bound (B&B) algorithm which can act better than CPLEX and also SADESP in driving policies for FS-RCPSP. Due to the following two reasons, branch-andbound was chosen for the current research: 1-FS-RCPSP belongs to the set of complex combinatorial problems and very often branch-and-bound is the only exact technique which can optimally solve these problems within an acceptable computational effort (Demeulemeester and Herroelen 2002), and 2-branch-and-bound can be truncated and used as a heuristic to find a near optimal solution within a reasonable time for the large-scale problems having many activities. Branch-and-bound has been used by many researchers to solve RCPSP and its variants. For example, Brucker et al. (1998) and Dorndorf et al. (2000) solve RCPSP with a B&B, Davari and Demeulemeester (2018) propose a B&B for the chance-constrained RCPSP, Kolisch et al. (1995), Sprecher et al. (1997), Sprecher and Drexl (1998), and Altintas and Azizoglu (2020) present a B&B for a multi-mode project scheduling problem with a single non-renewable resource, and so on. To the best of our knowledge, there is not any research to propose a B&B method for FS-RCPSP.
Author of the current study believes the contribution of this paper is threefold. First, the problem is modeled as a fuzzy stochastic version of RCPSP (FS-RCPSP) handling the two key sources of uncertainty in a unified way. Second, to the best of the author's knowledge, it is the first effort to propose a B&B algorithm to solve project scheduling problems in an uncertain environment. Third, to the best of the author's knowledge, this is the first time in the literature that the resource flow concept is employed in a B&B procedure for solving a variant of RCPSP.
The remainder of this paper is organized as follows. In the next section, some preliminaries on resource flow and fuzzy probability theories are presented. In Sect. 3, the FS-RCPSP problem description and its model formulations are provided. In Sect. 4, our proposed resource flow-based branch-and-bound algorithm for solving FS-RCPSP is presented. In Sect. 5, the results of computational experiments to test the potency of our method in solving the FS-RCPSP are reported. Finally, in Sect. 6, concluding remarks are drawn out and future works are presented.

Basic concepts and definitions
In this section, some general information is given on resource flow theory. Then, we review some key concepts in Buckley's (2005) approach to modeling uncertainty using fuzzy probability theory and fuzzy random variables.
Hereafter, we place a ''bar'' and a ''tilde'' over a letter to denote a fuzzy number and a fuzzy random number, respectively.

Resource flow theory
In any feasible schedule, upon the completion of an activity i its resources must be allocated to other activity/activities. The key idea of resource flow theory is to determine the number of resource units of type l that are transferred from the end of one activity to the beginning of another activity. Based on these resource links between the activities, a network, named resource flow network, can be constructed for each resource type. Therefore, this network depicts how units of a resource type are passed on between the activities in the schedule. In RCPSP, the problem of resource conflicts existing among the activities of the project can be solved by introducing appropriate resource flow networks.

Fuzzy probability theory
A discrete probability function P on all subsets of a finite set U ¼ fu 1 ; . . .; u s g is defined as Pðfu w gÞ ¼ a w , 0\a w \1 all w and P s w¼1 a w ¼ 1. This probability function is able to handle the uncertainties arising from natural variability (randomness) but cannot deal with the uncertainties due to imprecise information (fuzziness). In order to handle both kinds of uncertainty, Buckley (2005) substitutes a fuzzy number a w ¼ Pðfx w gÞ for each crisp number a w and assumes that 0\a w \1 for all w. The discrete set U together with the a w values is a discrete fuzzy probability distribution (Buckley 2005). Buckley (2005) considers the restriction ''there are a w 2 a wa ; a ¼ 1 so that P s w¼1 a w ¼ 1'' on the a w values and accordingly defines the fuzzy probability of C ¼ fu 1 ; . . .; u k g, 1 k\s, and the expected value of a discrete fuzzy probability distribution by their a-cut as follows: Let e U be a fuzzy random variable and let f ðu; hÞ be a fuzzy probability density, where u 2 IR and h ¼ h 1 ; . . .; h Q È É is for parameters h q , 1 q Q. The fuzzy probability of e U 2 ½z 1 ; z 2 and the expected value of fuzzy probability density f ðu; hÞ are defined by their a-cut as follows: ; for 0 a 1: In this research, the same as Alipouri et al. (2017Alipouri et al. ( , 2019, fuzzy random variables with fuzzy normal density functions are employed to represent the durations of activities. Moreover, all fuzzy parameters are represented by triangular fuzzy numbers (TFNs). A fuzzy normal density is shown as Nðl; r 2 Þ, wherein only the values of l and r 2 have become fuzzy compared to the crisp normal density Nðl; r 2 Þ. Buckley (2005) provides the result that the fuzzy mean of Nðl; r 2 Þ equals to l and its fuzzy variance is r 2 .

Problem description and formulation
The problem under study is FS-RCPSP, which is a variant of RCPSP considered in an uncertain environment in which randomness and fuzziness exist simultaneously. Here, we just consider the uncertainties existing in durations of activities, owing to the fact that extending the proposed algorithm to another algorithm which considers the uncertainties in other parameters is straightforward.
In FS-RCPSP, a single project consisting of n ? 2 activities is considered. The activities are numbered 0 to n þ 1, where 0th and (n þ 1)th activities are start and end dummy activities, respectively. Each activity j cannot be interrupted once in progress (i.e., preemption is not allowed) and has to be started after all its immediate predecessor activities i (i 2 IP j ) have finished (i.e., precedence constraint). We have L renewable resources (e.g., equipment and human resources), and each resource l 2 L has a limited capacity R l (1 l L) throughout the project duration. Each activity j requires r jl (0 j n þ 1, 1 l L) units of resource l once in progress. The sum of resource requirements for resource l at any time period t cannot exceed R l (i.e., resource constraint). The processing time of activity j is denoted asd j (0 j n þ 1) which is a fuzzy random variable. The start and finish time of each activity j are, respectively, shown bys j andf j (0 j n þ 1). The objective is to find precedence and resource feasible completion times for all activities which lead to the minimum expected makespan.
Based on the description above, FS-RCPSP can be formulated as follows: Minf nþ1 s:t: ð5Þ where E is the set of all precedence relations and B(t) is the set of activities being processed at the time period t. In this formulation, constraints (6) and (8) are, respectively, used to impose precedence and resource constraints. Due to the lack of a suitable approach for transferring set BðtÞ to a linear constraint and also due to the presence of fuzzy random parameters and variables, the above mathematical programming model cannot be solved directly. Therefore, Alipouri et al. (2017) x ij þ x ji 1 8ði; jÞ 2 V 2 ; i 6 ¼ j ð12Þ where x ij is a binary variable denoting that activity j is started immediately after the completion of activity i whenever x ij ¼ 1; otherwise, x ij ¼ 0. C l ij denotes the amount of resource l directly transferred from activity i to activity j. V and E are the set of all activities and the set of all precedence relations, respectively. The objective function (10) minimizes the expected completion time of the end dummy activity and consequently the expected completion time of the project, l f nþ1 ¼ ðf nþ1;1 ; f nþ1;2 ; f nþ1;3 Þ. To obtain this objective value, Alipouri et al. (2017) adopted the technique described in Buckley and Feuring (2000). First, the problem of minimizing the fuzzy number l f nþ1 was converted into a multi-objective problem; then, the multi-objective problem was converted into a single objective problem with the use of weights This is a flexible objective function for project managers since they can vary the values of weights to satisfy their different requirements. Equation (11) introduces the precedence relations between the activities. Constraint (12) is transiting constraint, and constraint (13) ensures that no cycles will exist in the network. Constraints (14)-(19) are employed for setting the completion times of activities. Constraints (20)-(22) are resource flow inequalities. By constraint (22), the resource flow values are limited to minðr il ; r jl Þ of arc ði; jÞ if the arc exists. Constraints (20) and (21) have been devised to ensure that the incoming flow on node i is equal to the outgoing flow from that node. In the next section, the proposed branch-and-bound procedure for solving this problem is presented.

The proposed branch-and-bound algorithm
In what follows, a branch-and-bound procedure is described to solve the FS-RCPSP. This algorithm can be employed to obtain optimal solutions and also can be truncated in order to find promising near optimal solutions. In Sect. 4.1, the branching strategy is presented. Section 4.2 is devoted to the description of five simple and effective rules for pruning the solution tree. The algorithmic description of the proposed B&B for solving FS-RCPSP is summarized in Sect. 4.3. Finally, in Sect. 4.4, an illustrative numerical example is given.

Branching and node selection
In this subsection, we describe how a solution tree is built by the proposed B&B. In other words, our strategy in creating new nodes of the tree and selecting a node for further branching is explained. Based on the literature, there are three main branching strategies for constructing a search tree.
• Best-first strategy, also known as frontier search or skiptracking strategy, which always selects the intermediate node with the lowest bound to branch from next. • Breadth-first strategy, which processes all the nodes at one level of the search tree before going to a higher level. • Depth-first strategy, also known as backtracking strategy, which explores a path by processing node by node down a branch until a node on the path is fathomed or a complete solution is reached. Then, if there exists at least one previously created, unfathomed, intermediate node which has not been fully explored, backtracking is accomplished.
In breadth-first strategy, the number of nodes at each level grows exponentially, leading to avoid using this strategy for complex problems. Also, the same as breadthfirst strategy, the best-first strategy is unable to find a feasible solution before ending the search and requires more storage capacity than the depth-first strategy. However, a depth-first strategy can find a feasible solution right away leading to have a rather good solution if the algorithm is truncated prematurely. All the advantages of depth-first strategy against other strategies encouraged us to employ this strategy for constructing the search tree.
While traversing the tree down to the leaves, a decision has to be made at each level of the tree about selecting a node for further branching. To this end, in this study the earliest start time (EST) concept is employed as explained later in this subsection.
In Sect. 3, it was explained that processing time of each activity i is represented as fuzzy normal density Nðl i ; r 2 i Þ and the values of l i and r 2 i are assumed to be triangular fuzzy numbers. Based on the fuzzy probability theory of Buckley (2005), the expected duration of each activity i with fuzzy normal density Nðl i ; r 2 i Þ is l i . Therefore, the expected makespan is also a triangular fuzzy number which can be calculated using the expected durations of activities, l d i ¼ ðd i;1 ; d i;2 ; d i;3 Þ. Thus, these expected durations will be the base of the calculations in the proposed B&B. Using the expected values of activity durations does not mean to ignore randomness, since these expected values are descriptive quantities calculated based on the probability distributions of fuzzy random variables. It is worthwhile noting that, while traversing on the search tree, if it is needed to compare two triangular fuzzy numbers, their corresponding crisp (Cr) values are compared. The corresponding crisp value of the triangular fuzzy number B ¼ ðb 1 ; b 2 ; b 3 Þ is calculated using the following equation: Each node of the search tree is associated with a partial schedule in which an expected fuzzy finish time is assigned to an activity of the project. Here, based on the concept of resource flow theory, resource precedence relations are defined between the activities that resource transferring has taken place between them. These relations will lead to have feasible partial and/or complete schedules in every node of the search tree (Artigues et al. 2003). To the best of author's knowledge, this is the first time in the literature that the resource flow concept is employed in a B&B procedure for solving a variant of RCPSP.
In the first step of constructing the search tree, the start dummy activity is scheduled at the root node and its expected start and finish times are set to (0,0,0). Moreover, all the available resources are assigned to this dummy activity, and also resource requirements of the end dummy activity are considered to be equal to the available amount of the resources. The branching process is started by adding eligible activities (i.e., unscheduled activities for which all predecessors have scheduled) of the first level, as children, to the root node, as parent. To this end, at this first level, the algorithm probes for all the activities that the start dummy activity is their only predecessor, and stores them in the eligible activities set of the first level, EAS 1 . Then, one node is created for each activity in EAS 1 and this activity is scheduled to start at expected time (0,0,0). The resources required by each of these newly scheduled activities are transferred from the start dummy activity. Among all these newly added nodes, one node has to be selected to branch from it in the next level of search tree. For this selection, an expected EST value and its corresponding crisp value (see Eq. (26)) are calculated for all the activities of the project. The competing activity with the lowest Cr EST value (using the lower activity number as a tie breaker) is selected for further branching.
where EST i;1 , EST i;2 , and EST i;3 are earliest start times of activity i which are, respectively, calculated by applying the critical path method on the most optimistic, most likely, and most pessimistic values of the expected durations of project activities. So far, we could explain the first iteration of the proposed B&B algorithm. From this iteration on until the end of search process, the procedure is somehow different, as some strategies for avoiding the node redundancy are employed by the algorithm to reduce the search space of the problem. After selecting a node (we call this node and its corresponding activity as SN and SA, respectively) from the set EAS m of level m to add to the tree, based on the resource flow network and precedence relations existing on the current path, activities which have been scheduled so far and are able to provide resources to the activity of the selected node are listed. The activities on this list are sorted based on the non-decreasing order of their expected finish time, and from the activity that has the lowest expected finish time it is started to transfer the resources into activity SA. Transferring of resources is continued until this activity has collected all its required resources. The resource transferor activities which are not the predecessor of activity SA are considered as the resource predecessor of this activity. Then, after taking into account the pruning rules presented in the next subsection, if activity SA is not fathomed, this activity is scheduled at the earliest possible expected start time equal to or greater than the expected time of the activity scheduled at the previous level.

Pruning rules
In order to speed up the above proposed branching procedure, it is necessary to introduce some pruning rules which are used to fathom the nodes that branching from them is proved to be unnecessary. Fathoming of a node can happen if it is proved that this node contains non-optimal or infeasible solutions, or by branching from the other intermediate nodes of the search tree a similar or better solution could be obtained (Demeulemeester and Herroelen 2002). In building our proposed search tree, we relied on the following five pruning rules which are simple and effective. These pruning rules are the extended version of the pruning rules which have been proposed in the literature for the precedence tree algorithm.
Rule 1 If a child node is scheduled to start in the same expected start time as its parent and it has ''smaller'' (arbitrarily one can use the word ''larger'') activity number, this child node is fathomed. This is because a same partial solution can be produced in another branch of the tree, this time by scheduling the smaller numbered activity before the one with the higher number.
Rule 2 Assume that activity SA 2 EAS m is going to be scheduled on expected start time l s SN at level m of the B&B tree. If this activity can be performed at a previous level with an expected start time less than l s SN and without violating resource and precedence constraints, the corresponding node of this activity (i.e., SN) is fathomed. This is because of the existence of a node at the previous level that branching from it will result in a same schedule that can be produced by shifting the identified activity/activities.
The other three pruning rules are based on the calculation of lower bounds. The lower bounds calculated for a node are compared with an upper bound. If at least one of the lower bounds of the node exceeds or equals the current upper bound of the search tree, that node is pruned. The initial value of the upper bound is set to the sum of the expected durations of all activities of the project. Each time a complete solution with a better expected makespan is found by the B&B algorithm, the upper bound is changed to be equal to this expected makespan.
Ideally, the best lower bound value for a given node is the value of the best accessible objective value to the problem at hand from this node. However, here obtaining this value is itself NP-hard leading us to relax the problem in order to calculate the lower bounds. Therefore, in this study, after the relaxation of resource constraints, an expected distance matrix, denoted as l DM , that its elements are triangular fuzzy numbers, is calculated as explained in Remark 1 and used in the calculation of lower bounds.
Remark 1 Let l dm ij ¼ ðdm ij;1 ; dm ij;2 ; dm ij;3 Þ be a triangular fuzzy number denoting the expected length of a longest path (expected distance) from activity i to activity j. The Floyed-Warshall algorithm, which has the time complexity of O(n 3 ) (see Lawler (1976)), is used three times to calculate the expected distance matrix l DM ¼ ðl dm ij Þ i;j2V . This algorithm is applied on the most optimistic values of expected durations of activities to calculate the dm ij;1 and on the most likely and most pessimistic values to calculate the dm ij;2 and dm ij;3 , respectively. Rule 3 If a child node is selected to start in the same expected start time as its parent node, a lower bound (LB1) for the child node can be calculated by adding the expected start time of the parent node to the remaining expected distance from the activity of the child node to the end dummy activity.
Rule 4 Let LN m be the set of unscheduled activities having numbers lower than the activity of node SN of the search tree. A lower bound (LB2) for this node is calculated as follows.
where l s SN denotes the expected start time of the activity of the node SN. This rule can be easily justified using Rule 1, where any lower numbered activity prohibited from being scheduled in an expected start time equal to or smaller than that of the currently scheduled, higher numbered activity. Rule 5 Let SN be a given node of the search tree with EAS m denoting the set of eligible activities to be scheduled after SN. A lower bound (LB3) for the node SN can be calculated as follows.

Algorithmic description
The algorithmic description of the proposed B&B for solving FS-RCPSP is given in this subsection. The B&B procedure consists of the following eight steps. Also, the flowchart of the proposed B&B is illustrated in Fig. 1.
Step 1 Read data related to the problem at hand: the number of non-dummy activities to be scheduled (n), expected processing time of each activity i (i = 0, …, n ? 1) as triangular fuzzy number ðd i;1 ; d i;2 ; d i;3 Þ, precedence relations, resource requirement of each activity, resource availabilities, truncating criteria.
Step 2 Calculate the expected distance matrix l DM using expected durations of activities and critical path method calculations as explained in Remark 1.
Step 3 Set the values of x 1 , x 2 , and x 3 to user-defined values. Set the initial upper bound as UB ¼ P n i¼1 l d i . Set the level of the B&B tree m ¼ 0. Calculate the expected EST value for each activity of the project.
Step 4 Schedule the start dummy activity to start and finish at the expected time (0,0,0) and assign all the available resources to this dummy activity. Assume that the resource requirement of dummy activity n?1 is equal to the availability of resources. Update the partial schedule PS={0}.
Step 5 Increment the level number by one and find the eligible activities of first level, EAS 1 . Schedule all the activities of EAS 1 to start at expected time (0,0,0) and transfer the resources required by these activities from start dummy.
Step 6 Select the node with the lowest expected EST value (using lower activity number as a tie breaker) to branch from it at the next level. Update PS with regard to the selected activity. Set m=m?1.
Step 7 WHILE there is at least one unexplored, unfathomed node or truncating criteria are not satisfied, Do Step 7.1 Find the eligible activities EAS m of level m, and select the activity with the lowest expected EST value (using lower activity number as a tie breaker), denoted as SA, for further branching.
Step 7.2 Based on the resource flow and precedence relations existing on the current path, list activities being able to send resources to SA and prioritize them considering their expected finish time such that lowest expected time results in highest priority.
Step 7.3 Starting from the activity with highest priority transfer resources to SA until SA gets all its required resources, consider the transferor activities which are not the predecessor of activity SA as the resource predecessor of this activity.
Step 7.4 Considering five pruning rules presented in section 4.2, if activity SA is not fathomed, go to step 7.5; otherwise go to step 7.6.
Step 7.5 Create a node for activity SA and schedule it at the earliest possible start time equal to or greater than the expected start time of the activity scheduled at the previous level. Set m=m?1, and set PS. If SA is the end dummy activity, update upper bound.
Step 7.6 If there exists at least one unexplored, unfathomed node, backtrack on the search tree to select another node to branch from next and then update m and PS with regard to the location of the selected node; otherwise go to step 8.
Step 8 END WHILE. It is necessary to mention that, based on Alcaraz and Maroto (2001), an optimal solution for RCPSP can only be achieved by the exact solution procedure in small projects, usually with fewer than 60 activities, which are not highly resource constrained. For projects with large-scale activities, there can be two scenarios to apply the proposed B&B: 1-using more powerful computers to solve the project to optimality, and 2-truncating the algorithm based on a user-defined criteria to find a near optimal solution.

Numerical example
In this section, our proposed resource flow-based B&B procedure is used to find an optimal solution for the project instance shown in Fig. 2. In this activity-on-the-node network, each node is the representative of an activity and the numbers above and below each node denote, respectively, the expected duration and resource requirement of that activity. The resource availability for this project equals 7.
Values of x 1 , x 2 , and x 3 are set to 0.1, 0.8, and 0.1, respectively. The initial upper bound is calculated as UB ¼ P 7 i¼1 l d i ¼ ð13; 27; 43Þ. Figure 3 depicts how the B&B tree can be built on this instance. In Fig. 3, the first line of every node shows the node number on the left and the activity number of that node on the right. Moreover, on the second line of every node the expected start and finish times for the activity scheduled at that node are indicated. Each pair of numbers on the last line of every node in Fig. 3 is assigned to a resource flow. The first and second numbers of a pair correspond, respectively, to the number of activity sending the resource to the node and the amount of the resource which is sent. In the search tree of Fig. 3, the fathomed nodes are depicted in gray or in dotted lines, such that the gray nodes are the nodes pruned by pruning rules 1 and/or 2 and the nodes with dotted lines are the nodes pruned because of the lower bound based pruning rules 3, 4, and/or 5.
The search tree initially (i.e., at level 0) consists of the root node in which start dummy activity is scheduled to start and finish at the expected time (0,0,0). The eligible activities of level 1 are EAS 1 = {1,2,3}, which get their required resources from start dummy activity and are scheduled to start at the expected time (0,0,0), see step 5 of the algorithm. The expected EST for all the activities in EAS 1 are equal, so activity 1 which has a lower number than the others is selected to branch from it at level 2. The eligible activities of level 2 are EAS 2 = {2,3,6}; among them activity 2 is selected as SA considering the expected EST values. In regard to the resource flow and precedence relations on the current path, activities 0 (having 5 units of resources) and 1 (having 2 units of resources) can send their resources to activity 2. Considering the expected finish times, since activity 0 has higher priority than activity 1 (see step 7.2 of the algorithm) for sending resources to activity 2, 3 units of resources are transferred from start dummy activity to activity 2. Activity 2 is not fathomed by any of the presented pruning rules, so a node is created for this activity with expected start and finish times equal to (0,0,0) and (3,5,7), respectively. The partial schedule has now three activities PS = {0,1,2}. At level 3, EAS 3 consists of three activities 3, 4, and 6. The priority of activity 3 in regard to the expected EST values is more than the others; therefore, activity 3 is considered as SA. Since this activity is not fathomed by any of the pruning rules, it is scheduled to start at expected time (1,2,3) getting its required resources from root node (2 units of resources) and activity 1 (one unit of resources). Activity 1 is a resource transferor to activity 3 but is not its predecessor; thus, activity 1 is considered as a resource predecessor of activity 3. It is continued to traverse the tree down until a complete feasible schedule is obtained by scheduling the end dummy activity to start and finish at expected time (8,19,31). The upper bound is now updated by setting its value to the expected finish time of end dummy activity, i.e., UB = (8,19,31). Then, since there exist unexplored, unfathomed nodes on the search tree, the backtracking is done and the node 12 on level 5 is selected as SN. However, this node is fathomed because its LB3 equals UB and another node, node 13 in Fig. 3, is selected for evaluation. This search process finishes after 58 nodes have evaluated, and the expected time (8,19,31) is resulted as the minimum value for the expected makespan.

Computational experiments
The algorithm proposed in Sect. 4.3 was implemented in MATLAB programming language on a laptop with Windows 7 operation system, Intel Core 2 Duo, and a CPU at 2.00 GHz. In Sect. 5.1, the tested benchmark problem instances are described. In Sect. 5.2, a sensitivity analysis is carried out to investigate the effects of different settings for the weights x 1 , x 2 , and x 3 on the performance of the B&B. Finally, in Sect. 5.3 the computational results are provided.

Sensitivity analysis
A sensitivity analysis was carried out to investigate the effects of different settings for the weights x 1 , x 2 , and x 3 on the performance of the B&B; 96 benchmark problems selected randomly from the 96 groups of problems are generated in Sect. 5.1. These were solved to optimality by the proposed B&B based on the 10 scenarios shown in Table 1 for the weights x 1 , x 2 , and x 3 .
For each of the scenarios in Table 1, the mean of the defuzzified makespans for both sets J30 and J60 is reported in Fig. 4. The well-known centroid of area (COA) method is employed to defuzzify fuzzy numbers. As shown in Fig. 4, an increase in the value of x 2 (from scenario 1 to 8) leads to a decrease in defuzzified makespans of both sets. However, from scenario 8 to scenario 10, increasing the value of x 2 increases the defuzzified makespans. Therefore, based on this sensitivity analysis, values 0.1, 0.8, and 0.1 are appropriate settings, respectively, for the weights x 1 , x 2 , and x 3 in order to calculate the fitness values in the proposed B&B.

Computational results
The literature review of FS-RCPSP indicated that so far only Nematian et al. (2010), Alipouri et al. (2017), and Alipouri et al. (2019) have proposed a method to solve this problem. Nematian et al. (2010) did not report any results from applying their MILP model on any benchmark problem and just illustrated their proposed procedure by an example with 7 activities. Hence, in this paper, the performance of our proposed B&B algorithm is only compared with the results obtained by Alipouri et al. (2017Alipouri et al. ( , 2019 to solve the 960 benchmark problem instances identified earlier. Comparing against CPLEX and SADESP provides quantitative results that are easily reproducible by other researchers. Based on the sensitivity analysis carried out in Sect. 5.2 and also in the same setting as Alipouri et al. (2017Alipouri et al. ( , 2019, in all experiments, values of x 1 , x 2 , and x 3 are considered to be 0.1, 0.8, and 0.1, respectively. Alipouri et al. (2019) showed based on a sensitivity analysis that this setting for the weights results in minimum defuzzified makespans. The results obtained by applying CPLEX, SADESP (with 100,000 computed schedules) and our proposed algorithm to solve the project instances generated from J30 and J60 sets are presented in Tables 2 and 3. These instances, based  on the input parameters of ProGen, are classified into 96 groups-48 groups for each set-each one containing 10 problems. The number of problems out of 10 that B&B was able to find an integer solution for them, the mean of objective values (Mean of Z), and the mean of computation times (Mean of CPU) are reported for each groups in both tables. In addition, the same results reported in Alipouri et al. (2017) and Alipouri et al. (2019) are presented in these tables under the CPLEX and SADESP columns for comparison. In both tables, for those groups which their results are differentiated in italics from other results, the CPLEX could only find integer solutions for some of the problems. The results of CPLEX for these groups are only for those problem instances for which it was able to find integer solutions. Bold numbers in the columns of Tables 2  and 3 show where either CPLEX, SADESP or B&B has performed better in regard to mean of Z and mean of CPU times. Since Alipouri et al. (2019) has compared the results of SADESP with CPLEX, this comparison is not discussed in this paper. As shown in Table 2, the CPLEX was unable to find an integer solution for all problem instances in 4 groups of J30 set, resulting in a corresponding 'N/A' entry. This was also the case for two and five instances in groups j30_41 and j30_45, respectively. In contrast, SADESP (after 100,000 trials) and B&B were able to find integer solutions for all 480 instances of J30. Compared to the CPLEX, proposed algorithm performed better in terms of the mean of Z in 32 of the 48 groups, and both algorithms were tied in 16 groups. When B&B's results for set J30 are compared to those of SADESP, corresponding figures are 32 and 16, respectively. Where B&B beats CPLEX, the average margin between the mean results was 1.24% with some groups showing larger margins than 2%, e.g., j30_5 (3.44%), j30_10 (3.07%), j30_14 (5.06%), j30_21 (2.82%), j30_30 (5.71%), j30_31 (3.42%), and j30_46 (3.11%). The average margin between the results returned by B&B and SADESP was 0.77% in favor of B&B. Where B&B beats SADESP, the average margin between the mean results was 1.02% with some groups showing larger margins than 2%, e.g., j30_6 (2.26%), j30_9 (4.50%), j30_13 (2.88%), j30_14 (3.76%), j30_29 (2.45%), j30_30 (2.40%), and j30_31 (2.58%). The average margin between the results returned by B&B and SADESP for J30 set was 0.68% in favor of B&B. Based on the results, our algorithm worked well (better or equal) on 100% of the problems in J30 set in terms of the mean of Z.
For J60 set, which has more challenging problems than J30 set, the CPLEX could not find an integer solution for all problem instances in 12 groups having 'N/A' entry in Table 3. In contrast, SADESP (after 100,000 trials) and proposed B&B were again able to find integer solutions for all 480 instances of J60. Compared to the CPLEX, in 32 of the 48 groups for J60 sets, B&B performed better in terms of the mean of Z, while both algorithms were tied in the remaining 16 groups. When B&B's results for set J60 are compared to those of SADESP, corresponding figures are 34 and 14, respectively. Where B&B beats CPLEX, the average margin between the mean results was 0.55% and all the margins were smaller than 2% except for j60_11(3.39%). The average margin between the results returned by the two algorithms was 0.23% in favor of B&B. Where B&B beats SADESP, the average margin between the mean results was 2.59% with 14 groups showing larger margins than 2%. The average margin between the results returned by B&B and SADESP for J60 set was 1.83% in favor of B&B. Based on the results, our algorithm worked well (i.e., better or equal) on 100% of the problems in J60 set in terms of the mean of Z.
CPU time of the algorithms should also be compared to clarify their costs. Since the specifications of the employed operating system for the current study are the same as Alipouri et al. (2017) and Alipouri et al. (2019), a fair comparison can be done easily just by comparing the CPU times reported in Tables 2 and 3. As seen in these tables, in 42 of the 48 groups for J30 and 40 of 48 groups for J60, the computation time taken by CPLEX was higher than B&B, while for the remaining 6 groups for J30 and 8 groups for J60, CPLEX performed better in terms of the mean of CPU times. Also, Tables 2 and 3 show that in 44 of the 48 groups for J30 and 45 of the 48 groups for J60 the mean of the mean of CPU time taken by SADESP (with 100,000 produced schedules) was higher than proposed B&B, while for the remaining 4 groups for J30 and 3 groups for J60, SADESP performed better in terms of the mean of CPU times. But should it be taken into consideration that, based  on the results of Alipouri et al. (2019) the computation time taken by SADESP scales proportionally with the number of schedules evaluated, and the distribution of computation times for SADESP was fairly uniform. Based on the above numerical experimental results on the project instances generated from J30 and J60 sets from PSPLIB, it can be seen that the proposed B&B outperforms SADESP and CPLEX. Meanwhile, the comparison of CPU times of algorithms indicates that in average the proposed B&B is faster than SADESP and CPLEX on 89% of the instances. In short, proposed B&B is more effective and efficient than SADESP and CPLEX in solving FS-RCPSP.

Conclusion
In this study, we addressed the fuzzy stochastic resourceconstrained project scheduling problem (FS-RCPSP), where the objective is to minimize the expected makespan of the project subject to precedence and resource constraints under the mixed uncertainty of fuzziness and randomness. The uncertainties were represented by fuzzy random variables and calculations carried out based on the expected value of the durations. We designed a resource flow-based B&B procedure incorporating with five pruning rules to find optimal solutions for this problem. The proposed B&B was successfully tested on 960 benchmark instances generated by ProGen project schedule generator. Our results indicated that the proposed approach can perform better than CPLEX solver and SADESP to solve the FS-RCPSP problems, which allows its use in practice to handle project schedule uncertainties.
Future research should consider testing the proposed B&B algorithm on a wider variety of benchmark problems in order to validate both the generality and robustness of this approach to fuzzy stochastic resource-constrained project scheduling. Currently, only the mean of the activity durations, modeled as fuzzy random variables, is used in the calculations of the proposed algorithm. Future studies can consider modeling other project parameters (e.g., the resource requirements of activities) as fuzzy random variables and include the variance of uncertain parameters along with the mean value in the model. Moreover, future research should try to introduce scenarios and to find a compromise between the solutions related to every scenario. Proposing more effective pruning rules can also be a good subject for future researches.
Funding The author has no relevant financial disclosures.

Declarations
Conflict of interest The author declares that he has no conflict of interest.
Ethical approval This article does not contain any studies with human participants or animals performed by any of the authors.