Ind. Eng. Chem. Res. 1998, 37, 1883-1892
1883
Robust Process Planning under Uncertainty Shabbir Ahmed Department of Mechanical & Industrial Engineering, University of Illinois at Urbana-Champaign, 1206 West Green Street, Urbana, Illinois 61801
Nikolaos V. Sahinidis* Department of Chemical Engineering, University of Illinois at Urbana-Champaign, 600 South Mathews Avenue, Urbana, Illinois 61801
The need to model uncertainty in process design and operations has long been recognized. A frequently taken approach, the two-stage paradigm, involves partitioning the problem variables into two stages: those that have to be decided before and those that can be decided after the uncertain parameters reveal themselves. The resulting two-stage stochastic optimization models minimize the sum of the costs of the first stage and the expected cost of the second stage. A potential limitation of this approach is that it does not account for the variability of the secondstage costs and might lead to solutions where the actual second-stage costs are unacceptably high. In order to resolve this difficulty, we introduce a robustness measure that penalizes secondstage costs that are above the expected cost. Incorporating this measure into stochastic programming formulations does not introduce nonlinearlities, thus making possible the solution of large-scale problems through linear programming techniques. The proposed framework is applied to the planning of chemical process networks for which we propose a number of specialized models and solution schemes. 1. Introduction Decision making in the design and operations of chemical processes is frequently based on parameters of which the values are uncertain. Uncertainties might arise in thermodynamic, kinetic, and other modeling parameters as well as market conditions such as prices of chemicals, product demands, and raw material availabilities. The need to model uncertainty has long been recognized as an essential one in process systems engineering, and many relevant approaches have been proposed including Avriel and Wilde (1969), Nishida et al. (1974), Grossmann and Sargent (1978), Grossmann and Morari (1983), Swaney and Grossmann (1985a,b), Wellons and Reklaitis (1989), Subrahmanyam et al. (1994), Liu and Sahinidis (1996), and Clay and Grossmann (1997). In the so-called two-stage approach, the decision variables of a design or operational problem under uncertainty are partitioned into two sets. The first-stage variables, which are often known as design variables, are those that have to be decided before the actual realization of the uncertain parameters. Subsequently, once the values of the design variables have been decided and the random events have presented themselves, further design or operational policy improvements can be made by selecting, at a certain cost, the values of the second-stage variables, also known as control or operating variables. This corrective action after the design phase is also known as recourse. Due to uncertainty, the second-stage cost is a random * Author to whom all correspondence is addressed. E-mail:
[email protected]. Phone: 217-244-1304. Fax: 217-333-5052.
variable. The objective, then, in the two-stage modeling approach to decision under uncertainty is to choose the design variables in such a way that the sum of firststage design costs and the expected value of the random second-stage recourse costs is minimized. The resulting two-stage stochastic programming models are, in general, difficult to solve. Recent works have demonstrated that it is possible to solve large-scale process planning problems under uncertainty through specialized algorithms for linear and mixed-integer linear programs (Clay and Grossmann, 1997; Liu and Sahinidis, 1996). A potential limitation of the two-stage approach described above is that it only accounts for the expected cost of the second-stage recourse without any conditions on the variability of these costs. It is likely that, under high-risk situations, a design with low expected secondstage cost might run into a specific realization of the uncertain parameters for which the actual recourse cost is unacceptably large. The importance of controlling variability of the costs as opposed to just optimizing its first moment is well recognized in portfolio management applications. To handle the tradeoff associated with the expected cost and its variability in stochastic programs, Mulvey et al. (1995) have recently proposed the concept of robustness. A decision is termed to be robust if the actual cost of the realized scenario remains “close” to the optimal expected cost for all scenarios. This is usually achieved by minimizing the variance of the scenario-dependent costs along with their expectation. Unfortunately, this approach requires the introduction of nonlinearities in the optimization formulation, making the already difficult to solve two-stage stochastic programs even more difficult.
S0888-5885(97)00694-5 CCC: $15.00 © 1998 American Chemical Society Published on Web 04/01/1998
1884 Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 Table 2. Two Solutions for the Motivating Example soln. no. 1 2
x1
y1A y2A y3A y4A
x2
85.42 14.58 58.33 41.67
0 0
0 0
0 0
0 0
y1B
y2B
y3B
y4B
0.01 115.01 14.58 115.01 0 50 41.67 50
Figure 1. Processing network for the motivating example.
Table 3. Expectation and Variability of the Costs for the Two Solutions
Table 1. Demand Scenarios for the Motivating Example
soln. no.
expected total cost
variance
maximum deviation
1 2
255.83 266.68
10747.60 833.33
103.34 16.66
scenario no.
ωkA
ωkB
pk
1 2 3 4
145 145 200 200
135 250 135 250
0.1 0.4 0.4 0.1
The main purpose of this paper is to propose alternative formulations in order to overcome the limitation of the standard two-stage approach described above. Contrary to earlier approaches, the proposed method does not introduce nonlinearities in the formulations, thus making the solution of large problems amenable to linear programming techniques. In section 2, we use an example to illustrate the need for robust optimization tools and review related literature. Section 3 proposes the a new variability index for robust optimization and a theoretical analysis that shows its advantages over the previously proposed robustness measure. In section 4, we apply the proposed ideas to develop robust formulations for the chemical process planning problem and develop a heuristic for solving one of the models. Finally, section 5 presents computational experience with several process planning examples while conclusions are drawn in section 6. 2. Recourse Robustness 2.1. Motivating Example. Refer to Figure 1, and consider the capacity planning decision problem involving processes P1 and P2 and products A and B. We would like to decide how much capacity to invest in the two different processing technologies P1 and P2. The installation cost per unit of capacity is $1 and $3 thousand for processes P1 and P2, respectively. Process P1 can produce 2 tons of product A and 1 ton of product B per unit of capacity installed. Process P2 can produce 1 ton of product A and 3.4 tons of product B per unit of capacity installed. A budget constraint does not allow installation of more than 100 units of total capacity for the two processes. Demands, ωA and ωB, for the products A and B must be satisfied. Any unmet demand has to be outsourced at a cost of $7 and $12 thousand/ ton for products A and B, respectively. There are four demand scenarios whose quantities, ωkA and ωkB, and corresponding probabilities, pk, for k ) 1, ..., 4, are shown in Table 1. The uncertainty in this problem is due to product demand which is not deterministically known. In terms of the two-stage modeling approach, the design variables for this problem are the capacities of the two processes. The control variables, on the other hand, are the amounts of production that are outsourced. First, one must decide how much capacity to install. Subsequently, demand presents itself. At that stage, the installed capacity may or may not be sufficient to satisfy demand. If demand exceeds the installed capacity, a penalty must be paid, as part of the production requirement must be outsourced. The cost of outsourcing is the recourse cost for this example.
To model this problem, let us use x1 and x2 to denote the capacity installed for processes P1 and P2, respectively. Let us also use ykA and ykB to denote amounts of A and B outsourced if demand scenario k materializes. The problem can then be formulated as the following two-stage stochastic linear program: 4
min x1 + 3x2 +
pk(7ykA + 2ykB) ∑ k)1
x1 + x2 e 100
s.t.
(1) (2)
ykA g ωkA - 2x1 - 2x2
k ) 1, ..., 4
(3)
ykB g ωkB - x1 - 3.4x2
k ) 1, ..., 4
(4)
x1, x2 g 0 ykA, ykB g 0
k ) 1, ..., 4
The objective function (1) minimizes the sum of the investment cost and the expected cost for outsourcing demand. Inequality (2) enforces the budget constraint while (3) and (4) ensure that demand is satisfied, if necessary through outsourcing. Two feasible solutions of the above two-stage optimization formulation are shown in Table 2. The first solution relies more heavily on process technology P1, whereas the second solution utilizes more capacity in process P2. Let us now consider two measures of variability of the second-stage costs: the variance and the maximum deviation. These are defined as follows: 4
variance )
∑ k)1
4
pk[(7ykA + 2ykB) -
pj(7yjA + 2yjB)]2 ∑ j)1
maximum deviation ) max{(7ykA + 2ykB) k
K
pj(7yjA + 2yjB)} ∑ j)1 The expected total cost (sum of investment cost and recourse cost), variance, and maximum deviation of the second-stage costs corresponding to the two solutions of Table 2 are presented in Table 3. It can be seen from this table that, although solution 2 is slightly more expensive in terms of the expected total cost, the variance associated with it is significantly smaller than the variance associated with solution 1. The fourth column of Table 3 presents the maximum cost deviation (over all possible scenarios) from the expected total cost. This maximum deviation corresponds to the scenario
Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 1885
for which the recourse cost deviates the most from the expected value. Since there is a positive probability for such a scenario to occur, the table suggests that there is a possibility that the actual cost of solution 1 for some scenario may be 103.33 (or 40%) greater than the expected total cost. In particular, for this example, this happens when scenario 4 occurs. On the other hand, the worst-case cost deviation for solution 2 is 16.67, a mere 6% of the expected cost. A risk-aversive planner might not be willing to tolerate the high variability of the costs associated with solution 1 and would thus opt for the design corresponding to solution 2. The above example illustrates the tradeoffs associated with minimizing the expected recourse cost and minimizing the recourse cost variance in a planning model. In the remainder of this section, we provide the general mathematical description of two-stage programming and discuss two techniques that deal with ensuring cost variability reduction. 2.2. Two-Stage Stochastic Programming. Assuming a finite number of scenarios {1, 2, ..., K} to model the uncertainty, the standard two-stage stochastic linear program with right-hand-side uncertainties is formulated as follows:
SP: min cTx + s.t.
K pkfTyk ∑k)1
Ax e b Tx + Dyk g ωk
(6) ∀k
xg0 yk g 0
(5)
(7) (8)
∀k
(9)
Above, x are the first-stage variables and yk are the second-stage or recourse variables corresponding to scenario k. The uncertain parameter realization and its associated probability under scenario k are denoted by ωk and pk, respectively. Recent applications of two-stage stochastic programming for solving large-scale chemical process planning problems include the linear and mixed-integer linear programming approaches of Clay and Grossmann (1997) and Liu and Sahinidis (1996). 2.3. Robust Optimization Framework. In model SP, there are no restrictions on the variability of the recourse costs. The robust optimization (RO) framework introduced by Mulvey et al. (1995) is a goal programming approach to balance the tradeoffs between expectation and variability of the recourse costs. Let ν(y1,y2,...,yK) be a measure of the variability of the recourse costs. Then, the RO approach is to modify the objective in SP by a weighted variability contribution as follows:
RO: K
T
min c x +
pkfTyk + Fν(y1,y2,...,yK) ∑ k)1
(10)
subject to constraints (6)-(9) where F is a goal programming weight. By changing
the weight F, the relative importance of the expectation and variability of the second-stage costs in the objective can be controlled to construct an “efficient frontier” of solutions. The classical approach to model the tradeoffs between expectation and variability is to use variance as the measure of variability. This gives rise to the wellknown mean-variance model of Markowitz (1959). Malcolm and Zenios (1994) presented an application of this approach to the problem of capacity expansion of power systems. With the variance of the scenariodependent costs included in the objective, the problem was formulated as a large-scale nonlinear program. Escudero et al. (1993) proposed a robust formulation of a production planning problem and recognized the potential problems associated with the nonlinearities. Another application using variance is presented by Lee and Park (1996) for the expansion of chemical processes. Other measures of variability have also been suggested. Paraskevopoulos et al. (1991) developed a capacity planning model for the PVC industry. Robustness was enforced through the incorporation of a penalty function on the sensitivity of the objective to the uncertain parameters. For uncertainty associated with data forecasts, this penalty turns out to be a function of the variance-covariance matrix of the forecast errors. Thus, the developed model is nonlinear. A similar approach was taken by Watanabe and Ellis (1993) to enforce robustness in the problem of allocation of emission reductions to control acid rain. Sensitivity measures were incorporated into a two-stage stochastic programming model resulting in a highly nonlinear formulation. 2.4. Restricted Recourse Framework. The restricted recourse approach of enforcing robustness was recently introduced by Vladimirou and Zenios (1996). Their method differs from robust optimization in that, here, the variability of the second-stage solution is reduced by adjusting the bounds on the second-stage variables. The original formulation of Vladimirou and Zenios (1996) has been developed by the inclusion of a constraint on the mean Euclidean deviation of the second-stage solution vector from the expected solution vector. In this way, the dispersion of the second-stage solutions is restricted to a prescribed level. Instead of restricting the dispersion of the solution vectors, the approach to restrict the dispersion of the random recourse costs can also be considered. Since the variability of the recourse costs is a measure of its dispersion, the restricted recourse formulation of SP can then be stated as follows:
RR: K
min cTx +
pkfTyk ∑ k)1
subject to constraints (6)-(9) and ν(y1,y2,...,yK) e
(11)
where is a prescribed tolerance level for the variability of the second-stage costs.
1886 Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998
In the work of Vladimirou and Zenios (1996), constraint (11) is not considered directly. This allows the problem formulation to be the same as that in SP. The authors describe methods to iteratively enforce constraint (11) by adjusting the bounds on the second-stage variables. It should be noted that solutions obtained by restricting the bounds of the recourse variables are not necessarily optimal to RR. The method of Vladimirou and Zenios (1996) generates sequences of feasible solutions that have smaller cost variability than that of the optimal solution of SP.
RO-UPM: K
min cTx +
∑ k)1
K
pkfTyk + F
pk∆k ∑ k)1
(14)
subject to Ax e b Tx + Dy g ωk
(15) ∀k
(16)
K
∆k g fTyk -
3. Proposed Variability Metric An important issue in enforcing recourse robustness in both the robust optimization and restricted recourse approaches is the choice of the variability metric ν(y1,y2,...,yK). In most applications, the natural choice is variance. There are two main arguments against the use of variance in these formulations: (1) Variance is a symmetric risk measure, penalizing costs both above and below the expected recourse cost equally. As is the case for the chemical process planning problem we are interested in solving and, in general, whenever cost minimization is involved, it is more appealing to use an asymmetric risk measure that would penalize only costs above the expected value. Scenarios whose cost is below the expected recourse cost should not be penalized, as their occurrence leads to lower actual costs than the ones expected. (2) The use of variance results in the inclusion of nonlinearities in the formulation. For capacity expansion problems that involve integer decisions, inclusion of the variance would transform the original mixedinteger linear program into a mixed-integer nonlinear one, making it computationally very difficult to solve. We propose to use the upper partial mean as our measure of variability of the recourse costs. The upper partial mean (UPM) of the recourse costs is defined as follows: K
∆ h (y1,y2,...,yK) )
pk∆k(y1,y2,...,yK) ∑ k)1
(12)
where K
∆k(y1,y2,...,yK) ) max{0, fTyk -
pkfTyk} ∑ j)1
(13)
For scenario k, ∆k is the positive deviation of that scenario’s cost from the expected recourse cost. In this way, ∆ h is defined as the expectation of the positive deviations over all scenarios. Clearly, this index is asymmetric, as it measures only costs that are higher than the expected one. Our definition of the upper partial mean was motivated by the paper of Eppen et al. (1989), where the notion of downside risk is proposed. The key advantage of UPM over downside risk is that UPM does not require the a priori specification of a target level for variance and is therefore more flexible. The upper partial mean can be used as the variance measure in the robust optimization and restricted recourse formulations through linear constraints by introducing additional variables (∆k) as follows:
∑ pk′fTyk′ k′)1
∀k
(17)
∀k
(18)
pkfTyk ∑ k)1
(19)
x, yk, ∆k g 0 RR-UPM: K
min cTx +
subject to constraints (15)-(18) and K
pk∆k e ∑ k)1
(20)
Note that constraints (17) and the nonnegativity of ∆k together with the minimization in the objective satisfy the definition of upper partial mean. Theorem 1. For a given g 0, (x*, y*, ∆*) is an optimal solution of RR-UPM if and only if there exists F g 0 such that (x*, y*, ∆*) is an optimal solution of ROUPM. Proof. Penalizing constraint (20) results in the Lagrangian relaxation of problem RR-UPM:
LR: K
min cTx +
∑ k)1
K
pkfTyk + F(
pk∆k - ) ∑ k)1
(21)
subject to (15)-(18) where F is the corresponding Lagrange multiplier. Since constraint (20) is convex, from strong duality, we have that (x*, y*, ∆*) is an optimal solution of RR-UPM if and only if there exists F g 0 such that (x*, y*, ∆*) is an optimal solution of LR. Clearly, a solution is optimal to LR if and only if it is also optimal to RO-UPM. 0 From the above theorem it can be seen that, for appropriate choices of F and , the formulations ROUPM and RR-UPM are entirely equivalent. In addition to providing the necessary intuition for the development of algorithms, this equivalence is key to establishing stochastic dominance of the solutions of the two models. A solution (x1, y1) is said to be stochastically dominated by another solution (x2, y2) if, for every scenario, the cost of (x2, y2) is at least as small as that of (x1, y1) and, for at least one scenario, the cost of (x2, y2) is strictly smaller. Theorem 2. The optimal solution of RO-UPM and RR-UPM is not stochastically dominated by any other solution. Proof. Let (x1, y1) be an optimal solution to RR-UPM. Let (x2, y2) be another solution to RR-UPM that sto-
Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 1887
chastically dominates (x1, y1). From the definition of stochastic dominance, there exists kˆ such that T
c x2 +
fTyk2 T
T
e c x1 +
c x2 +
fTyk2ˆ
∀ k * kˆ
fTyk1 T
< c x1 +
fTyk1ˆ
(22) (23)
PPP:
∑ ∑ (RitXit + βityit) pk{∑ ∑ δitWkit + ∑ ∑ (ΓjtPkjt - γjtSkjt)} ∑ k∈K t∈u i∈P t∈u j∈C
max NPV ) -
t∈u i∈P
(25)
subject to Multiplying the inequalities (22) and (23) by their corresponding scenario probabilities and summing over all scenarios, we obtain K
cTx2 +
Qit ) Qit-1 + Xit
K
pkfTyk2 < cTx1 + ∑ pkfTyk1 ∑ k)1 k)1
∀ i, t
yitXLit e Xit e yitXU it
This contradicts the assumption that (x1, y1) is an optimal solution of RR-UPM. Thus, an optimal solution is not stochastically dominated by any other solution. From the equivalence of the formulations RR-UPM and RO-UPM, the theorem follows for RO-UPM as well. 0 The theorem argues that the solutions of formulations RO-UPM and RR-UPM are not stochastically dominated by any other solutions. The importance of this property is that one only needs to consider solutions on the expected cost vs UPM efficient frontier. These solutions are not dominated by others irrespectively of the actual realization of the uncertain parameters. Situations in which this property is not necessarily obeyed when the variance is used as the robustness measure are discussed in Eppen et al. (1989) and references therein. 4. Robust Chemical Process Planning The previous section has established the theoretical properties of using the upper partial mean to enforce robustness in stochastic optimization formulations. In this section, we apply the proposed RO-UPM and RRUPM frameworks to the problem of capacity expansion in the chemical process industries. We first present the problem statement and review the two-stage stochastic programming formulation. The RO-UPM and RR-UPM extensions of the SP formulation are then developed. A heuristic for solving the RR-UPM problem is also developed. 4.1. Process Planning under Uncertainty. Given a network of existing and potential processing plants and chemicals, and forecasts for operating and raw material costs, selling prices, demands, and availabilities of the chemicals, the process planning problem consists of determining the operating and capacity expansion policy to maximize the net present value of the operation over a time horizon. The deterministic version of this problem has been formulated as a mixedinteger linear program by Sahinidis et al. (1989). For a survey of literature on this problem and solution methodologies refer to Ahmed and Sahinidis (1998). Assuming linear mass balances for the chemicals across the processes and a finite number of scenarios to model the uncertainty in the availabilities and demands of the chemicals, Liu and Sahinidis (1996) developed a two-stage stochastic programming formulation of the process planning problem. Using the notation provided in the Nomenclature section, their formulation is as follows:
∀ i, t
(27)
∀i
(28)
∑ yit e NEXPi
(24)
(26)
t∈u
∑ (RitXit + βityit) e CIt
∀t
(29)
i∈P
Wkit e Qit Pkjt +
∀ i, t, k
µijWkit ) Skjt + ∑ ηijWkit ∑ i∈P i∈P
(30) ∀ j, t, k
(31)
aLkjt e Pkjt e aUkjt
∀ j, t, k
(32)
dLkjt e Skjt e dUkjt
∀ j, t, k
(33)
Qit, Wkit g 0 yit ∈ {0, 1}
∀ i, t, k
(34)
∀ i, t
(35)
The objective function (25) maximizes the expected net present value, which is the difference between the sales revenues of the final products and the investment, operating, and raw material costs. Equation (26) is a variable bounding constraint for capacity expansions. A zero value of yit forces the capacity expansion of process i in period t to become zero. If the binary variable equals 1, then the capacity expansion is performed within prescribed bounds. Constraint (27) in the above formulation defines the total capacity available at period t as a sum of the capacity available in period t - 1 and the capacity expansion at the beginning of period t. The parameter Qi0 represents the initial capacity, i.e., at t ) 0. Constraints (28) and (29) express bounds on the number of allowable expansions and capital investment, respectively. The condition that the operating level of any process cannot exceed its installed capacity is modeled by constraint (30). Equation 31 expresses mass balances for chemicals across processes. For each chemical, in each time period, the total amount purchased plus the total amount produced from all processes must equal the total amount sold and the total amount consumed by all processes. This balance must be satisfied independently for each time period since it is assumed that no inventories are carried from one period to another. The scenario-dependent bounds on the availabilities and demands of the chemicals are expressed through constraints (32) and (33). Finally, eq 35 expresses the binary nature of variable yit. The above formulation is a two-stage stochastic program of the type SP. The relationship between SP and PPP becomes obvious by recasting PPP as a
1888 Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998
minimization problem and observing that the first-stage problem is defined by variables x ) (Xit, Qit, yit) and constraints (26)-(29). The second-stage problem is then defined by the scenario-dependent variables yk ) (Wkit, Pkjt, Skjt) and constraints (30)-(33). In line with the two-stage philosophy, the capacity expansion decisions are to be made first. Subsequently, once product demands and raw material availabilities have presented themselves, the operating levels of processes will be adjusted so as to maximize revenue. Note that other choices of design and control variables are possible and that the robust optimization framework would extend to any particular choice of the two variable sets. In formulation PPP, the second-stage problem is separable by time periods as well as scenarios. Moreover, the integer variables appear only in the first-stage problem. In this way, decomposition methods such as the sampling-based method of Liu and Sahinidis (1996) are very effective. In such methods, the first-stage problem is typically solved by a branch and bound scheme. Many of the strong bounding techniques developed for the deterministic problem can be used to expedite the convergence of the first-stage problem as reviewed by Ahmed and Sahinidis (1998). 4.2. Robust Optimization Formulation. To develop a robust optimization model of the process planning problem, we introduce the non-negative variables ∆k to measure the positive deviations of the second-stage costs from the expected value. The following constraints are introduced to the second-stage problem:
∆k g
δitWkit + ∑ ∑ (ΓjtPkjt - γjtSkjt) ∑ ∑ t∈u i∈P t∈u j∈C k′ p {∑ ∑ δitWk′ ∑ it + k′∈K t∈u i∈P k′ (ΓjtPk′ ∑ ∑ jt - γjtSjt )}, t∈u j∈C
∀ k (36)
The objective function of PPP will also be modified by recasting the problem as a minimization problem and including a weighted contribution of the upper partial mean. The formulation is as follows:
4.3. Restricted Recourse Formulation. The RR formulation of PPP is obtained by adding a constraint to restrict the upper partial mean to a prescribed level as follows:
PPP-RR:
∑ ∑ (RitXit + βityit) + (ΓjtPkjt - γjtSkjt)} ∑ p {t∈u ∑ i∈P ∑ δitWkit + t∈u ∑∑ k∈K j∈C
min Z )
t∈u i∈P k
subject to constraints (26)-(36) and
∑ pk∆k e
where is the prescribed upper partial mean tolerance. Note that PPP-RR has the same limitations as PPPRO because of the nonseparability of constraints (36) and (39), and the additional variables ∆k. However, without these constraints and variables, the problem is essentially the same as the two-stage stochastic programming model PPP and can be efficiently solved by decomposition. To take advantage of this property, in the next subsection, we do not consider constraints (36) and (39) and the additional variables ∆k explicitly in the problem. The required condition of restricting the upper partial mean to within a prescribed tolerance is achieved by iteratively adjusting the bounds on the second-stage variables. 4.4. Heuristic for the Restricted Recourse Formulation. For PPP, the upper partial mean (∆ h ) of the second-stage costs is defined as
pk[∑ ∑ δitWkit + ∑ ∑(ΓjtPkjt - γjtSkjt) ∑ k∈K t∈u i∈P t∈u j∈C k′ k′ k′ + p { δ W + (ΓjtPk′ ∑ ∑ ∑ ∑ ∑ it it jt - γjtSjt )}] k′∈K t∈u i∈P t∈u j∈C
∆ h )
(40)
˜ kjt, S ˜ kjt) the optimal second-stage Let us denote by (W ˜ kit, P solution for scenario k and define
˜ kit} W ˆ it ) max {W
∀ i, t
(41)
P ˆ jt ) max {P ˜ kjt}
∀ j, t
(42)
S ˆ jt ) min {S ˜ kjt}
∀ j, t
(43)
W h it )
pkW ˜ kit ∑ k∈K
∀ i, t
(44)
P h jt )
pkP ˜ kjt ∑ k∈K
∀ j, t
(45)
S h it )
pkS ˜ kit ∑ k∈K
∀ j, t
(46)
k∈K
min Z )
k∈K
t∈u i∈P k
k∈K
and
subject to constraints (26)-(35) and (36) Note that, due to the nonseparability of constraint (36), we can no longer apply the decomposition methods developed for the stochastic programming model. This, along with the larger number of variables introduced, is the major limitation of the RO framework in the context of process planning. This limitation, in turn, makes it imperative that we avoid the nonlinearities that are introduced by standard robust optimization approaches that use variance in the objective. The model proposed above is a large-scale mixed-integer linear program and can be solved by branch and bound. The bounding methods developed for the deterministic problem can be used to speed up convergence.
(39)
k∈K
PPP-RO:
∑ ∑ (RitXit + βityit) + p {∑ ∑ δitWkit + ∑ ∑ (ΓjtPkjt - γjtSkjt)} + ∑ k∈K t∈u i∈P t∈u j∈C F ∑ pk∆k (37) k∈K
(38)
With the above definitions in (40), we can construct an upper bound on the upper partial mean as follows:
∆ h e
[Γjt(P ˆ jt - P h jt) ∑ ∑ δit(Wˆ it - Wh it) + t∈u ∑∑ j∈C
t∈u i∈P
γjt(S ˆ jt - S h jt)] (47)
Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 1889
In the above bound, the terms in parentheses are merely the upper half-ranges of the second-stage solutions. We can, therefore, lower this bound by restricting these ranges. This can be accomplished by decreasing the upper bounds on the variables Wkit and Pkjt and increasing the lower bound on the variables Skit. However, these variables are tied together through the mass balance equation (31). Thus, to avoid infeasibilities, we only consider restricting the upper bounds on Wkit and Pkjt. The framework of this heuristic procedure is as follows. The indices on the second-stage variables have been dropped for brevity. Step 0. Let be the upper partial mean tolerance. Let UWk and UPk be the current upper bounds on the variables Wk and Pk. Select a parameter λ ∈(0, 1). Step 1. Solve PPP. If the problem is infeasible, stop; no further robustness of second-stage costs can be ˜ k, S ˜ k) be the optimal secondachieved. Else, let (W ˜ k, P stage solutions for scenario k. Step 2. Compute the upper partial mean (∆ h ) for the current solution using (40). If ∆ h e , stop: the current solution is robust within tolerance. Else, go to step 3. Step 3. Calculate W ˆ,P ˆ and W h,P h using (41) and (42) and (44) and (45). Compute the upper half-ranges as follows:
Figure 2. Processing network of example 1.
Figure 3. Processing network of example 2.
RW ) W ˆ -W h ˆ -P h RP ) P Calculate the upper bounds on variables Wk and Pk as follows:
h + λRW} UWk ) min {UWk, W
Figure 4. Processing network of example 3.
UPk ) min {UPk, P h + λRP} Modify the upper bounds in PPP and go to step 1. In the above procedure, the upper bounds on the second-stage variables are iteratively tightened. The half-ranges of W and P are reduced by at least a factor of λ, resulting in a decrease of the upper bound on upper partial mean (47). An important concern is the choice of the parameter λ. Smaller values will lead to faster termination of the procedure but may result in infeasibilities or inferior solutions. The heuristic has the appealing feature that the problem that is solved in step 1 is PPP, a two-stage stochastic program, for which efficient solution strategies have already been developed. A potential drawback of this approach is that, in each iteration, the dispersion of second-stage solutions is restricted with respect to the mean solution of the previous iteration. This condition might be restrictive. Vladimirou and Zenios (1996) describe techniques that might be used to counter this problem but at the expense of a higher computational burden. It should be emphasized that the solution produced by the above heuristic may not be optimal to PPP-RR. This is because the procedure solves problems that have tighter bounds on the second-stage variables than those in PPP-RR. However, note that the objective of the restricted recourse framework is to produce solutions to the two-stage stochastic program PPP such that the variability of the second-stage costs is within a prescribed level. This is accomplished by the above heuristic.
Figure 5. Processing network of example 4.
5. Computational Results In this section, we present computational results from the application of the proposed frameworks to four chemical process planning example problems. The aim of our experiments is to demonstrate how the proposed formulations can serve to determine solutions according to the decision maker’s preference toward risk. The quality of the heuristic solutions for the restricted recourse formulation is also investigated. For all four example problems, the robust optimization formulation was solved using the MIP solver of OSL (IBM, 1991). The heuristic procedure of section 4.4 was coded in GAMS (Brook et al., 1988) using OSL to solve the subproblems. All computations were carried out on an IBM RISC system/6000 model-43P machine with 64 MB of memory. 5.1. Test Problems. The processing networks of the four example problems are shown in Figures 2-5 where processes and chemicals are denoted by squares and
1890 Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 Table 4. Dimensions of the Example Problems size of deterministic equivalent binary continuous problem |P | |C | |u | |K | variables variables constraints 1 2 3 4
3 4 6 5
4 5 6 6
2 2 2 3
100 40 20 30
6 8 12 15
2412 1216 784 1650
3218 1623 1033 2199
Figure 7. Revenue, NPV, and investment vs UPM for example 2.
Figure 6. Revenue, NPV, and investment vs UPM for example 1.
circles, respectively. The first processing network is taken from Ahmed (1997), while the rest are taken from Liu (1995). The dimensions of the networks and number of scenarios are provided in Table 4. This table presents the number of processes, |P |, number of chemicals, |C |, number of time periods, |u |, number of scenarios, |K |, and the size of the deterministic equivalent of the RO formulation. For each problem, we considered all product demands and availabilities as uncertain, thus giving rise to a total of |C | × |u | uncertain parameters. The data for the uncertain parameters were sampled from normal distributions. It is clear from Table 4 that, despite the modest size of the networks, the large number of uncertain parameters leads to large-scale two-stage programming formulations. It should be noted that the RR formulation includes only |u | additional constraints for the restriction of the upper partial mean. 5.2. Effect of Enforcing Robustness. For each problem, a spectrum of solutions was generated by progressively enforcing robustness by increasing the goal programming weight (F) for the RO formulation and decreasing the upper partial mean tolerance () for the RR formulation. Figures 6-9 show the variation of the expected second-stage revenue, expected total NPV, and the firststage investment as a function of the upper partial mean, for both the RO and RR formulations. Although the first-stage investment and second-stage revenues do not show any particular trend by themselves, the expected NPV almost always decreased as robustness was enforced. Figure 10 shows the percentage reduction in the upper partial mean corresponding to the percentage loss in the expected NPV from the stochastic programming solution. From these plots it can be seen that up to 38% of the expected NPV may have to be sacrificed to achieve
Figure 8. Revenue, NPV and investment vs UPM for example 3.
Figure 9. Revenue, NPV, and investment vs UPM for Example 4.
about 95% reduction in the upper partial mean (example 4). However, the slopes of the plots for examples 1-3 are quite steep close to the origin, implying that a substantial degree of robustness can be achieved at little loss to the expected NPV. In fact, the results of example 1 reveal that about 42% reduction in the upper partial mean is possible at the expense of less than 1% of the expected NPV. Such alternative robust solutions could not have been identified from the stochastic programming formulation.
Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998 1891
Figure 12. Iterations for different values of the range restriction parameter. Table 5. CPU Seconds for the RR Heuristic Figure 10. Cost of robustness.
problem
CPU s
problem
CPU s
1 2
20.08 77.73
3 4
8.25 25.11
6. Conclusions
Figure 11. Expected NPV for different values of the range restriction parameter.
The tradeoff between the expected total profit and the upper partial mean can be clearly observed in all cases. These tradeoff curves (“efficient frontiers”) can serve to help the decision maker to select a solution according to his preference toward risk. 5.3. Quality of the RR Heuristic Solutions. Since the RR formulation was solved using a heuristic procedure, the RO solutions can be expected to dominate those of RR. This is clearly observed in Figures 6-9. The presented results were obtained using a range restriction parameter (λ) value of 0.8. The heuristic solutions for the RR formulation exhibit tradeoff behavior similar to that of RO solutions. For larger values of the upper partial mean the heuristic solution approaches the actual frontiers. The difference between the solutions is mainly contributed by the difference between the second-stage revenues. This is because the approach taken by the heuristic is to restrict the second-stage decisions to reduce the variability. The effect of the range restriction parameter (λ) on the RR heuristic is shown in Figures 11 and 12. These curves are for example 1. As low values of λ reduce the dispersion of the second-stage variables at a faster rate, the algorithm terminates faster (Figure 12). However, because of the tighter restriction of the recourse variables, the expected total profit is lower for smaller λ values (Figure 11). Typical computational times for the heuristic are shown in Table 5. These CPU times are obtained using an upper partial mean tolerance () that is 50% of the UPM of the stochastic programming solution. The range restriction parameter (λ) had a value of 0.8.
This paper was motivated by the need to develop a robust optimization framework for the problem of longrange planning in the process industries. The standard stochastic programming formulation of the problem does not address the variability of the uncertain recourse costs across the uncertain parameter scenarios. The need for enforcing the robustness of these costs is particularly important to a risk-aversive planner in a high-variability environment. We extended the stochastic programming formulation to account for robustness of the recourse costs through the use of an appropriate variability criterion. In particular, the upper partial mean has been proposed as the measure of variability for its intuitive appeal and to avoid stochastically dominated solutions and nonlinear formulations. The developed models can be used to generate a spectrum of solutions with varying degrees of recourse robustness. These models provide the decision maker with a tool to analyze the tradeoff associated with the expected profit and its variability. To overcome the difficulty associated with solving the robust models which include nonseparable terms, we developed a heuristic procedure for the restricted recourse formulation. This method iteratively enforces recourse robustness while solving the standard stochastic program in each step. The heuristic generates similar but more conservative tradeoff frontiers for the profit and its upper partial mean. Future research in this area should be directed at the development of efficient exact algorithms for the solution of robust optimization and restricted recourse models. The main difficulty with these problems is the nonseparability of the variability terms, which prevents the direct application of sampling methods. The possibility of using a combination of Lagrangian relaxation and Benders decomposition seems particularly attractive in this context. Finally, it should be clear that the applicability of the proposed variability criterion is not restricted to chemical process planning. For example, its application to other process design and process operations problems could be explored.
1892 Ind. Eng. Chem. Res., Vol. 37, No. 5, 1998
Acknowledgment
ν( ) ) measure of the variability of the recourse costs
Partial financial support by the National Science Foundation under CAREER Award DMII 95-02722 to N.V.S. is gratefully acknowledged.
Literature Cited
Nomenclature Indices i ) for the set P of processes that constitute the network (i ∈ P ) j ) for the set C of chemicals that interconnect the processes (j ∈ C ) k ) for the set K of scenarios (k ∈ K ) t ) for the set u of time periods of the planning horizon (t ∈u ) Variables ∆k ) nonnegative variable to model the upper deviations of the recourse costs under scenario k Pkjt ) units of chemical j purchased at the beginning of period t under scenario k Qit ) total capacity of process i in period t. The capacity of a process is expressed in terms of its main product Skjt ) units of chemical j sold at the end of period t under scenario k Wkit ) operating level of process i in period t under scenario k x ) vector of first-stage variables Xit ) units of expansion of process i at the beginning of period t yk ) vector of second-stage variables under scenario k yit ) 0-1 integer variable; if process i is expanded during period t, then yit ) 1; else yit ) 0 Parameters Rit ) per unit expansion cost for process i at the beginning of period t βit ) fixed cost of establishing or expanding process i at the beginning of period t γjt, Γjt ) forecasted selling and buying prices of chemical j in period t δit ) unit operating cost for process i during period t ) variability tolerance for the RR formulation ηij ) input proportionality constant for chemical j in process i λ ) range restriction parameter for the RR heuristic µij ) output proportionality constant for chemical j in process i F ) goal programming weight for the RO formulation ωk ) second-stage right-hand-side vector for scenario k aLkjt, aUkjt ) lower and upper bounds for the availability of chemical j in period t under scenario k A ) first-stage constraint matrix b ) first-stage right-hand-side vector c ) first-stage cost vector CIt ) investment budget for period t dLkjt, dUkjt ) lower and upper bounds for the demand of chemical j in period t under scenario k D ) second-stage recourse matrix f ) second-stage cost vector NEXPi ) allowable number of expansions for process i pk ) probability of occurrence of scenario k T ) second-stage technology matrix XLit, XU it ) lower and upper bounds for the capacity expansion of process i in period t Functions ∆ h ( ) ) upper partial mean of the recourse costs ∆k( ) ) upper deviation of recourse cost under scenario k
Ahmed, S. Robust Process Planning under Uncertainty. M.S. Thesis, University of Illinois at Urbana-Champaign, Urbana, IL, May 1997. Ahmed, S.; Sahinidis, N. V. Techniques in long range planning in chemical manufacturing systems. In Computer Aided and Integrated Manufacturing Systems: Techniques and Applications; Leondes, C. T., Ed.; Gordon & Breach Publishers: New York, 1998; in press. Avriel, M.; Wilde, D. J. Engineering Design under Uncertainty. Ind. Eng. Chem. Process Des. Dev. 1969, 32, 124. Brook, A.; Kendrick, D.; Meeraus, A. GAMS-sA User’s Guide; The Scientific Press: Redwood City, CA, 1988. Clay, R. L.; Grossmann, I. E. A Disaggregation Algorithm for the Optimization of Stochastic Planning Models. Comput. Chem. Eng. 1997, 21, 751. Eppen, G. D.; Martin, R. K.; Schrage, L. A Scenario Approach to Capacity Planning. Oper. Res. 1989, 37, 517. Escudero, J.; Kamesam, P. V.; King, A.; Wets, R. J.-B. Production Planning via Scenario Modeling. Ann. Oper. Res. 1993, 43, 311. Grossmann, I. E.; Sargent, R. W. H. Optimum Design of Chemical Plants with Uncertain Parameters. AIChE J. 1978, 24, 1021. Grossmann, I. E.; Morari, M. Operability, Resiliency and FlexibilitysProcess Design Objectives for a Changing World. 2nd International Conference on Foundations of ComputerAided Process Design, 1983. IBM. Optimization Subroutine Library, Guide and Reference; IBM Corporation: Kingston, NY, 1991. Lee, H.; Park, S. Long-Range Robust Investment Model for Capacity Expansion of Chemical Processing Networks under Uncertain Demand Forecast Scenarios. AIChE Annual Conference, Chicago, IL, 1996; Paper 142b. Liu, M.-L. Optimization Tools for Process Planning. Ph.D. Thesis, University of Illinois at Urbana-Champaign, Urbana, IL, 1995. Liu, M. L.; Sahinidis, N. V. Optimization in Process Planning under Uncertainty. Ind. Eng. Chem. Res. 1996, 35, 4154. Malcolm, S. A.; Zenios, S. A. Robust Optimization of Power Systems Capacity Expansion under Uncertainty. J. Oper. Res. Soc. 1994, 45, 1040. Markowitz, H. M. Portfolio Selection: Efficient Diversification of Investments; John Wiley & Sons: New York, 1959. Mulvey, J. M.; Vanderbei, R. J.; Zenios, S. A. Robust Optimization of Large-Scale Systems. Oper. Res. 1995, 43, 264. Nishida, N. A.; Ichikawa, A.; Tazaki, E. Synthesis of Optimal Process Systems under Uncertainty. Ind. Eng. Chem. Process Des. Dev. 1974, 13, 209. Paraskevopoulos, D.; Karakitsos, E.; Rustem, B. Robust Capacity Planning under Uncertainty. Manage. Sci. 1991, 37, 787. Sahinidis, N. V.; Grossmann, I. E.; Fornari, R. E.; Chathrathi, M. Optimization Model for Long Range Planning in the Chemical Industry. Comp. Chem. Eng. 1989, 13, 1049. Subrahmanyam, S.; Pekny, J. F.; Reklaitis, G. V. Design of Batch Chemical Plants under Market Uncertainty. Ind. Eng. Chem. 1994, 33, 2688. Swaney, R. E.; Grossmann, I. E. An Index for Operational Flexibility in Chemical Process Design. Part 1: Formulation and Theory. AIChE J. 1985a, 31, 621. Swaney, R. E.; Grossmann, I. E. An Index for Operational Flexibility in Chemical Process Design. Part 2: Computational Algorithms. AIChE J. 1985b, 31, 631. Vladimirou, H.; Zenios, S. A. Stochastic Linear Programs with Restricted Recourse. Eur. J. Oper. Res. 1997, 101, 177. Watanabe, T.; Ellis, H. Robustness in Stochastic Programming Models. Appl. Math. Model. 1993, 17, 547. Wellons, H. S.; Reklaitis, G. V. The Design of Multiproduct Batch Plants under Uncertainty. Comput. Chem. Eng. 1989, 13, 115.
Received for review September 30, 1997 Revised manuscript received February 3, 1998 Accepted February 4, 1998 IE970694T