Robust execution strategies for project scheduling with unreliable resources and stochastic durations

Abstract

The resource-constrained project scheduling problem with minimum and maximum time lags (RCPSP/max) is a general model for resource scheduling in many real-world problems (such as manufacturing and construction engineering). We consider RCPSP/max problems where the durations of activities are stochastic and resources can have unforeseen breakdowns. Given a level of allowable risk, $\alpha$, our mechanisms aim to compute the minimum robust makespan execution strategy. Robust makespan for an execution strategy is any makespan value that has a risk less than $\alpha$. The risk for a makespan value, $M$ given an execution strategy, is the probability that a schedule instantiated from the execution strategy will not finish before $M$ given the uncertainty over durations and resources. We make three key contributions: (a) firstly, we provide an analytical evaluation of resource breakdowns and repairs on executions of activities; (b) we then incorporate such information into a local search framework and generate execution strategies that can absorb resource and durational uncertainties; and (c) finally, to improve robustness of resulting strategies, we propose resource breakdown aware chaining procedure with three different metrics. This chaining procedure computes resource allocations by predicting the effect of breakdowns on robustness of generated strategies. Experiments show effectiveness of our proposed methods in providing more robust execution strategies under uncertainty.

Appendices

Appendix 1: Proof of Equation 11

From Eq. 8, we can express $\tilde{z}^{+}$ and $\tilde{z}^{-}$ by the following set of equations:

$$\begin{aligned} \begin{array}{l} \tilde{z}^{+}=(\tilde{z}+|\tilde{z}|) / 2, \ \ \ \tilde{z}^{-}=(|\tilde{z}|-\tilde{z}) / 2. \end{array} \end{aligned}$$

(27)

As $\tilde{z}\sim N\{0,\sigma \}$, $|\tilde{z}|$ follows half-normal distribution with the expected value and variance given by,

$$\begin{aligned} \begin{array}{l} E[|\tilde{z}|] = 2\sigma / \sqrt{2\pi }, \ \ var[|\tilde{z}|] = \sigma ^2 (1- 2 / \pi ). \end{array} \end{aligned}$$

(28)

From definitions of $\tilde{z}^{+}$, $\tilde{z}^{-}$ and Eq.  28, we have,

$$\begin{aligned} \begin{array}{l} E[\tilde{z}^{+}]=E[\tilde{z}^{-}]= \sigma / \sqrt{2\pi }, \\ Var[\tilde{z}^{+}]= E\left[ (\tilde{z}^{+})^{2}]-(E[\tilde{z}^{+}]\right) ^2 = (\pi -1)\sigma ^{2} / {2\pi }. \end{array} \end{aligned}$$

(29)

Similarly, $Var[\tilde{z}^{-}] = (\pi -1)\sigma ^{2} / {2\pi }$.

Appendix 2: Proof of Equation 21

To determine values of $E[\tilde{\delta _i}]$ and $Var[\tilde{\delta _i}]$, we need to calculate mean and variance values of $N_i$ and $R_i$. The detailed mathematical derivation will be shown in the following subsections.

1.1 Number of interruptions

In this subsection, we

analyze random variable $N_i$ describing the total number of interruptions due to resource breakdowns experienced by activity throughout its execution.

Lemma 9.1

Hillier and Lieberman (2004) Let $X_1, X_2,\ldots X_n$ be independent exponential random variables with parameters $\lambda _1, \lambda _2,\ldots \lambda _n$, respectively. Then the minimum of these random variables follows an exponential distribution; that is,

$$\begin{aligned} \mathrm{min}\{X_1, X_2,\ldots X_n\} \sim EXP\left( \sum \limits _{i=1}^n \lambda _i\right) . \end{aligned}$$

During execution, activity $a_i$ would be interrupted as soon as one of these resource units allocated to $a_i$ breaks down. Suppose $T_i$ is a random variable representing the minimum time to failure over all resource units allocated to activity $a_i$. We can represent $T_i$ as follows.

$$\begin{aligned} \begin{array}{l} T_i = \mathrm{min}\{X_{11}, \ldots , X_{r_{i1}1}, \ldots , X_{1K}, \ldots , X_{r_{iK}K}\}. \end{array} \end{aligned}$$

(30)

