Bounding in Multi-Stage Stochastic Programming

Problems

Olga Fiedler^{∗}

and

András Prékopa^{†}

31 May 1995

Abstract

We assume that the random variables corresponding to the subsequent periods are all discrete with finite supports. We also assume that the problem is written in the λ-representation form. Starting from the last period, and proceeding in the backward direction, we create first a dual feasible basis which provides us with a lower bound. Next, a primal feasible basis is created which gives us an upper bound. These steps are repeated several times. If the bounds are not satisfactorily, then a few dual steps are performed. The proposed bounding technique is very simple in the case of a multiperiod simple recourse problem. In this case, the dual steps are executed effectively by careful selections of the incomimg and outgoing vectors.

∗Freie Universität Berlin, FB Mathematik u. Informatik, Arnimalle 2-6, 14195 Berlin, Germany †RUTCOR, Rutgers Center for Operations Research, New Brunswick, New Jersey 08903-5062, USA

0. Introduction

A simple Bounding Procedure to obtain fast bounds for the optimum value of the Simple Recourse Stochastic Programming Problem, where the violation of the random constraints are penalized by piecewise linear, convex functions, was proposed Prekopa 1990 (see [8]). Later, Prekopa and Li used this bounding technique in a PERT optimization problem (see [9]) and then generalized it for the case, where the penalty-function is a multivariate polyhedral convex function (see [10]). The main idea is the following: assuming that the random variables are all discrete with finite supports and that the problem is written in the λrepresentation form, one can easily construct an initial dual feasible basis and compute the corresponding lower bound for the optimum value; then, a primal feasible basis is constructed which gives an upper bound. In this paper the Bounding Procedure will be generalized for Multi-Stage Stochastic Programming Problems. If we formulate the multi-stage problem using the ideas from dynamic programming, we can observe that the matrix of the equality constraints (written in the λ-representation form in each stage), has the same block-structure as those in [8] and [10]. Now, starting from the last stage, and proceeding in the backward direction, one can create a dual feasible basis, combining the bases of all stages. Then, starting from the first stage and proceeding in the forward direction, one constructs a primal feasible basis. This provides us with fast lower and upper bounds for the optimum value. The procedure is repeated several times. If the bounds are not satisfactorily close, then a few primal or dual steps are performed. Bounding procedures are very important in single as well as multi-stage stochastic programming problems. In the solving algorithms for single stage problems, such as the dual type method (DTM) of Prékopa (see [8]) and Improved DTM of Fiedler, Prékopa and Fábián (see [4]), one can use them to get fast information about the optimum value or to gain a good starting position to solve the problem. Even more important are bounding procedures in multi-stage problems. In the last decade several authors, such as Birge, Gassman, Higle, Sen, Louveaux, Ruszczynski and Wets, developed various solution techniques for multi-stage stochastic programming problems (see e.g. [1], [5], [6], [2], [3], [7]). These methods are capable to solve problems at most 4-5 stages, because the size of the stochastic programming problem exponentially increases with the number of stages. This is called the “curse of dimensionality” in dynamic programming. Thus, bounding procedures may be the only tools to obtain information about the optimum value in case of a large number of stages, where the solution of the problem cannot be carried out.

1. Formulation of the problem

We start from the following underlying deterministic (N+1)-period linear programming problem with staircase structure:

(1.1)  \begin{aligned}min \{ & q_{0}^{T} x^{0} & + & q_{1}^{T} x^{1} & + & q_{2}^{T} x^{2}+ & … + & q_{N}^{T}x^{N} \} \\ s.t. & A x^{0} & & & & & & & = b \\ & T^{1}x^{0} & + & W^{1}x^{1} & & & & & = ξ^{1} \\ & & & T^{2}x^{1} & + & W^{2}x^{2} & & & = ξ^{2} \\ & & & & & ... & & & ⋮ \\ & & & & & & T^{N}x^{N−1} + & W^{N}x^{N} & = ξ^{N} \\ & x^{0} ≥ 0 , & & x^{1} ≥ 0 , & & x^{2} ≥ 0 , & …, & x_{N} ≥ 0 ; \end{aligned}

