Estimating the Probability of Meeting a Deadline in Hierarchical Plans

1 Introduction

Numerous planning tools produce plans that call for executing tasks non-linearly. Usually, such plans are represented as a tree, where the leaves indicate primitive tasks, and other nodes represent compound tasks consisting of executing their sub-tasks either in parallel (also called "concurrent" tasks [Gabaldon, 2002]) or in sequence. [Erol et al., 1994; Nau et al., 1998; 2003; Kelly et al., 2008].

Given such a hierarchical plan representation, it is frequently of interest to evaluate its desirability in terms of resource consumption, such as fuel, cost, or time. The answer to such questions can be used to decide which of a set of plans, all valid as far as achieving the goal(s) are concerned, is better given a user-specified utility function. Another reason to compute these distributions is to support runtime monitoring of resources, generating alerts to the execution software or human operator if resource consumption in practice has a high probability of surpassing a given threshold.

While most tools aim at good average performance of the plan, in which case one may ignore the full distribution and consider only the expected resource consumption [Bonfietti et al., 2014], our paper focuses on providing guarantees for the probability of meeting deadlines. This type of analysis is needed, e.g., in Service-Level-Agreements (SLA) where guarantees of the form: "response time less than 1mSec in at least 95% of the cases" are common [Buyya et al., 2011] Section 8 discusses additional related work.

We assume that an hierarchical plan is given in the form of a tree, with uncertain resource consumption of the primitive actions in the network, provided as a probability distribution.

The problem is to compute a property of interest of the distribution for the entire task network. In this paper, we focus mainly on the issue of computing the probability P (t < T ) of satisfying a deadline T (i.e. that the makespan t of the plan is less than a given value). Since in the above-mentioned applications for these computations, one needs results in real-time (for monitoring) or multiple such computations (in comparing candidate plans), efficient computation here is crucial, and is more important than in, e.g., off-line planning.

We show that computing P(t < T) is NP-hard (see Section 6) even for the simple sum of independent random variables (r.v.s) , the first contribution of this paper. A deterministic polynomial-time approximation scheme for this problem is proposed, the second contribution of this paper. Error bounds are analyzed and are shown to be tight. For discrete r.v.s with finite support, finding the distribution of the maximum can be done in low-order polynomial time. However, when compounded with errors generated due to approximation in subtrees, handling this case requires careful analysis of the resulting error. The approximations developed for both sequence and parallel nodes are combined into an overall algorithm for task trees, with an analysis of the resulting error bounds, yielding a polynomial-time (additive error) approximation scheme for computing the probability of satisfying a deadline for the complete network, another contribution of this paper.

We also briefly consider computing expected makespan. Since for discrete r.v.s, in parallel nodes one can compute an exact distribution efficiently, it is easy to compute an expected makespan in this case as well as for sequence nodes. Despite that, we show that for trees with both parallel and sequence nodes, computing the expected makespan is NP-hard.

Experiments are provided in order to examine the quality of approximation in practice when compared to the theoretical error bounds. A simple sampling scheme is also provided as a yardstick, even though the sampling does not come with error guarantees, but only bounds in probability. Finally, we examine our results in light of related work in the fields of planning and scheduling, as well as probabilistic reasoning.

2 Problem statement

We are given a hierarchical plan represented as a task tree consisting of three types of nodes: primitive actions as leaves, sequence nodes, and parallel nodes. Primitive action nodes contain distributions over their resource consumption. Although other resources can be handled with our approach (e.g. memory: tasks running in parallel use the sum of the memory space of each of the tasks; if they run in sequence, only the maximum thereof is needed), we will assume henceforth, in order to be more concrete, that the only resource of interest is time. A sequence node vs denotes a task that has been decomposed into subtasks, represented by the children of vs, which must be executed in sequence in order to execute vs. We assume that a subtask of vs begins as soon as its predecessor in the sequence terminates. Task vs terminates when its last subtask terminates. A parallel node vp also denotes a decomposed task, but subtasks begin execution in parallel immediately when task vp begins execution; vp terminates as soon as all of the children of vp terminate.