We then rewrite $T_i$ as:

$$\begin{aligned} T_i = \mathrm{min}\{M_1, M_2, \ldots , M_K\}, \end{aligned}$$

(31)

where $M_k$ is a random variable representing the minimum time to breakdown of all resource units of resource type k needed by activity $a_i$, and $M_k = \mathrm{min}\{X_{1k},X_{2k},\ldots X_{r_{ik}k}\}$.

Since $X_{jk} \sim EXP(\lambda _{jk})$, using Lemma 9.1, we have,

$$\begin{aligned} \begin{array}{l} M_k \sim EXP(\sum _{j=1}^{r_{ik}} \lambda _{jk}). \end{array} \end{aligned}$$

(32)

Similarly, following Eqs. 31 and 32 and Lemma 9.1, we can see that $T_i$ is also exponentially distributed with the parameter $\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}$, so that

$$\begin{aligned} \begin{array}{l} E[T_i] = 1/\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}. \end{array} \end{aligned}$$

(33)

Without consideration of resource variability, duration of activity $i$ is represented as the sum of its mean $d^0_i$ and natural variability $\tilde{z}_i$ with zero mean, i.e.,

$$\begin{aligned} \begin{array}{l} E[\tilde{z}_i] = 0. \end{array} \end{aligned}$$

(34)

Let $\{z_i\}$ denote the realization set of durational variability $\tilde{z}_i$. Thus, for a certain realization $z_i$, duration of activity $a_i$ is a deterministic value, i.e., $d^0 + z_i$.

Lemma 9.2

Hillier and Lieberman (2004) Suppose the time between consecutive occurrences of arrivals follows an exponential distribution with parameter $\lambda$. Let X(t) be the number of occurrences by time t, then X(t) follows a Poisson distribution with parameter $\lambda *t$.

Using Lemma 9.2 by Hillier and Lieberman Hillier and Lieberman (2004), the number of interruptions of activity $a_i$ during its execution $d_i^0 + z_i$, follows a Poisson distribution with the mean number of occurrences given by $(d_i^0+z_i)/E[T_i]$. In other words, mean and variance values for the number of interruptions for the realization are:

$$\begin{aligned} E[N_i|z_i]=Var[N_i|z_i]=(d^0 + z_i)\sum \nolimits _{k=1}^{K}\sum \nolimits _{j=1}^{r_{ik}}\lambda _{jk}. \nonumber \\ \end{aligned}$$

(35)

The above shows the analysis on the number of interruptions for fixed duration, now we extend that with consideration of stochastic duration. According to Law of Total Expectation (Billingsley 1995), the expected value of $N_i$ is the same as the conditional expected value of $N_i$ given $\tilde{z}_i$, i.e.,

$$\begin{aligned} \begin{array}{l} E[N_i]= E_{\tilde{z}_i}[E_{N_i|\tilde{z}_i}[N_i|\tilde{z}_i]]. \end{array} \end{aligned}$$

(36)

Based on the representation of the mean number of interruptions for the deterministic case in Eq. 35 and Eq.  34, mean of $N_i$ can be calculated as:

$$\begin{aligned} E[N_i]= & {} E_{\tilde{z}_i}[E_{N_i|\tilde{z}_i}[N_i|\tilde{z}_i]]\nonumber \\= & {} E_{\tilde{z}_i}\left[ (d^0 + \tilde{z}_i)\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}\right] \nonumber \\= & {} d^0\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}. \end{aligned}$$

(37)

According to Law of Total Variance (Billingsley 1995), variance of the number of interruptions of activity $a_i$ can be represented as:

$$\begin{aligned} \begin{array}{l} Var[N_i]= E[Var[N_i|\tilde{z}_i]]+Var[E[N_i|\tilde{z}_i]]. \end{array} \end{aligned}$$

(38)

According to our analysis for the variance representation of the deterministic case in Eq. 35, we have,

$$\begin{aligned} E[Var[N_i|\tilde{z}_i]]= & {} E\left[ (d^0 + \tilde{z}_i)\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}\right] \nonumber \\= & {} d^0\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}, \end{aligned}$$

(39)

and