here ξ^{t} ∈ R^{n_{t}}, W^{t} ∈ R^{n_{t}} x R^{m_{t}}, t = 1, …, N , b ∈ R^{n_{0}}, x^{0} ∈ R^{m_{0}} and the other data are of suitable sizes, compatible with the formulation (1.1). We assume that the RHS-vectors ξ^{t} , t = 1, …, N are discrete random variables, each with a finite number of possible values ξ_{k}^{t} and with corresponding path probabilities p_{k}^{t} (cf. [5]). Furthermore, we suppose that the technology matrix in stage t is equal to the recourse matrix or its negative in the previous stage t − 1 for t = 2, …, N, , i.e.

T^{t} = ±W^{t−1} (t = 2, …, N).

One formulates the multi-period stochastic programming problem, based on problem (1.1), in the form of a single large scale linear programming problem as follows:

(1.2)  \begin{aligned}min \{ & q_{0}^{T} x^{0} + ∑_{k=1}^{M_{1}}p_{k}^{1}q_{1}^{T} x_{k}^{1} + ∑_{k=1}^{M_{2}} p_{k}^{2}q_{2}^{T} x_{k}^{2} + … + ∑_{k=1}^{M_{N}} p_{k}^{N}q_{N}^{T}x_{k}^{N} \} \\ s.t. & A x^{0} = b \\ & T_{k}^{1}x^{0} + W_{k}^{1}x_{k}^{1} = ξ_{k}^{1}, k = 1, …, M_{1} , \\ & ±W_{α(k)}^{t−1}x_{α(k)}^{t−1} + W_{k}^{t}x_{k}^{t} = ξ_{k}^{t}, k = 1, …, M_{t} , t = 2, …, N \\ & x^{0} ≥ 0 , x_{k}^{t} ≥ 0 , k = 1, …, ∑_{t=1}^{N} M_{t} , t = 1, …, N , \end{aligned}

where W_{k}^{t} is the recourse matrix for node k at stage t, α(k) is the immediate ancestor of node k at stage t, M_{t} is the number of possible values at stage t, t = 1, …, N.

2. The basic idea

We reformulate the problem by the use of the λ-representation and the construct dual and primal feasible bases. First a dual feasible basis (DFB) is constructed which provides us with a lower bound. Then, we construct a primal feasible basis (PFB) and obtain an upper bound. To create the DFB for (1.2) we adopt the idea for the construction of the initial DFB in the DTM, developed by Prékopa (see [8]) for the simple recourse problems, and then extended by Prékopa and Lee (see [10]) for linearly constrained optimization problems with convex polyhedral objective functions. The λ-representation gives rise to a specially block-constrained LP for which the DFB can easily be constructed based on the fundamental property of DFBs for such problems (see Th. 2.1 in [8]). We start with the application of λ-linearization method. For each decision variable x_{k}^{t}, k = 1, …, ∑_{t=1}^{N} M_{t} , t = 1, …, N, we introduce the following functions:

(2.1)  \begin{matrix}f_{k}^{t}(z_{k}^{t}) = & min_{x_{k}^{t}} p_{k}^{t}q_{t}^{T} x_{k}^{t} \\ s.t. & W_{k}^{t}x_{k}^{t} = z_{k}^{t}, x_{k}^{t} ≥ 0, z_{k}^{t} ∈ R^{n_{k}^{t}}. \end{matrix}

Since the optimum value of a linear program, that depends on its RHS parameters z_{k}^{t} ∈ R^{n_{k}^{t}}, is a convex polyhedral function, there exists a subdivision of the space R^{n_{k}^{t}} into the convex polyhedra (simplices) S_{kh}^{t}, with pairwise disjoint interiors, such that f_{k}^{t} is linear on each S_{kh}^{t} and continuous and convex on their union ⋃_{h}S_{kh}^{t}. For fixed k and t let z_{k1}^{t}, …, z_{kH_{k}^{t}}^{t} designate the set of all vertices of all S_{kh}^{t}. Then, the λ-representation of the function f_{k}^{t}(·) is the following:

(2.2)  \begin{aligned}f_{k}^{t}(z_{k}^{t}) = min_{λ} & ∑_{h=1}^{H_{k}^{t}} f_{k}^{t}(z_{kh}^{t})λ_{kh}^{t} \\ s.t. & ∑_{h} z_{kh}^{t}λ_{kh}^{t} = z_{k}^{t} \\ & ∑_{h} λ_{kh}^{t} = 1, λ_{kh}^{t} ≥ 0, \end{aligned}