Resource consumption is uncertain, and described as probability distributions in the leaf nodes. We assume that the distributions are independent (but not necessarily identically distributed). We also assume initially that the r.v.s are discrete and have finite support (i.e. number of values for which the probability is non-zero). As the resource of interest is assumed to be completion time, let each leaf node v have a completion-time distribution Pv, in some cases represented as a cumulative distribution function form (CDF) Fv.

The main computational problem tackled in this paper is: Given a task tree τ and a deadline T , compute the probability that τ satisfies the deadline T (i.e. terminate in time t ≤ T). We show that this problem is NP-hard and provide an approximation algorithm. The above deadline problem reflects a step utility function: a constant positive utility U for all t less than or equal to a deadline time T, and 0 for all t > T. We also briefly consider a linear utility function, requiring computation of the expected completion time of τ , and show that this expectation problem is also NP-hard.

3 Sequence nodes

Since the makespan of a sequence node is the sum of durations of its components, the deadline problem on sequence nodes entails computation of (part of) the CDF of a sum of r.v.s which is an NP-hard problem (shown in Section 6). Thus, there is a need for an approximation algorithm for sequence nodes. The main ingredient in our approximation scheme is the Trim operator specified as follows:

Definition 1 (The Trim operator). For a discrete r.v. X and a parameter ε > 0, consider the sequence of elements in the support of X defined recursively by: x1 = min support(X) and, if the set C = {x > xi : Pr(xi < X ≤ x) > ε} is not empty, xi+1 = min{x > xi : Pr(xi < X ≤ x) > ε}. Let l be the length of this sequence, i.e., let l be the first index for which C is empty. For notational convenience, define xl+1 = ∞. Now, define X′ = Trim(X, ε) to be the random variable specified by:

For example, if X is a r.v. such that Pr(X=1)=0.1, Pr(X=2)=0.1 and Pr(X=4)=0.8, the r.v. X′=Trim(X, 0.5) is given by Pr(X′=1)=0.2 and Pr(X′=4)=0.8. Intuitively, the Trim operator removes consecutive domain values whose accumulated probability is less than ε and adds their probability mass to the element in the support that precedes them.

Using the Trim operator, we are now ready to introduce the main operator of this section:

Sequence takes a set of r.v.s and computes a r.v. that represents an approximation of their sum by applying the Trim operator after adding each of the variables. The parameter ε, (see below) specifies the accuracy of approximation.

The Sequence operator can be implemented by the procedure outlined in Algorithm 1. The algorithm computes the distribution ∑ni=0 Xi using convolution (the Convolve() operator in line 3) in a straightforward manner. Computing the convolution is itself straightforward and not further discussed here. However, since the support of the resulting distribution may be exponential in the number of convolution operations, the algorithm must trim the distribution representation to avoid this exponential blow-up. This decreases the support size, while introducing error. The trick is to keep the support size under control, while making sure that the error does not increase beyond a desired tolerance. Note that the size of the support can also be decreased by simple "binning" schemes, but these do not provide the desired guarantees. In the algorithm, the PDF of a r.v. Xi is represented by the list DXi , which consists of (x, p′ ) pairs, where x ∈ support(Xi) and p′ is the probability Pr(Xi=x), the latter denoted by DXi [x]. We assume that DXi is kept sorted in increasing order of x.

Algorithm 1: Sequence (X1, . . . ,Xn , ε)

1 D = ((0, 1)) // Dummy random var.: 0 with prob. 1 2 for i = 1 to n do 3 D=Convolve(D,DXi ) 4 D=Trim(D, ε) 5 return D 6 Procedure Trim(D,ε) 7 D′ = (), p = 0 8 d0 = dprev =minsupport(D) 9 foreach d ∈ support(D) \ {d0} in ascending order do 10 if p + D[d] ≤ ε then 11 p = p + D[d] 12 else 13 Append (dprev, D[dprev] + p) to D′ 14 dprev = d, p = 0 15 Append (dprev,D[dprev] +p) to D′ 16 return D′

We proceed to show that Algorithm 1 indeed approximates the sum of the r.v.s, and analyze its accuracy/efficiency tradeoff. A notion of approximation relevant to deadlines is:

Note that this definition is asymmetric, because, as shown below, our algorithm never underestimates the exact value. For the proof of the Lemma 1 below, we establish the following technical claim (can be proven by induction on n):

Claim 1. Let (ai)ni=1 and (bi)ni=1 be sequences of real numbers such that ∑ki=1 ai ≥ 0 for all 1 ≤ k ≤ n and b1 ≥ b2 ≥ · · · ≥ bn ≥ 0, then ∑ni=1 aibi ≥ 0

We now bound the approximation error of sums of r.v.s:

Lemma 1. For discrete random variables X1,X′1 ,X2,X′2 and ε1, ε2 ∈ [0, 1], if X′1 ≈ε1 X1 and X′2 ≈ε2 X2, then X′1 +X′2 ≈ε1+ε2 X1 +X2.

Proof. Pr(X′1 +X′2 ≤ T) -Pr(X1 +X2 ≤ T)

Finally, we show that the difference is also nonnegative:

The first term here is non-negative by Claim 1, the second is nonnegative because it is a sum of nonnegative numbers.

We now show that Trim(X, ε) is an ε-approximation of X: Lemma 2. Trim(X, ε) ≈ε X

Proof. Let X′ = Trim(X, ε). Let x1<· · ·<xm be the support of X′ and let l = max{i : xi ≤ T}. We have,

because, after Trim, the probabilities of elements that were removed from the support are assigned to the element that precedes them. From Equation (1) we get:

The inequality Pr(T < X < xl+1) ≤ ε follows from the observation that, for all i, Pr(xi < X < xi+1) < ε, because p is never greater than ε in Algorithm 1.

To bound the amount of memory needed for our approximation algorithm, the next lemma bounds the size of the support of the trimmed r.v.:

Lemma 3. | support(Trim(X, ε))| ≤ 1/ε

Proof. Let X′ = Trim(X, ε) and let x1 < · · · < xm be the support of X′ . And, for notational convenience, let xm+1 = ∞. Let pi = ∑xi<x<xi+1 Pr(X=x). Then, 1 = ∑mi=1 Pr(X′=xi) = Pr(X=x1) + ∑m-1 i=1 (pi + Pr(X=xi+1)) + pm. According to algorithm 1, lines 11-12, pi ≤ ε and pi +Pr(X=xi+1) > ε for all 0 ≤ i < m. Therefore, 1 = ∑mi=1 P r(X′=xi) > P r(X′=x1) + (m - 1) · ε + pm. Using the fact Pr(X′=x1) > 0, we get: (m-1) · ε < 1, therefore m ≤ 1/ε.

These lemmas highlight the main idea behind our approximation algorithm: the Trim operator trades off approximation error for a reduced size of the support. The fact that this trade-off is linear allows us to get a linear approximation error in polynomial time, as shown below:

Theorem 1. If X′ i ≈εi Xi for all i ∈ {1, . . . , n} and ˆX = Sequence(X′1 , . . . ,X′n , ε) then ˆX ≈e ∑ni=1 Xi, where e = ∑ni=1 εi + nε.

Proof. (outline) For n iterations, from Lemma 1, we get an accumulated error of ε1 + · · · + εn. From Lemma 2, we get an additional error of at most nε due to trimming.

Theorem 2. Assuming that m ≤ 1/ε, the procedure Sequence(X′1 , . . . ,X′n , ε) can be computed in time O((nm/ε) log(m/ε)) using O(m/ε) memory, where m is the size of the largest support of any of the X′ i s.

Proof. From Lemma 3, the size of list D in Algorithm 1 is at most m/ε just after the convolution, after which it is trimmed, so the space complexity is O(m/ε). Convolve thus takes time O((m/ε) log(m/ε)), where the logarithmic factor is required internally for sorting. Since the runtime of the Trim operator is linear, and the outer loop iterates n times, the overall run-time of the algorithm is O((nm/ε) log(m/ε)).

Example 1. The error bound provided in Theorem 1 is tight, i.e. Sequence(X1, . . . ,Xn, ε/n) may results in error ε: Let 0≤ε<1 and n∈N such that 1-ε>ε/n. Consider, for very small δ > 0, the r.v. X1 defined by:

and, for i = 2, . . . , n, let the r.v.s Xi be such that Pr(Xi=0)=1-δ, Pr(Xi=n2 )=δ, and zero otherwise.

4 Parallel nodes

Unlike sequence composition, the deadline problem for parallel composition is easy to compute, since the execution time of a parallel composition is the maximum of the durations:

where the last equality follows from independence of the r.v.s. We denote the construction of the CDF using Equation (2) by Parallel(X1, . . . ,Xn). If the r.v.s are all discrete with finite support, Parallel(X1, . . . ,Xn) incurs linear space, and computation time O(nmlog(n)).

If the task tree consists only of parallel nodes, one can compute the exact CDF, with the same overall runtime. However, when the task tree contain both sequence and parallel nodes we may get only approximate CDFs as input, and now the above straightforward computation can compound the errors. When the input CDFs are themselves approximations, we bound the resulting error:

Lemma 4. For discrete r.v.s X′1 , . . .X′n , X1, . . . ,Xn, if for all i = 1, . . . , n, X′ i ≈εi Xi and 0 ≤ εi ≤ 1 n(Kn+1) for some K > 0, then, for any ε ≥ εi, we have: maxi∈[1:n] X′ i ≈e maxi∈[1:n] Xi where e = ∑ni=1 εi + ε/K.

Proof. Pr(maxi∈[1:n] X′ i ≤ T )-P r(maxi∈[1:n] Xi ≤ T )

Since P r(X′ i ≤ T) > Pr(Xi ≤ T) for each i, this expression is nonnegative.

5 Task trees: mixed sequence/parallel

Given a task tree τ and a accuracy requirement 0 < ε < 1, we generate a distribution for a r.v. X′τ approximating the true duration distribution Xτ for the task tree. We introduce the algorithm and prove that the algorithm indeed returns an ε-approximation of the completion time of the plan. For a node v, let τv be the sub tree with v as root and let childv be the set of children of v. We use the notation |τ| to denote the total number of nodes in τ .

Algorithm 2, that implements the operator Network, is a straightforward postorder traversal of the task tree. The only remaining issue is handling the error, in an amortized approach, as seen in the proof of the following theorem.

Theorem 3. Given a task tree τ , let Xτ be a r.v. representing the true distribution of the completion time for the network. Then Network(τ, ε) ≈ε Xτ .

Algorithm 2: Network(τ , ε)

1 Let v be the root of τ // Hence, τv = τ 2 nv = |childv| 3 if v is a Primitive node then 4 return the distribution of v 5 if v is a Sequence node then 6 for c ∈ childv do 7 Xc = Network (τc, |τc|ε |τv| ) 8 return Sequence ({Xc}c∈childv , ε nv |τv | ) 9 if v is a Parallel node then 10 for c ∈ childv do 11 Xc = Network (τc, min( |τc|ε |τv| , 1 nv (|τv |nv+1) )) 12 return Parallel ({Xc}c∈childv )

Proof. By induction on |τ|. Base: |τ| = 1, the node must be primitive, and Network will just return the distribution unchanged which is obviously an ε-approximation of itself. Suppose the claim is true for 1 ≤ |τ| < n. Let τ be a task tree of size n and let v be the root of τ. If v is a Sequence node, by the induction hypothesis that Xc ≈|τc|ε/|τv| Xτc , and by Theorem 1, the maximum accumulated error is ∑c∈childv |τc|ε/|τv| + ε/|τv| = (n - 1)ε/|τv| + ε/|τv| = ε for v, therefore, Sequence({Xc}c∈childv , ε/n) ≈ε Xτ as required. If v is a Parallel node, by the induction hypothesis that Xc ≈ec Xτc , where ec = min( |τc|ε |τv| , 1 nv (|τv |nv+1) ) so ∑c∈childv ec ≤ ∑c∈childv |τc|ε |τv| ≤ ε - ε/|τv|. Then, by Lemma 4, using K = |τv| and n = nv, we get that Parallel({Xc}c∈childv ) ≈ε Xτ as required.