$$\begin{aligned} Var[E[N_i|\tilde{z}_i]]= & {} Var[(d^0 + \tilde{z}_i)\sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}]\nonumber \\= & {} \sigma _i^2\left( \sum _{k=1}^{K}\sum _{j=1}^{r_{ik}}\lambda _{jk}\right) ^2. \end{aligned}$$

(40)

Thus, the variance value of the number of interruptions can be represented by:

$$\begin{aligned} \begin{array}{l} Var[N_i]= d^0\sum \limits _{k=1}^{K}\sum \limits _{j=1}^{r_{ik}} \lambda _{jk}+\sigma _i^2\left( \sum \limits _{k=1}^{K} \sum \limits _{j=1}^{r_{ik}}\lambda _{jk}\right) ^2. \end{array} \end{aligned}$$

(41)

1.2 Total vacancy time

Lemma 9.3

(Lambrechts et al. 2011) Let X and Y be independent random variables that are both exponentially distributed, respectively with parameter $\lambda$ and $\mu$. The probability that X will be smaller than Y is then:

$$\begin{aligned} P(X < Y) = \lambda / (\lambda + \mu ). \end{aligned}$$

Let $P_{jk}$ denote the probability that the interruption for activity $a_i$ is caused by the resource unit j of resource type k. That is, among all resource units used by activity $a_i$, this resource unit has the smallest time to breakdown. Then $P_{jk}$ can be calculated as:

$$\begin{aligned} P_{jk} = P(X_{jk} < \min \limits _{l=1,\ldots K;i=1,\ldots r_{il};(i,l)\ne (j.k)} X_{il}). \end{aligned}$$

(42)

From Eq. 42 and Lemma 9.3, we can rewrite $P_{jk}$ as:

$$\begin{aligned} P_{jk}= & {} \lambda _{jk} / (\lambda _{jk} + \sum \limits _{l=1,\ldots K;i=1,\ldots r_{il};(i,l)\ne (j.k)}\lambda _{il})\nonumber \\= & {} \lambda _{jk} / \sum \limits _{l=1}^{K}\sum \limits _{i=1}^{r_{il}}\lambda _{il}. \end{aligned}$$

(43)

Since we assume that only one unit of resource is allowed to break down at a time, we can calculate the mean and variance values of $R_i$ describing the vacancy time for $a_i$ to wait once interrupted as:

$$\begin{aligned} E[R_i]= & {} \sum \limits _{k=1}^K\sum \limits _{j=1}^{r_{ik}}P_{jk}E[Y_{jk}].\nonumber \\ Var[R_i]= & {} E[R_i^2] - E^2[R_i]\nonumber \\= & {} \sum \limits _{k=1}^K\sum \limits _{j=1}^{r_{ik}}P_{jk}E[Y_{jk}^2] - E^2[R_i]\nonumber \\= & {} \sum \limits _{k=1}^K\sum \limits _{j=1}^{r_{ik}}P_{jk}\left( Var[Y_{jk}] + E^2[Y_{jk}]\right) - E^2[R_i].\nonumber \\ \end{aligned}$$

(44)

Here, random variable $Y_{jk}$ represents the time needed to repair the resource unit j of resource type k once broken. Since we assume that $Y_{jk}\sim EXP(\mu _{jk})$, then we can obtain:

$$\begin{aligned} E[Y_{jk}] = 1/\mu _{jk}, \quad Var[Y_{jk}] = 1/(\mu _{jk})^2. \end{aligned}$$

(45)

Therefore, from Eqs. 19, 20, 37, 41, 44 and 45, mean and variance values of the stochastic part of processing time can be derived as:

$$\begin{aligned}&E[\tilde{\delta _i}] = d_i^0\sum \limits _{k=1}^{K}\sum \limits _{j=1}^{r_{ik}} \frac{\lambda _{jk}}{\mu _{jk}},\nonumber \\&Var[\tilde{\delta _i}]= d_i^0\sum \limits _{k=1}^{K}\sum \limits _{j=1}^{r_{ik}} \frac{2\lambda _{jk}}{(\mu _{jk})^2}+\sigma _i^2\left( \sum \limits _{k=1}^{K} \sum \limits _{j=1}^{r_{ik}}\frac{\lambda _{jk}}{\mu _{jk}}\right) ^2.\nonumber \\ \end{aligned}$$

(46)