where the summation over h in the last two rows is the same as in the first one. Let us specialize z_{k}^{t} from (1.2) as follows:

z_{k}^{t} = ξ_{k}^{t} ± W_{α(k)}^{t−1}x_{α(k)}^{t−1}, k = 1, …, ∑_{t=1}^{N} M_{t} , t = 1, …, N,

where α(k) is the immediate ancestor of node k at stage t and W_{α(k)}^{0} := T_{k}^{1}. If we substitute the functions f_{k}^{t}(z_{k}^{t}) by their λ-representations from (2.2), we can reformulate the multi-stage stochastic programming problem (1.2) in the following manner: (2.3)

min_{x^{0},λ} {q_{0}^{T} x^{0}+ ∑_{k=1}^{M_{1}} ∑_{h=1}^{H_{k}^{1}}f_{kh}^{1}λ_{kh}^{1} + ∑_{k=1}^{M_{2}} ∑_{h=1}^{H_{k}^{2}}f_{kh}^{2}λ_{kh}^{2} +  …  +∑_{k=1}^{M_{N}} ∑_{h=1}^{H_{k}^{N}}f_{kh}^{N}λ_{kh}^{N}}
Ax^{0}  = b
T_{k}^{1}x^{0} +  ∑_{h=1}^{H_{k}^{1}} z_{kh}^{1}λ_{kh}^{1}  (k = 1, …, M_{1})  = ξ_{k}^{1}
∑_{h=1}^{H_{k}^{1}} λ_{kh}^{1}  = 1
∑_{h=1}^{H_{k}^{1}} z_{α(k)h}^{1}λ_{α(k)h}^{1} +  ∑_{h=1}^{H_{k}^{2}} z_{kh}^{2}λ_{kh}^{2}  (k = 1, …, M_{2})  = ξ_{k}^{2}
∑_{h=1}^{H_{k}^{2}} λ_{kh}^{2}  = 1
...  ⋮
∑_{h=1}^{H_{k}^{N−1}} z_{α(k)h}^{N−1}λ_{α(k)h}^{N−1} + ∑_{h=1}^{H_{k}^{N}} z_{kh}^{N}λ_{kh}^{N} = ξ_{k}^{N}
∑_{h=1}^{H_{k}^{N}} λ_{kh}^{N}  = 1
(k = 1, …, M_{N})

The matrix of problem (2.3) has a special block-structure which is illustrated in Table 1 in case of three stages, i.e., N = 2, and the number of scenarios at the 3rd stage is equal to the number of descendants at the 2nd stage, i.e., M_{1} = M_{2} = M.

It is easy to see that in our special situation, where T^{t} = ±W^{t−1} all scenarios {ξ_{k}^{t}, k = 1, …, M_{t}} at stage t (t = 2, …, N) belong to at most M_{1} different spaces R^{n_{j}}, j = 1, …, M_{1}, where M_{1} is the number of all possible values of the random RHS-vector ξ^{1}. Let β^{1}(k, t) denotes the 2nd-stage ancestor of node (k, t). Then, if ξ_{β^{1}(k,t)}^{1} = ξ_{j}^{1} ∈ R^{n_{j}}, j ∈ {1, …, M_{1}}, then ξ_{k}^{t} ∈ R^{n_{j}}, too. Let us introduce the following notations:

(2.4)  \begin{matrix}K^{t} := \{1, …, M_{t}\} \\ I_{j}^{t} := \{k | k ∈ K^{t}, ξ_{k}^{t} ∈ R^{n_{j}}\} , j = 1, …, M_{1} \\ K^{t} = ⋃_{j=1}^{M_{1}}I_{j}^{t} \end{matrix} \} t = 1, …, N.

Table 1 Matrix of equality constraints in case of 3-stage problem:

q_{o} f_{11}^{1} … f_{1H_{1}^{1}}^{1} … f_{M1}^{1} … f_{MH_{M}^{1}}^{1} f_{11}^{2} … f_{1H_{1}^{2}}^{2} … f_{M1}^{2} … f_{MH_{M}^{2}}^{2}
A  b
T_{1}^{1} z_{11}^{1} … z_{1H_{1}^{1}}^{1}  ξ_{1}^{1}
0  1… 1  1
⋮  ...  ⋮
T_{M}^{1}  z_{M1}^{1} … z_{MH_{M}^{1}}^{1}  ξ_{M}^{1}
0  1 … 1  1
z_{11}^{1} … z_{1H_{1}^{1}}^{1}  z_{11}^{2} … z_{1H_{1}^{2}}^{2}  ξ_{1}^{2}
0 … 0  1 … 1  1
...  ...  ⋮
z_{M1}^{1} … z_{MH_{M}^{1}}^{1}  z_{M1}^{2} … z_{MH_{M}^{2}}^{2} ξ_{M}^{2}
0 …0  1 … 1  1

The dual vector corresponding to any basis B̂ for problem (2.3) will be partitioned as

(2.5)  \{ y ; \begin{matrix}v_{k}^{1} \\ w_{k}^{1} \end{matrix}  k = 1, …, M_{1} ; ⋮ ; \begin{matrix}v_{k}^{N} \\ w_{k}^{N} \end{matrix}  k = 1, …, M_{N},

where y ∈ R^{n_{0}}, (v_{k}^{1}, w_{k}^{1}) ∈ R^{n_{k}} × R^{1}, k ∈ K^{1}, (v_{k}^{t}, w_{k}^{t}) ∈ R^{n_{j}} × R^{1}, k ∈ K^{t}, and v_{k}^{t} ∈ R^{n_{j}} for k ∈ I_{j}^{t} , i.e. if the 2nd-stage ancestor of the scenario k at stage t ξ_{β^{1}(k,t)}^{1} belongs to the space R^{n_{j}}, j ∈ K^{1}, t = 2, …, N.

3. An Algorithm

Step 0: Set it := 1. For all t = 1, …, N, j = 1, …, M_{1} and for each block k ∈ K^{t} at stage t such that k ∈ I_{j}^{t} choose any (n_{j}+1) vectors of z_{kh}^{t} with subscripts {(k, h_{1}), …, (k, h_{n_{j}+1})} so that the selected vectors constitute the set of the vertices of one of the (subdividing) R^{n_{j}}−dimensional simplices S_{kh}^{t}.

Step 1: Set t := N and f̄_{kh}^{N} := f_{kh}^{N} for all h = 1, …, H_{k}^{N}, k = 1, …, M_{N}. Go to Step 3.

Step 2: Compute the modified costs-coefficients at stage t (corresponding to the vectors selected in Step 1 as follows:

(2.6)  f̄_{kh_{i}}^{t} = f_{kh_{i}}^{t} − v_{k}^{t+1}z_{kh_{i}}^{t} , i = 1, …, n_{j} + 1, k = 1, …, M_{t},

Step 3: Solve the system of linear equations:

(2.7)  \begin{matrix}(z_{kh_{i}}^{t})^{T} v_{k}^{t} + w_{k}^{t} = f̄_{kh_{i}}^{t}, i ∈ \{1, …, n_{j} + 1\}, \\ v_{k}^{t} ∈ R^{n_{j}}, w_{k}^{t} ∈ R^{1}, k ∈ I_{j}^{t} ( ⋃_{j=1}^{M_{1}}I_{j}^{t} = K^{t}) . \end{matrix}

Step 4: If t = 1, then go to Step 5. Otherwise, set t = t − 1, and go to Step 2.

Step 5: Solve the following linear programming problem by any method which produces a pair of primal and dual solutions but not necessarily an optimal basis:

(2.8)  \begin{matrix}min_{x_{1}}\{( q_{0}^{T} − & ∑_{j=1}^{M_{1}}(v_{k}^{1})^{T} T_{k}^{1})x^{0}\} \\ s.t. & Ax^{0} = b , x^{0} ≥ 0. \end{matrix}

Step 6: Compute a lower bound \underline{V} for the optimum value V of problem (2.3):

(2.9)  \underline{V} := b^{T} y + ∑_{t=1}^{N} ∑_{k=1}^{M_{t}} [(v_{k}^{t})^{T} ξ_{t} + 1^{T} w_{k}^{t}] ≤ V ,

where y is the optimal dual vector for the problem (2.8), v_{k}^{t}, w_{k}^{t} (k ∈ K^{t}, t = 1, …, N) are the solutions of (2.7) and the symbol 1 denotes a vector with all components equal to 1.

Step 7: If t > N, then go to Step 8. Otherwise, for all t = 1, …, N, j = 1, …, M_{1} and for each block k ∈ K^{t} at stage t, such that k ∈ I_{j}^{t}, choose any (n_{j}+1) vectors of z_{kh}^{t} with the subscripts (k, h̃_{1}), …, (k, h̃_{n_{j}+1}) , so that the equations

(2.10)  \{ \begin{matrix}∑_{i=1}^{n_{j}+1} z_{kh̃_{i}}^{t}λ_{kh̃_{i}}^{t} = ξ_{k}^{t} ± W_{k}^{t−1} , \\ ∑_{i=1}^{n_{j}+1} λ_{kh̃_{i}}^{t} = 1, λ_{kh̃_{i}}^{t} ≥ 0, k = 1, …, M_{t} \end{matrix}

are satisfied. Note, that the selected vectors z_{kh̃_{i}}, i = 1, …, n_{j} + 1, constitute the set of the vertices of one of the subdividing R^{n_{j}} simplices S_{kh}^{t}; this simplex can differ from the simplex identified in Step 1 for the corresponding stage t. Set t := t + 1, and go to Step 7.

Step 8: Compute an upper bound V̄ for the optimum value V of the problem (2.3):

(2.11)  V̄ := q_{0}^{T} x_{opt}^{0} + ∑_{t=1}^{N} ∑_{k=1}^{M_{t}} ∑_{i=1}^{n_{j}+1} f_{k̃h_{i}}^{t} λ_{kh̃_{i}}^{t} ,

where x_{opt}^{0} is the primal optimal solution of problem (2.8) and λ_{kh̃_{i}}^{t} are defined in Step 7.

□

Remark 2.1 If we combine now the vectors selected in Step 0 with the vectors in the working basis B we obtain a dual feasible basis for the problem (2.3). The dual feasibility of B̂ can easily be shown using similar considerations as those in Lemma 2.1 in [8] and Remark 1.1 in [4]. On the other hand, if we combine the vectors selected in Step 7 with the vectors in the working basis B, we obtain a primal feasible basis for the problem (2.3).

□

Remark 2.2 The described bounding procedure can be iterated as follows: Step 9: Compute

ΔV := V̄ − \underline{V} .

If ΔV is sufficiently small, then STOP. Otherwise, set it := it + 1, replace the subscripts selected in Step 0 for the subscripts defined in Step 7, i.e. set (k, h_{i}) = (k, h̃_{i}) , and go to Step 2.

However, there is no garantee that the new bounds will be better than the previous ones. Therefore, recommend to repeat the Bounding Procedure subsequently several times and accept those bounds for which the duality gap ΔV is the smallest. If ΔV is unsatifactorily large, perform a few dual steps.

□

References

[1] Birge J.R. (1985). Decomposition and partitioning methods for multi-stage stochastic linear programs, Operations Research, 33, 989–1007

[2] Birge J.R., Louveuax F.V. (1988). A multicut algorithm for two-stage linear programs, European Journal of OR 34, 484–392

[3] Birge J.R., Wets J.B. (1986). Designing approximation schemes for stochastic optimization problems, in particular for stochastic programs with recourse, Mathematical Programming 27, 54–102

[4] Fiedler O., Prékopa A., Fábián C.I. (1995). On a dual method for a specially structured linear programming problem, RUTCOR Research Report (to appear)

[5] Gassman H. (1990). MSLiP: A computer code for the multi-stage stochastic linear programming problem, Mathematical Programming 47, 407–423

[6] Higle J., Sen. S. (1991). Stochastic decomposition: an algorithm for two-stage linear programs with recourse, Mathematics of OR 16, 650-669

[7] Ruszczýnski A. (1993). Parallel decomposition of multistage stochastic programming problems, Mathematical Programming 58, 201–228

[8] Prékopa A. (1990). On a dual method... Dual method for the solution of a one-stage stochastic programming problem with random RHS obeying a discrete probability distribution, ZOR-Methods and Models of Operations Research 34, 441-461.

[9] Prékopa A., Li W. (1992). On an optimization problem concerning the stochastic PERT problem, RUTCOR Research Report 18–92

[10] Prékopa A., Li W. (1995). Solution of and bounding in a linearly constrained optimization problem with convex, polyhedral objective function, Mathematical Programming (to appear)