Theorem 4. Let N be the size of the task tree τ, and M the size of the maximal support of each of the primitive tasks. If 0 ≤ ε ≤ 1 N(N 2+1) and M < N/ε, the Network approximation algorithm runs in time O((N 5 /ε2 ) log(N 3/ε2 )), using O(N 3 /ε2 ) memory.

Proof. The run-time and space bounds can be derived from the bounds on Sequence and on Parallel, as follows. In the Network algorithm, the trim accuracy parameter is less than or equal to ε/N. The support size (called m in Theorem 2) of the variables input to Sequence are O(N 2 /ε). Therefore, the complexity of the Sequence algorithm is O((N 4 /ε2 ) log(N 3/ε2 )) and the complexity of the Parallel operator is O((N 3 /ε) log(N)). The time and space for sequence dominate, so the total time complexity is N times the complexity of Sequence and the space complexity is that of Sequence.

If the constraining assumptions on M and ε in Theorem 4 are lifted, the complexity is still polynomial: replace one instance of 1/ε by max(m, 1/ε), and the other by max(1/ε, N(N 2 +1)) in the runtime complexity expression.

6 Complexity results

We show that the deadline problem is NP-hard, even for a task tree consisting only of primitive tasks and one sequence node, i.e. linear plans.

Lemma 5. Let Y = {Y1, . . . , Yn} be a set of discrete realvalued r.v.s specified by probability mass functions with finite supports, T ∈ Z, and p ∈ [0, 1]. Then, deciding whether P r(∑ni=0 Yi < T ) > p is NP-Hard.

Proof. By reduction from SubsetSum [Garey and Johnson, 1990, problem number SP13]. Recall that SubsetSum is: given a set S = {s1, . . . , sn} of integers, and integer target value T, is there a subset of S whose sum is exactly T? Given an instance of SubsetSum, create the two-valued r.v.s Y1, . . . , Yn with Pr(Yi = si) = 1/2 and Pr(Yi = 0) = 1/2. By construction, there exists a subset of S summing to T if and only if P r(∑ni=0 Yi = T ) > 0.

Suppose that algorithm A(Y, T, p) can decide Pr(∑ni=0 Yi < T) > p in polynomial time. Then, since the r.v.s Yi are two-valued uniform r.v.s, the only possible values of p are integer multiples of 1/2n , and we can compute p = Pr(∑ni=0 Yi < T) using a binary search on p using n calls to A. To determine whether Pr(∑ni=0 Yi = T) > 0, simply use this scheme twice, since Pr(∑ni=0 Yi = T) > 0 is true if and only if P r(∑ni=0 Yi < T ) < P r(∑ni=0 Yi < T + 1).

Theorem 5. Finding the probability that a task tree satisfies a deadline T is NP-hard.

Proof. Given a task tree consisting of leaf nodes, all being children of a single sequence node, its makespan is the sum of the completion times of the leaves. The theorem follows immediately from Lemma 5.

Finally, we consider the linear utility function, i.e. the problem of computing an expected makespan of a task network. Note that although for linear plans the deadline problem is NP-hard, the expectation problem is trivial because the expectation of the sum of r.v.s Xi is equal to the sum of the expectations of the Xis. For parallel nodes, it is easy to compute the CDF and therefore also easy to compute the expected value. Despite that, for task networks consisting of both sequence nodes and parallel nodes, these methods cannot be effectively combined, and in fact, we have:

Theorem 6. Computing the expected completion time of a task network is NP-hard.

Proof. By reduction from subset sum. Construct r.v.s ("primitive tasks") Yi as in the proof of Lemma 5, and denote by X the r.v. ∑ni=1 Yi. Construct one parallel node with two children, one being the a sequence node having the completion time distribution defined by X, the other being a primitive task that has a completion time Tj with probability 1. (We will use more than one such case, which differ only in the value of Tj, hence the subscript j). Denote by Mj the r.v. that represents the completion time distribution of the parallel node, using this construction, with the respective Tj . Now consider computing the expectation of theMj for the following cases: T1 = T+1/2 and T2 = T+1/4. Thus we have, for j ∈ {1, 2}, by construction and the definition of expectation:

where the second equality follows from the Yi all being integer-valued r.v.s (and therefore X is also integer valued). Subtracting these expectations, we have E[M1] - E[M2] = 14 P(X ≤ T). Therefore, using the computed expected values, we can compute P(X ≤ T), and thus also P(X = T), in polynomial time.

7 Empirical Evaluation

We examine our approximation bounds in practice, and compare the results to exact computation of the CDF and to a simple stochastic sampling scheme. Three types of task trees are used in this evaluation: task trees used as execution plans for the ROBIL team entry in the DARPA robotics challenge (DRC simulation phase, http://in.bgu.ac.il/en/Pages/news/dar pa.aspx), linear plans (seq), and plans for the Logistics domain (from IPC2 http://ipc.icaps-conference.org/). The primitive task distributions were uniform distributions discretized to M values. For every entry of M in Table 1 the first line is the runtime in seconds, the second line presents the estimation error.

In the Logistics domain, packages are to be transported by trucks or airplanes. Hierarchical plans were generated by the JSHOP2 planner [Nau et al., 2003] for this domain and consisted of one parallel node (packages delivered in parallel), with children all being sequential plans. The duration distribution of all primitive tasks is uniform but the support parameters were determined by the type of the task, in some tasks the distribution is fixed (such as for load and unload) and in others the distribution depends on the velocity of the vehicle and on the distance to be travelled.

After running our approximation algorithm we also ran a variant that uses an inverted version of the Trim operator, providing a lower bound of the CDF, as well as the upper bound generated by Algorithm 2. Running both variants allows us to bound the actual error, costing only a doubling of the run-time. Despite the fact that our error bound is theoretically tight, in practice and with actual distributions, according to Table 1, the resulting error in the algorithm was usually much better than the theoretical ε bound.

We ran the exact algorithm, our approximation algorithm with ε ∈ {0.1, 0.01, 0.001}, and a simple simulation with 103 to 107 samples (number of samples is denoted by s in the table), on networks from the DRC implementation, sequence nodes with 10, 20, and 50 children (number of nodes denoted by N in the table), and 20 Logistics domain plans, and several values of M (M, N are as in Theorem 4). Results for a typical indicative subset (regretfully reduced due to page limits) are shown in table 1.

The exact algorithm times out in some cases. Both our approximation algorithm and the sampling algorithm handle all these cases, as our algorithm's runtime is polynomial in N, M, and 1/ε as is the sampling algorithm's (time linear in number of samples).

The advantage of the approximation algorithm is mainly in providing bounds with certainty as opposed to the bounds inprobability provided by sampling. Additionally, as predicted by theory, accuracy of the approximation algorithm improves linearly with 1/ε (and almost linear in runtime), whereas accuracy of sampling improves only as a square root of the number of samples. Thus, even in cases where sampling initially outperformed the approximation algorithm, increasing the required accuracy for both algorithms, eventually the approximation algorithm overtook the sampling algorithm.

8 Discussion and RelatedWork

Numerous issues remain unresolved, briefly discussed below. Trivial improvements to the Trim operator are possible, such as the inverse version of the operator used to generate a lower bound for the empirical results. Other candidate improvements are not performing trimming (or even stopping a trimming operation) if the current support size is below 1/ε, which may increase accuracy but also the runtime. Another point is that in the combined algorithm, space and time complexity can be reduced by adding some Trim operations, especially after processing a parallel node, which is not done in our version. This may reduce accuracy, a trade-off yet to be examined. Another option is, when given a specific threshold, trying for higher accuracy in just the region of the threshold, but how to do that is non-trivial. For sampling schemes such methods are known, including adaptive sampling [Bucher, 1988; Lipton et al., 1990], stratified sampling, and other schemes. It may be possible to apply such schemes to deterministic algorithms as well - an interesting issue for future work.

Extension to continuous distributions: our algorithm can handle them by pre-running a version of the Trim operator on the primitive task distribution. Since one cannot iterate over support values in a continuous distribution, start with the smallest support value (even if it is -∞), and find the value at which the CDF increases by ε. This requires access to the inverse of the CDF, which is available, either exactly or approximately, for many types of distributions.

We showed that the expectation problem is also NP-hard. A natural question is on approximation algorithms for the expectation problem, but the answer here is not so obvious. Sampling algorithms may run into trouble if the target distribution contains major outliers, i.e. values very far from other values but with extremely low probability. Our approximation algorithm can also be used as-is to estimate the CDF and then to approximate the expectation, but we do not expect it to perform well because our current Trim operator only limits the amount of probability mass moved at each location to ε, but does not limit the "distance" over which it is moved. The latter may be arbitrarily bad for estimating the expectation. Possibly adding simple binning schemes to the Trim operator in addition to limiting the moved probability mass to ε may work, another issue for future research.

Related work on computing makespan distributions includes [Hong, 2013], which examines sum of Bernoulli distributed r.v.s. Other work examines both deterministic [Mercier, 2007] and Monte-Carlo techniques [Bucher, 1988; Lipton et al., 1990]. Distribution of maximum of r.v.s was studied in [Devroye, 1980], with a focus mostly on continuous distributions.

Complexity of finding the probability that the makespan is under a given threshold in task networks was shown to be NP-hard in [Hagstrom, 1988], even when the completion time of each task has a Bernoulli distribution. Nevertheless, our results are orthogonal as the source of the complexity in [Hagstrom, 1988] is in the graph structure, whereas in our setting the complexity is due to the size of the support. In fact for linear plans (an NP-hard case in our setting), the probability of meeting the deadline can be computed in low-order polynomial time for Bernoulli distributions, using straightforward dynamic programming. Makespan distributions in series parallel networks in the i.i.d. case was examined in [Gutjahr and Pflug, 1992], without considering algorithmic issues. There is also a significant body of work on estimating the makespan of plans and schedules [Herroelen and Leus, 2005; Fu et al., 2010; Beck and Wilson, 2007], within a context of a planner or scheduler. The analysis in these papers is based on averaging or on limit theorems, and does not provide a guaranteed approximation scheme. Temporal planing and in particular TPNs (temporal plan network) are presented in [Kim et al., 2001], the model is similar to ours, but the focus is on lower/upper bounds, rather than probability distributions. Hierarchical constraint-based plans in MAPGEN [Ai-Chang et al., 2004] allow for more general dependencies than seriesparallel, providing additional expressive power but making the deadline problem even harder.

Computing the distribution of the makespan in trees is a seemingly trivial problem in probabilistic reasoning [Pearl, 1988]. Given the task network, it is straightforward to represent the distribution using a Bayes network (BN) that has one node per task, and where the children of a node v in the task network are represented by BN nodes that are parents of the BN node representing v. This results in a treeshaped BN, where it is well known that probabilistic reasoning can be done in time linear in the number of nodes, e.g. by belief propagation (message passing) [Pearl, 1988; Kim and Pearl, 1983]. The difficulty is in the potentially exponential size of variable domains, which our algorithm, essentially a limited form of approximate belief propagation, avoids by trimming.

Looking at makespan distribution computation as probabilistic reasoning leads to interesting issues for future research, such as how to handle task completion times that have dependencies, represented as a BN. Since reasoning in BNs is NP-hard even for binary-valued variables [Dagum and Luby, 1993; Cooper, 1990], this is hard in general. But for cases where the BN toplogy is tractable, such as for BNs with bounded treewidth [Bodlaender, 2006], or directed-path singly connected BNs [Shimony and Domshlak, 2003], a deterministic polynomial-time approximation scheme for the makespan distribution may be achievable. The research literature contains numerous randomized approximation schemes that handle depenencies [Pearl, 1988; Yuan and Druzdzel, 2006], especially for the case with no evidence. In fact, our implementation of the sampling scheme in ROBIL handled dependent durations. It is unclear whether such sampling schemes can be adapted to handle dependencies and arbitrary evidence, such as: "the completion time of compound task X in the network is known to be exactly 1 hour from now". Finally, one might consider additional commonly used utility functions, such as a "soft" deadline: the utility is a constant U before the deadline T, decreasing linearly to 0 until T + G for some "grace" duration G, and 0 thereafter.

Acknowledgments. This research was supported by the ROBIL project, by the EU, by the ISF, and by the Lynne and William Frankel Center for Computer Science.