This paper focuses on estimating maximum resource loads in ResourceConstrained Project Scheduling Problem (RCPSP). We examine a variant of vector sum problem with fractions: considering preemptions allowed, determine what part of each of n jobs should be accomplished in order to minimize sum of squares of non-consumed amounts of resources, taking into account the resource constraints along with minimizing the number of preemptions. We prove that in case of 2 resources, the optimal solution contains only 2 or less preemptions, and present two polynomial algorithms of finding such solution with time complexities of O(n2) operations. We also present an investigation of the general case of arbitrary number of resources.
Resource-Constrained Project Scheduling Problem (RCPSP) is one of the core problems of scheduling theory. It has been widely explored throughout the last years. RCPSP is NP-hard in the strong sense [Garey et al., 1975], which means no exact approach to it may find a solution in polynomial time. The exact approaches to RCPSP include, but are not limited to branch-and-bound methods, linear programming, constraint programming and dynamic programming methods. For all of these methods, knowing the lower makespan bounds allows to vastly reduce search time. In [Arkhipov et al., 2017], a polynomial algorithm of evaluating lower makespan bounds is presented. The algorithm is based on estimating maximum resource load for a given pair of resources, which is the main goal of the present paper. 2
The problem of estimating maximum resource load can be formulated in the following way. There is a set of jobs (tasks) N and a set of non-renewable resources R, and available amounts rj ∈ R+ of resources j ∈ R. Each task may be either completed, rejected or executed partially. In order to complete task i ∈ N , amount rij ∈ R+ of each resource j ∈ R is required.
The completion ratio xi defines the executed part of each task i ∈ N . Completion ratio value xi = 1 corresponds to full completion of task i ∈ N , xi = 0 correspond to rejection, and all the other values in (0, 1) correspond to partial execution (preemption). The required resources are consumed in proportion to completion ratio of this job xi ∈ [0, 1], i.e. amount rij xi of each resource j ∈ R is needed.
The goal is to determine the vector X = {xi} of completion ratios such that the aggregate amount of consumed resources is maximized ∑j∈R ∑i∈N xirij → max and the resource constraints ∑i∈N xirij ≤ rj ∀j ∈ R are satisfied. 3
In this section we will examine the case of two resources. As will be shown later, this problem can be solved in polynomial time. Let us denote the available amounts of resources as A, B and resource consumption amounts as ai, bi. We will also associate each task i ∈ N with a vector s⃗i = (ai, bi) in a two-dimensional rectangular Cartesian coordinate system. For the sake of convenience, here and below we will consider N = {1, . . . , n}, where n is the total number of tasks. 3.1
Let us reformulate the problem in a formal manner: given a pair of positive real numbers A, B ∈ R+ and a list of tasks described by a set N = {1, . . . , n} of indices and set S = {s⃗i} ⊂ R+2, i ∈ N of vectors s⃗i = (ai, bi) on a two-dimensional rectangular Cartesian coordinate system, find vector X⃗ = {xi} ∈ Rn subject to constraints A − ∑i∈N aixi ≥ 0, B − ∑i∈N bixi ≥ 0, xi ∈ [0, 1] with objective function ∑i∈N aixi + ∑i∈N bixi → max.
Consider the classical ”0-1” knapsack problem, which is known to be N P -hard, and can be solved in pseudopolynomial time [Pisinger, 1997]. If boolean decision variables are replaced with real decision variables (knapsack with fractions), the problem becomes solvable in polynomial time, as was shown by Dantzig [Dantzig, 1957]. The problem of estimating resource loads can also be regarded as a variant of subset sum (or vector sum) problem, which is essentially a special case of the knapsack problem such that the cost of any item is equal to the weight of this item. Therefore, one-dimensional subset sum problem with fractions is also solvable in polynomial time. In [Baburin et al., 2007], a variant of vector sum problem is considered, and polynomial algorithms are presented for the case of arbitrary dimensionality. However, the problem considered in [Baburin et al., 2007] differs from the one considered in our article in that it does not involve any resource constraints, and the decision variables are boolean. All the problems mentioned above are related to the vector partitioning problem, which has been covered in [Onn et al., 1999]. Thus, the formulation of estimating resource loads problem takes the form of two-dimensional variant of subset sum problem with fractions. One can assume that since the one-dimensional subset sum problem with fractions is solvable in polynomial time, a polynomial approach to the two-dimensional variant should exist as well, on which we will elaborate further.
As was mentioned before, the objective of the problem is to find a list of completion ratios that maximizes resource consumption, taking the resource constraints into account. Obviously, if the amounts of available resources exceed or are equal to their maximum possible consumption (∑in=1 ai ≤ A, ∑in=1 bi ≤ B), the solution would be X = {xi}|∀ i xi = 1, i.e. each task should be fully completed. However, this is rarely the case. 3.2
Consider two permutations γup and γdown of vectors si acquired by sorting these vectors according to non-increase and non-decrease of the ratio abii . Let us note the order of vectors according to these permutations as si up and ∀i ∈ 1, n − 1 bi up bi+up1 ai up ≥ ai+up1 For γdown the sign of this inequality is reversed. The division operation here is rather pro forma: obviously, the case of dividing by zero must be handled with a separate conditional operator. If, for any two vectors, ratios of their components ( bi ) happen to be equal, the order of these vectors in permutations γup and γdown for now is ai considered to be of no importance.
Note: obviously, sorting changes order of the vectors, which means the original indexation like si, xi, i ∈ N would not make any sense. However, we will continue to use this style of indexation just for convenience’s sake. When implementing the algorithms, it is useful to keep the original index stored as an object’s (vector’s) property.
Let us place vector s1up starting at the origin of our coordinate system, then place the second vector s2up starting at the end of first one, and so forth. We will name the resultant polygonal line as γup, too. Let us also do the same operations with elements of γdown, resulting in polygonal line γdown. Clearly, any permutation of n vectors corresponds to some polygonal curve built with vectors s⃗i, and, if provided a set of coefficients X⃗ = {xi}, together this permutation and this set of coefficients correspond to a polygonal curve built with vectors (xis⃗i).
Since in γup, γdown vectors are sorted according to non-increase and non-decrease of tangent of the angle between each vector and the abscissa axis of our coordinate system, these permutations also correspond to sorting according to the value of this angle. Thus, the polygonal line γup is convex upward, and γdown is convex downward, and together these lines form a convex polygon on the coordinate plane. We will name this polygon the solvability domain of considered instance of the problem.
Lemma: Any polygonal curve γ′ comprised of vectors si ∈ S is located completely within the solvability domain of considered instance of the problem.
Proof : We will view the considered polygonal curves as plots of piece-wise linear functions. For example the curve γ corresponds to function b (a). Of course, if some of the vectors has a zero A component (ai = 0), the resulting curve cannot be viewed as a function of A. However, these cases too can be handled by additional condition operators.
Let γ be such a polygonal curve that ∀γ′, ∀a ∈ [0, ∑n i=1 ai] b (a) ≥ b ′ (a). Then γ is identical to γup (neglecting the cases in which ∃ i, j : abi = abj ), i.e. it is convex upward in each of its breaking points. To prove i j this, suppose the contrary: let γ be convex downward in its j-th breaking point, i.e. ∃ j ∈ 1, n − 1 : aj < abjj++11 . b j Consider a polygonal line ψ that differs from γ in that its j-th and (j + 1)-th vectors are placed in reverse order: sj = sj+1, sj+1 = sj , ∀i ∈ N \ {j, j + 1} si = si . Then ψ is convex upward in its j-th breaking point, and ∀ a ∈ (∑ij=−11 ai , ∑ij=+11 ai ) b (a) < b (a), which contradicts the supposition that ∀γ′, ∀a ∈ [0, ∑in=1 ai] b (a) ≥ b ′ (a). Thereby γ is identical to γup. Similarly, by considering another polygonal curve γ, this time with the condition ∀γ′, ∀a ∈ [0, ∑n i=1 ai] b down (a) ≤ b ′ (a). Thereby ∀ γ′, ∀a ∈ [0, ∑n i=1 ai] b (a) ≤ b ′ (a), we can prove that ∀ a ∈ [0, ∑n
i=1 ai] b up (a) ≥ b ′ (a) ≥ b down (a).
Corollary: Obviously, for any polygonal curve γ′ composed of vectors xis⃗i such that X⃗ = {xi} does not violate the constraints xi ∈ [0, 1]
∀ γ′, ∀X⃗ | xi ∈ [0, 1] ∀i ∈ N a similar inequality holds true
n ∀a ∈ [0, ∑ aixi] b up (a) ≥ b ′ (a) ≥ b down (a)
i=1 i.e. such polygonal curves are also located within the solvability domain, which leads to the following: ∀ X = {xi}| xi ∈ [0, 1] ∀ i ∈ N point (∑n i=1 bixi) is located within the solvability domain.
i=1 aixi, ∑n Note: X should not necessarily be a solution to regarded instance of the problem.
This allows us to formulate the necessary and sufficient condition of possibility to consume all the available resources: full consumption of the available resources A and B is possible if and only if point (A, B) is located within the solvability domain. To verify the latter, consider the boundary of the solvability domain, which is a pair polygonal curves γup and γdown, to be defined by a pair of piece-wise linear functions b up (a), b down (a) and check if b up (A) ≥ B ≥ b down (A). For an empty task set S = ∅ the problem is considered solvable only if A = 0, B = 0. 3.3
In this section we’ll provide two algorithms of finding the solution that leads to full consumption of available resources; the time complexities are O(n2).
First we’ll examine the problem of full resource consumption, which is to check if ∃ X = {x1, . . . , xn} : A − ∑in=1 aixi + B − ∑in=1 bixi = 0, taking into consideration all the constraints mentioned before. The case in which full resource consumption is not possible will be examined after both of the algorithms (see Auxiliary Algorithm). 3.3.1
The first algorithm is aiming at consequently reducing the current instance of full resource consumption problem to a solvable subproblem that has a smaller number of tasks. This continues until the cardinality of current set of tasks S reaches either 2, 1 or 0. Then the remaining problem is identical to solving a simple system of linear equations.
Consider a fixed set of tasks N = {1, . . . , n} associated with vector set S = {s1, . . . , sn} and a resource vector (A, B). Suppose this instance of the full resource consumption problem has a solution X⃗ = {x1, . . . , xn}.
Notice that if ∃ i ∈ 1, n : xi = 0, then the solution X⃗ of this problem is identical to the solution X⃗ ′ of a subproblem for the task set S′ = S \ si and resource vector (A, B), except that in X⃗ ′ the component xi is skipped: ∀ j < i x′j = xj , ∀ j > i x′j = xj+1. Similar applies if ∃ i ∈ 1, n : xi = 1, then the solution X⃗ of this problem is identical to the solution X⃗ ′ except that the resource vector in the subproblem would be (A−ai, B −bi).
Let us introduce a boolean function u(s) : S → {0, 1}. u(s) = 1 means vector task s has been marked as used. The first algorithm is as follows.
Step 1. Sort the tasks of the original set S by non-increase of ratio abii . If ai = 0, place the corresponding vector at the end of the list.
Step 2. Check if full resource consumption is possible (check if the point (A, B) lies within the solvability domain). If not – interrupt the algorithm.
Step 3. Mark all the task vectors si ∈ S as not used: ∀ si ∈ S u(si) := 0.
Cycle 1: until there is an unused task (∃ si ∈ S : | u(si) = 0): 1) Find an unused task vector si ∈ S | u(si) = 0.
Mark si as used: u(si) := 1. Form a subset of tasks Si′ by removing task si from S: Si′ = S \ si. 2) Check if full resource consumption is possible for Si′ and (A, B). If it is, assign xi := 0, S := S′. i
Step 4. Cycle 2: until there is an unused task (∃ si ∈ S : u(si) = 0): 1) Select an unused task vector si ∈ S | u(si) = 0.
Mark si as used: u(si) := 1. Form a subset of tasks Si′ = S \ si and a modified resource vector (A′, B′) where A′ = A − ai, B′ = B − bi. 2) Check full resource consumption is possible for Si′ and (A′, B′). If it is, assign xi := 1, S := S′. i
Step 5. If set S consists of two tasks j, k – solve the following system { xj aj + xkak = A xj bj + xkbk = B If S consists of only one task j, solve any of the equations of the corresponding system. If S is empty, no action is needed.
End.
Proposition: Time complexity of Algorithm 1 is at most O(n2), where n is the number of tasks in the initial task set.
Proof : Examining each step of the Algorithm 1 individually, one can make sure that the most time-consuming steps are Step 3 and Step 4, at which the solvability of a selected subproblem is tested, taking O(n) operations for each subproblem instance. In total at most n subproblems are examined during both steps, taking O(n2) or less operations. Thus, the time complexity of Algorithm 1 is O(n2). 3.3.2
The second algorithm is aimed at finding a pair of indices j, k such that solution α, β of the system k−1 ∑ai + βak = A + i=1
k−1 ∑bi + βbk = B + i=1 j−1 ∑ai + (1 − α)aj i=1 j−1 ∑bi + (1 − α)bj i=1 αaj + αbj + k−1 ∑ ai + βak = A i=j+1 k−1 ∑ bi + βbk = B i=j+1 satisfies the constraints 0 ≤ α, β ≤ 1 (the indexation is provided as if the vectors si were already sorted). Then, the solution of the whole problem is constructed from this pair of indices j, k and coefficients α, β.
Let us modify the system:
Examining each pair j, k, which would take at least O(n2) time, is not necessary: note that for each fixed j, the coefficients β form a monotoneous non-increasing sequence depending on k (we’ll provide a way of finding such pair of α, β even if (1) has more than one solution).
Suppose that j ∈ N is fixed, and k1, k2 ∈ N are such that k1 < k2. Let α1, β1 be the solution of (1) for the pair of indices j, k1, and α2, β2 the solution of (1) for j, k2. Then, since all the variables and coefficients in (1) are non-negative, β1 ≥ β2. Therefore, a binary search can be implemented to slightly speed up the algorithm. The algorithm is constructed as follows.
Step 1. Sort the tasks of set S by non-increase of ratio abii . If ai = 0, place the corresponding vector at the end of the list.
Step 2. Check if the problem of full resource consumption is solvable (check if the point (A, B) lies within the solvability domain). If the problem is unsolvable – interrupt the algorithm.
For each j, k ∈ N calculate and memorize values of sums k−1 ∑ ai, i=j+1 k−1 ∑ bi. i=j+1 Step 4. Set j = 1. Cycle:
For a fixed j ∈ N , try to find (via binary search) the index k such that the solution α, β of (1) meets the constraint 0 ≤ β ≤ 1. When solving (1), use the values of sums calculated previously.
It may happen that for some pairs j, k system (1) has more than one solution, i.e. its determinant turns out to be 0. In such cases, we simply assign one of the variables α, β equal to 0 and solve the remaining equation to find the other variable.
If the obtained values α, β fit within the constraints 0 ≤ α, β ≤ 1, then interrupt the cycle. Otherwise, set j := j + 1 and repeat the cycle.
Note: the situation in which none of the constraints were met and j reaches n is impossible, since it would mean that full resource consumption is not possible and the algorithm must have been interrupted at Step 2. End. set.
Proposition: Time complexity of Algorithm 2 is O(n2), where n is the number of tasks in the initial task Proof : Examining each step of the Algorithm 2 individually, one can make sure that the most time-consuming step is Step 3, which can be done in O(n2). Thus, the time complexity of Algorithm 2 is O(n2). 3.3.3
To account for cases in which full resource consumption is not possible (e.g., ∑i = 1naixi < A), another algorithm is needed.
Step 1. Sort vectors si according to non-decrease of ratio abii . If ai = 0, place the corresponding vector at the end of the list.
Step 2. Assign a := 0, b := 0, i := 1 Step 3. Cycle: while (a + ai ≤ A)&(b + bi ≤ B)&(i < n) do (xi := 1, a := a + ai, b := b + bi, i := i + 1). Step 4. If (i < n), assign xi := min ( Bb−i b , Aa−ia ), i := i + 1 (only the minimum of these two values will meet the constraint 0 ≤ xi ≤ 1).
Cycle: while (i ≤ n) do (xi := 0, i := i + 1).
End.
Proposition: Time complexity of Auxiliary Algorithm is at most O(n log n), where n is the number of tasks in the initial task set.
Proof : Examining each step of the Auxiliary Algorithm individually, one can make sure that the most timeconsuming step is Step 1, at which sorting the tasks takes O(n log n) time using Hoare’s quicksort. Thus, the time complexity of Auxiliary Algorithm is O(n log n). 3.3.4
The algorithm to solve the problem of estimating maximum resource load for two resources (the main algorithm) consists of either Algorithm 1 or Algorithm 2 supplemented with the Auxiliary Algorithm. Overall time complexity, therefore, can be either O(n2) in case ∑n i=1 bixi < B or O(n log n) otherwise.
i=1 aixi < A, ∑n 4
The problem is close to a linear relaxation of multi-dimensional knapsack problem, with resource constraints being regarded as knapsack capacity constraints. The multi-dimensional knapsack problem has been widely described in [Kellerer et al., 2004].
We’ll call complete solution such a solution of the problem that allows for full resource consumption. We’ll provide a polynomial algorithm that allows to find a complete solution (if such exists and is unique) for the case r ≥ n.
Let us denote the matrix of required resource amounts rij as R, vector column of sought-for completion ratios xi as X, and the vector column of available resource amounts ai as A. Finding a complete solution of the original problem is then equivalent to solving a system of linear algebraic equations RX = A subject to constraints 0 ≤ xi ≤ 1, i = 1, . . . , n.
We’ll describe the algorithm first and then discuss the ideas behind it.
Algorithm for the case r ≥ n Step 1. Check if r ≥ n. If not – interrupt the algorithm.
Use the Gauss-Jordan elimination method on the augmented matrix (R|A) to obtain the reduced row echelon form R∗ of matrix R and the transformed right-hand side A∗.
Check if rank R∗ = rank (R∗|A∗). If not – interrupt the algorithm.
Step 2. Check if rank R = rank R∗ = n. If not, interrupt the algorithm.
Since these matrices are if reduced row echelon form, to calculate ranks of these matrices one needs simply to find the numbers of non-zero rows.
For all i = 1, . . . , n assign xi := ai∗.
Step 3. Check if all the completion ratios obtained xi fit within the constraints 0 ≤ xi ≤ 1. If not – interrupt the algorithm.
End.
Since Gauss-Jordan method uses only elementary row operations, the rank of the matrix is not affected by it: rank R = rank R∗. Steps 1 and 2 are based upon the Rouche-Capelli (Kronecker-Capelli) theorem. According to it, the system RX = A is consistent and has a unique solution if and only if rank R = rank (R|A) = n, which is equivalent to rank R∗ = rank (R∗|A∗) = n. This condition also implies that r ≥ n, which is best to be ascertained at the very beginning of the algorithm in order to save computation time. If the solution to system RX = A is unique, then the augmented matrix (R∗|A∗) produced by Gauss-Jordan method defines a set of n equalities ∀ i ∈ 1, . . . , n xi = a∗. Step 3 is needed to ascertain that the solution provided by Gauss-Jordan method satisfies i the resource constraints.
Proposition: Time complexity of the provided algorithm is at most O(rn2), where n is the number of tasks in the initial task set and r is the number of resources.
Proof : Examining each step of the algorithm individually, one can make sure that the most time-consuming step is Step 1, at which Gauss-Jordan elimination method is implemented to obtain the reduced row-echelon matrix R∗. The Gauss-Jordan elimination method requires n pivots, each of which requires O(nr) calculations [Nemhauser et al., 1999]. Thus, the time complexity of the algorithm is O(n2).
Note that the algorithm does not provide a solution for any instance of the problem. However, it does demonstrate that the problem is solvable polynomially in at least some of all possible cases. 5
In this paper we’ve examined the problem of estimating maximum resource load for Resource-Constrained Project Scheduling Problem. We’ve discussed in detail the case of two resources and presented two polynomial algorithms for solving this case, with time complexity of at most O(n2). The former allows for relatively easy further modifications, while the latter is easier to implement and was used in designing a pseudopolynomial RCPSP lower bound estimation algorithm. A series of tests was done using PSPLIB benchmark to compare the designed pseudopolynomial algorithm with the best known lower bound estimation algorithms (Tab. 1).
The results demonstrate that for 66% of instances the lower bound provided by our algorithm is not worse than the best known value. For all instances, our value of lower bound is at most 31, 5% worse than the best known lower bound. Moreover, for 5 instances the best known lower bound was improved. We have also discussed the general case of arbitrary number of resources and described how it can be solved in O(rn2) if r ≥ n, full resource consumption is possible and the solution to the corresponding linear equation system is unique. This work was supported by the Russian Science Foundation (grant 17-19-01665).