Global optimization algorithm for heat exchanger networks - American

bilinear terms which provide a tight lower bound to the global optimum. ... within a spatial branch and bound method for which branching rules are giv...
0 downloads 0 Views 1MB Size
Znd. Eng. Chem. Res. 1993,32, 487-499

487

Global Optimization Algorithm for Heat Exchanger Networks Ignacio Quesada a n d Ignacio E. Grossmann' Department of Chemical Engineering, Carnegie Mellon University, Pittsburgh, Pennsylvania 15213 This paper deals with the global optimization of heat exchanger networks with fixed topology. It is shown that if linear area cost functions are assumed, as well as arithmetic mean driving force temperature differences in networks with isothermal mixing, the corresponding nonlinear programming (NLP) optimization problem involves linear constraints and a sum of linear fractional functions in the objective which are nonconvex. A rigorous algorithm is proposed that is based on a convex NLP underestimator that involves linear and nonlinear estimators for fractional and bilinear terms which provide a tight lower bound to the global optimum. This NLP problem is used within a spatial branch and bound method for which branching rules are given. Basic properties of the proposed method are presented, and ita application is illustrated with several example problems. The results show that the proposed method only requires few nodes in the branch and bound eearch. Introduction The synthesis of heat exchangers networks (HEN) has received considerable attention over the past two decades (see Gundersen and Naess (1988)). In the recent past a number of mathematicalprogrammingmethods have been proposed (Floudas et al., 1986;Yee and Grossmann, 1990) that can automatically synthesize networks and account for the different tradeoffs involved in the synthesis (energy, area, and number of units). Also these methods can handle special constraints such as forbidden matches and no stream splits. A major limitation, however, of these methods is that HEN problems involve inherent nonconvexities which can lead to multiple local optimal solutions. Therefore, when traditional local search techniques are used to solve these problems, the solution will depend on the initial point and no guarantee can be given of global optimality. Global optimization of heat exchanger networks has been addressed with stochastic and deterministic methods. As an example for the former type of methods, Dolan et al. (1989) have applied simulated annealing to the synthesis of heat exchanger networks. Although this method has proved to be successful in terms of finding improved solutions compared to the ones reported in the literature, rigorous global optimality cannot be guaranteed and the computationalexpense can be high. Deterministic methods for HEN'S with fixed topology have been proposed by Westerberg and Shah (1978) and by Floudas and Ciric (1989). The former authors developed a Lagrangian decomposition scheme which, however, requires the global optimization of subproblems. The latter authors applied generalized Benders decomposition by treating the flows as complicating variableswith which the NLP subproblems and the master problems are convex. The limitation of this method is that global optimality cannot be guaranteed since the master problem may not predict valid lower bounds. In this paper the optimal design of heat exchanger networks with fixed topology is considered. It is assumed that the network corresponds to a particular configuration within the superstructure given by Yee and Grossmann (1990)in which a sequence of stages with isothermalmixing points is considered. Even under some further simplifying assumptions, the optimal design problem may have multiple local solutions. In particular, by assuming linear costs for the areas and arithmetic mean temperature

* Author to whom correspondence should be addressed.

differences for the driving forces, it will be shown that the NLP can be formulated with linear constraints but with nonconvexities in the objective function that are given by the sum of linear fractional terms. Recently there has been an important effort in addressing the global optimization of nonconvex nonlinear programming problems (e.g., see Horst (1990)). Falk and Palocsay (1992)proposed a rigorous algorithm for solving these NLP problems with fractional terms in the objective. The algorithm consists of a sequence of linear programming (LP) bounding problems in which cuts from the contours of the objective function are added. The drawback of this method is that convergence can be slow. The NLP model considered in this paper is also closely related to bilinear programming problems which can be solved to global optimality with several methods, including the ones proposed by AlKhayyal and Falk (1983), Swaney (1990), Floudas and Visweswaran (1990), and Sherali and Alameddine (1992). All these methods have relied on the use of linear underestimators. Rather than applying these more general solution techniques, it is the objective of this work to develop a specialized method that can better exploit the structure of the HEN problem considered in this paper. In this paper, a global optimization algorithm is proposed that is based on a new underestimator problem that involves both linear and nonlinear estimator functions which provide an exact approximation in the boundary of the feasible region. The proposed underestimator problem, which gives rise to a convex NLP problem that predicts tight lower bounds of the original problem, is coupled with a spatial branch and bound search procedure for finding the global optimum solution. Several examples are presented to illustrate the properties and performance of the algorithm. As will be shown, the results are very encouraging since only a small number of subproblems are required in the branch and bound search. Mathematical Model Two major simplifications have been assumed in the optimization model of heat exchanger networks that provide a mathematical structure that can be exploited for developing a rigorous method for the global optimization. The area cost is given by a linear function, and the driving force for the heat exchangers is calculated by the arithmetic mean temperature difference. Although these assumptionsare somewhatrestrictive, the proposed model should still be useful for practical applications. The arithmetic mean temperature difference is often very close to the logarithmic mean and underestimates the area. The

o ~ ~ ~ - ~ ~ ~ ~ I ~ ~ 0I 1993 ~ ~American ~ ~ - Chemical o ~ ~ Society ~ ~ o ~ . o o / o

488 Ind. Eng. Chem. Res., Vol. 32, No. 3, 1993 cot $10, ooo

c, T3

4

Q

*

initial point. If C C C* set C* = C and store the solution as the incumbent solution ( *). If CL~+'C C L ~invert +~ and Ok+2 in F. Go to step 2.

Convergence As for the convergence of the algorithm, it should be noted that Al-Khayyal and Falk (1983) and Sherali and Alameddine (1992) presented branch and bound algorithms with partition rules that are similar to the one used here. The convergence proof given below is in the same spirit of the one given by Sherali and Alameddine (1992). Property 5. The algorithm will either terminate in a finite number of partitions at a global optimal solution or generate a sequence of bounds that converge to the global solution. Proof. Given the branch and bound procedure, there are two possibilities. In the first one, at a given node the lower bound CLof the underestimator NLPL is identical to the original objective function, in which case the algorithm terminates in a finite number of partitions. In the second possibility an infinite sequence of partitions is generated. This in turn implies that there is a subregion that is being infinitely partitioned. Let the sequence of solutions be denoted by superscript k' and z = [Q, A, AT, XI. By the termination criteria it is known that CUh' - CLkl> 0

(37) Since the upper bound is at least as strong as the evaluation of the actual objective function for the current solution zk, C(zk) - C,(zk') 1 CQk'- C/ > 0 (38) there must exist an exchanger, m,for which its feasible region is infinitely partitioned. By the partition rule 1,

($-Ai)

s

(g

- A,,,),

i = 1,...,n

(39)

Summing over all the exchangers i, it follows that

n(

2

- A,,,) 1 C(zk')- CL(zk')2 CUlk'- C/ > O

(40)

The variables for exchanger m have some bounds defining an interval. Since the partition is of the same nature as that the one used by Al-Khayyal and Falk (1983), the variables in the sequence must converge to one of the bounds. Moreover, the series has to converge to a point since when one of the bounds of a variable are not changing this one is selected as the partition variable in the algorithm. When one of the variables is at the boundary the representation is exact, Qm/ATm= A,. Therefore, 0 L C(zk')- CL(Zk')L CUVk'- C F > 0

(41)

which means that equality between the lower bound CL and the original cost function C must hold. Since by Property 4 C L ~is' a lower bound to the global optimal solution, it corresponds to the global solution. rn

Examples In this section several examples are presented that illustrate the performance of the algorithm. The size and computational results of the examples are summarized in Table 11. The NLP time is the time used for solving the NLPL problems, and the LP time is the time for solving linear problems to obtain the initial bounds and subsequent

Ind. Eng. Chem. Res., Vol. 32, No. 3, 1993 495 Table 11. Size and Computational Results of Example Problems size'

initial objective

problem

original

NLPL

CL

C*

global optimum

no. of nodes

LP timeb

NLP timeb

example 1 example 2 example 3 example 4 example 5 example 5 (c = 0.2%) example 6 example 7

(12,13) (12,13) (7,6) (11,9) (11,9)

(16,31) (16,31) (9,W (14,23) (14,23)

32300 18640 9.49 6427 6269

36163 19168 9.49 6427 6482

36163 19168 9.49 6427 6398

5 5 1

(20,21) (26,30)

(27,52) (34,68)

242.78 536.76

252.8 536.76

245.6 536.76

3.20 3.20 0.65 1.50 2.70 (3.90) 5.23 5.03

0.90 1.05 0.10 0.22 1.25 (2.20) 2.10 0.43

0

1

5 (9) 3 1

(rn, n): rn = number of variables; n = number of constraints. CPU seconds IBM/R6000-530.

Table 111. Initial Lower and Upper Bounds for Example 1 area (m*) heat exchanger LB UB 1 11.11 73.52 2 4.24 47.43 3 0.106 44.74 4 0.162 14.70

heat load (kW) LB UB 978.125 55.55 21.875 944.44 21.875 944.44 55.55 978.125

1

AT (K) LB 50 51.5 21.11 66.5625

UB 133.03 199.11 205.625 343.33

Nonconvex

45 \ \

Table IV. Area for First Convex Underestimator in Example 1 heat convex exact incumbent exchanger solution solution solution 1 59.5 67.4 73.521 2 14.6 27.0 4.248 3 1.4 1.4 0.0106 4 6.0 6.0 14.698

updates. Note that in three of the examples only one node was required, which means that the problem was solved with the underestimator NLPLwithout requiring branch and bound enumeration. These problems were solved with MINOS 5.2 using the modeling language GAMS on an IBM/R6000-530. Note that all the examples required less than 10 s of CPU time. Moreover, the time for solving the LP bounding problems can be further reduced by using a warm start. Example 1. Consider the motivating example introduced earlier in this paper. If bounds are first computed for the variables, this yields the values shown in Table 111. The solution that had the lowest cost among these calculations has a cost of C* = $36 163 and the areas (Ai/ Vi) are shown in Table IV as the incumbent solution. The initial NLPL problem is then constructed using the nonlinear and linear estimator functions. In this example it is possible to obtain projections for the upper bound of the driving force over its heat load for exchangers 3 and 4 that can be used togenerate underestimators of the form (24). These projections are given by

AT3 I210 - 0.2Q3

(42)

AT4 I360 - 0.3Q4

(43)

The importance of these projections is that they provide an exact approximation for exchangers 3 and 4 as it can be seen by the solution of the first NLPLin Table IV. The solution is CLO = $32 300. The actual objective function, the objective function of a linear underestimator problem, and the objective function of the convex NLPLare plotted versus 81 in Figure 7. As can be seen, when only linear overestimators are used in problem NLPL a lower bound of CLP= $27 700 is obtained. Therefore, using nonlinear underestimators, the gap between the lower bound and the incumbent solution is reduced from 23.4?6 to 10.6 '3%. Since there is a gap between the incumbent solution C* and CLO,a partition is made. The largest difference in the

40

Convex

I\

\

I ,'I

200

4 0 0 6 0 0 8 0 0

I.&*Qp.l

Figure 7. Plot of first convex underestimator in example 1.

approximations corresponds to exchanger 1. Since the incumbent solution lies a t an extreme point, it cannot be used to partition the space. Instead the convex solution is used to generate two subregions (A2 1 27 m2 and A2 I 27 m2 for Q1 and Q2, respectively). For each one of these subregions the bounds on Q2 and AT2 are recalculated taking into account the partition constraint. The solutions of these two new convex problems are C L =~ $39 360 and C L= ~ $36 070. The first is greater than the incumbent solution (C* = $36 163), so no further partition of this subregion is required. For the second subregion the lower bound is below the incumbent solution (0.2?6) and a new partition is done similar to the first one. The only different is that now the second exchanger is selected for partition. The solutions of the new subregions are C L =~ $36 270 (A1 I73 m2) and C L = ~ $36 163 (Al 1 73 m2). The first of these subregions can be rejected while the first is equal to the incumbent solution, and therefore the global solution with AI = 73.521 m2, A2 = 4.248 m2,A3 = 0.106 m2, A4 = 14.698 m2, and C = $36 163 has been obtained. The approximations of the objective function in the different subregions can be seen in Figure 8. The problem required a total of 4.1 s to find the global optimal solution as seen in Table 11. Example 2. The same network as in example 1 is considered with the data given in Table V with ci = $10oO/ m2,EMAT = 10, and the overall heat-transfer coefficients U1= U,= 0.1 kW/(K m2)and Us = Ua = 1.5 kW/(K m2)). In this case there are two local minima and their objective functions are close (C= $19 564 and C = $19 168). The algorithm obtains a tight lower bound in the first iteration (C= $18 640) and behaves in a similar fashion as in the first example. The global optimum with A1 = 9.647 m2, A2 = 7.756 m2, A3 = 0.577 m2, A4 = 1.188 m2, and C = $19168 is obtained after the solution of five NLP underestimator problems (see Table 11). Example 3. The relevance of the projected underes-

496 Ind. Eng. Chem. Res., Vol. 32, No. 3, 1993

cost $10. ooo

Nonconvu

function

Undenrtlmrtor I

I

1

(a)

Figure 8. Partitions for example 1: (a) first partition; (b) second partition. 500

450

I

Total Area Imz)

la, Fcp=l kW/K

Fcp=l kW/K

300

T2

TI

Figure 9. Network for example 3.

I

Table V. Data for Example 2 stream

Fcp (kWIK)

CI

1 1 1 1 1

cz c3

Hi

Hz

temperature (K) inlet outlet 250 370 352.5 580 512.5

400

In this case the cost function is given by the total area and the formulation is given by Q2

Qi

minA=-+O.lAT,

Q1

150 - T 2

T - 300

Q2 = 400 - T

Q1 =

Q2

I

t

AT2 =

(b)

Figure 10. (a) Original objective function for example 3. (b) Objective function with linear and nonlinear underestimatom.

projectionsare generated for the upper bound of the driving forces in terms of their heat loads (22), the following inequalities are obtained:

ATl I150 - Q1

0.1AT2

+ Ti

Total Area (m?

380 362.5

9).

ATl =

QM

(a)

timators (24) is illustrated by the following example. Consider the problem in Figure 9. It consists of a cold stream that goes through two heat exchangers in series. For the two hot streams the inlet temperatures are given and the outlet temperatures are not specified (see Figure

subject to

*

I

i8.5

IOO-T+T2 2

(44) AT2 5 100 (45) In this way, from (24), a new underestimator function can be generated for the first exchanger:

450 - T I

= 500-

7'2

300 5 T 5400, T , I450, T , I500

Although the formulation of this problem is nonconvex, there is only one optimal solution. Moreover, if the objective function is plotted versus the heat load of the first heat exchanger (81)it is a convex function (see Figure loa). When the algorithm is applied using only the nonlinear and underestimator and linear overestimator functions (1)-(4) to approximate the nonconvex terms, it is not possible to obtain an exact approximation of the convex objective of this problem which is shown in Figure lob. If

Once this projected underestimator (46) is added to the NLP underestimator problem, an exact approximation is obtained. This occurs because the original formulation (P3) can be reformulated to obtain the following convex problem: min A =

Qi

0.1(150- Q1)

+-

Q2

0.1(100)

(47)

When the objective function of NLPLO is plotted, it matches exactly the original function in Figure loa. In this way the solution of the underestimator problem has

Ind. Eng. Chem. Res., Vol. 32, No. 3, 1993 497 278

HI

3

3

H.

TI T2 Figure 11. Network for example 4.

T3

520

Table VI. Data for Example 4 stream c1

Hi Hz H3

Fcp (kW/K) 10 10 10

10

temperature (K) inlet outlet 300 400 450 500 550

Temperature stream Inlet ouuet

c

293

HI HI

353 383 395 343 405 288

c: z

483

Table VII. Data for Example 5 stream c1

Hi

Hz H3

Fcp (kW/K) 10 10 10 10

temperature (K) inlet outlet 300 400 490 500 550

a total area of ALO = 9.49 m2 and it is possible to prove global optimality in one iteration, by simply evaluating the actual objective function of this feasible solution and obtaining an incumbent solution of A* = 9.49 m2with A1 = 2.34 m2 and A2 = 7.15 m2. As seen in Table 11, the solution of this problem required only0.75 8. This example then shows that the proposed method has the capability of recognizing convexity in a problem and obtaining its solution efficiently. Example 4. This example consists of one cold stream and three hot streams in a network of three heat exchangers in series (Figure 11). The data are given in Table VI. This network is similar to the one presented in example 2 with the difference that now the problem does not project in a convex form. The objective function is to minimize the cost where Ci = $1000/m2 and Vi = 1 kW/(K m2). The algorithm is applied, and it is possible to generate two extra projected underestimators (24) with the projections ATl I150 - O.lQl

(48)

AT2 I200 - 0.1Q2 (49) The solution of the first convex underestimator problem is C = $6427. When the real objective function for this feasible Solution is evaluated, the incumbent solution is C* = $6427, and hence the global solution with A1 = 0.0 m2,A2 = 1.547 m2,and A3 = 4.880 m2 is obtained in only one iteration, As seen in Table 11, this problem required 1.72 8. Example 5. The same network as in the previous example is considered, with the new data given in Table VII. In this case the global solution is not an extreme point like in the previous one and the approximations are not exact in the first iteration. The solution of this first convex problem is CLO= $6269, and the incumbent solution is C*= $6482. The approximations for the first and third heat exchangers are exact, so the second heat exchanger is selected. Here the linear overestimators are redundant and the partitions are made over the heat load. T w o subregions are obtained with the constraints Q2 2 26.8 and Q2 I26.8 with convex solutions of C L =~ $6408 and C L =~ $6278, respectively. The incumbent solution now is C* = $6408, so the first subregion does not require more

/H

Figure 13. Network for example 7.

partitioning. For the remaining subregion a partition is made, taking now the driving force of the second heat exchanger. The subregions are AT2 I165.89 and AT2 2 165.89 with solutions C L ~= $6355 and C L ~= $6373, respectively. No better upper bound is obtained, and since the lower bound and upper bound are close, the original problem is solved given a new incumbent solution of C* = $6398. The subregionsare within 0.6 % of the incumbent solution, and the algorithm stops with A1 = 0.625 m2, A2 = 1.220 m2, and A3 = 4.553 m2. If a smaller c is used (c = 0.2%),four more NLPL’s are solved. Thus, with the larger tolerance 3.95 sand five NLP’swere required, while with the smaller tolerance nine NLP’s were required with 5.1 s as seen in Table 11. Example 6. The next example is the one presented by Colberg and Morari (1990) and Yee et al. (1990). Here a fixed configurationis considered and the network is shown in Figure 12. The objective is to minimize the area. T w o local solutions are listed in Table VIII. Applying the algorithm, it is possible to generate projections for three of the seven heat exchangers. In the initial underestimator a lower bound of 242.78 m2 is obtained. The actual objective function for this solution is 252.8 m2. The largest difference corresponds to A c ~ Hand ~ , the level chosen for the partition is 24.696. For A c ~ HI ~24.696 the convex solution has an objective of A = 246.39 m2, and for the other subregion the solution is A = 245.1 m2. A better incumbent solution is obtained with an objective of A = 245.6 m2, and the solution is within 0.2% of the global solution. Only three convex NLP underestimator subproblems were required to converge in a total of 7.33 s as seen in Table 11. Example 7. The following example consists of a network reported in Grossmann and Floudas (1987). The

498 Ind. Eng. Chem. Res., Vol. 32, No. 3, 1993 Table VIII. Areas for Local Solutions in Example 6 solution AHIC~ AHPCI AHZCZ 1 38.27 325.11 73.18 2 30.53 45.41 107.9 ~~

~

Aci 6.55 6.70

AHI 1.59 0.00

Ac2 4.15 2.74

AHZ 42.76 52.34

total 491.61 245.62

Table IX. Data for Example 7 stream

FCD(kW/K) 40 20 20 25 50 60 20

temperature (K) inlet outlet 325 400 1350 450 1360 400 298 430 540 540 380 310 410 290 340 285

configuration is illustrated in Figure 13,the data are given in Table IX, and a EMAT = 5 is considered. The objective function is to minimize the total area of the network and Vi = 0.5 kW/(K mz). For this problem it was possible to identify the following six projection terms: AT, I90 - 0.0225Q1 AT, L 90 - 0.0225Q1 AT4 I46.24 - 0.008Q4 AT, I63.66 - O.OllQ, AT6 I13 - 0.005Q6 AT6 1 13 - 0.005Q6 (50) The solution of the first underestimator problem is CL, = 536.76 mz, and the evaluation of the original function has the same value. Therefore, the optimal solution with A1 = 53.33 mz, AZ= 100 mz,A3 = 34.28 mz,A4 = 87.96 m2, Ab = 146.94 m2, As = 52.17 mz, A7 = 26.67 m2, and A6 = 35.41 mz is obtained after 5.46 s (see Table 11).

Conclusions This paper has presented a global optimization algorithm for heat exchanger networks that can be formulated in terms of linear fractional functions in the objective function and linear constraints. The key elements in the proposed algorithms are the proposed convex nonlinear underestimators for fractional terms which complement the linear underestimators for bilinear terms by McCormick (1976). As has been shown with thenumerical results, the resulting NLP underestimator problem provides rigorous tight lower bounds to the global optimum with which the computational effort in the spatial branch and bound method is greatly reduced. In fact, for all cases except one, only a maximum of five nodes had to be enumerated in the example problems. Finally, it should be noted that the method proposed in this paper has been generalized to nonlinear programs that involve convex, linear fractional, and bilinear terms in the objective function and constraints (Quesada and Grossmann, 1992). Also, work is under way to be able to handle concave cost functions and logarithmic mean temperature differences in the heat exchanger network, as well as the optimization of the configurations. Acknowledgment The authors would like to acknowledge financial support from the Engineering Design Research Center at Carnegie Mellon University. Appendix A. Linear and Nonlinear Under- and Overestimators A concave overestimating function of a product of functions is given by (McCormick, 1976)

C,(Y) 1 g(Y)

(A31

Considering the function AiATi, where f ( x ) = Ai and g(z) = ATi, the concave functions Cf(x)and C,(Y) are functions f ( x ) and g(Y), respectively. The concave overestimating function, which is linear, is given by

AiATi Imin[A:ATi

+ AiATy - AFAT?, AyATi +

AiATk - AyATpl (A4) In a similar way as in (Al), the convex underestimating function of a product of functions is given by f ( x ) g(Y) 1m=[fuc,(Y)

+ gucf(x) - fugu, P'c,b) + &cf(x) - ?'&I (A51

where f", gut f L , and gL are positive bounds over the functions f ( x ) and g(Y) such that

f. If ( x ) IP

(A61

&Ig(Y) I g " (A71 and c f ( x ) and c&) are convex functions such that for all x and y in some convex set

If ( x )

(A81

c,(Y) SgcY)

(A91

Cf(X)

On the basis of (A5) one can generate nonlinear underestimators for the rational terms in the objective function of problem (P) as follows. In particular, consider the function Qi/ATi,wheref(x) = Qiandgb) = l/ATi. The functions and bounds are given by c f ( x ) = Qi

(A101

1 c,(Y) = A Ti

(All)

Qk IQi IQiu

(A121

From (A5) the convex underestimator function is given by Qi

Qi

A Ti

Appendix B. On the Exact Approximation of the Linear and Nonlinear Estimators in the Boundary The estimator functions A4 and A14 have the property that they match the function when one of the variables is

Ind. Eng. Chem.Res., Vol. 32, No. 3, 1993 499 at a bound. This is because the individual convex and concave approximation functions in (A2), (A3), (A8), and (A9) are the functions themselves. Equation A4 for the linear overestimator reduces as follows: if Ai = A?

A?ATi Imin[A)ATi + AyAT? - AiLAT?, AyATil = min[AyATi (A? - A:)(ATiU - ATi),A?ATil = AYAT, (Bl)

+

if Ai = A:

A)ATi I min[A:ATi, A)ATi + (A? - AF)(ATi - AT?)] = A t A T i (B2)

if^^^ = AT? AiAT? Imin[A,AT?, AiATy + (A? - Ai)(AT? - AT))] = AiAT? (B3) if

AT^ = AT;

AiATF Imin[A,AT) + (A, - A?)(ATiU - AT?), AiAT?l = AiAT) (B4) Equation A14 for the convex nonlinear underestimator reduces as follows: i f a T i = AT?

Similarly when ATi = ATiL

For Qi = Qiu

Finally for Qi = QiL

Literature Cited Al-Khayyal, F. A.; Falk, J. E. Jointly constrained biconvex programming. Math. Oper. Res. 1983,8,273-286. Charnes, A.; Cooper, W. W. Programming with linear fractional functions. Nav. Res. Logistics Q. 1962,9, 181-186. Colberg, R. D.; Morari, M. Area and capital cost targeta for heat exchanger network synthesis with constrained matches and unequal heat transfer coefficients. Comput. Chem. Eng. 1990,14, 1-22. Dolan, W. B.; Cummings, P. T.; LeVan, M. D. Process optimization via simulated annealing: application to network design. AIChE J. 1989,35,725-736. Falk, J. E.; Palocsay, S. W. Optimizing the sum of linear fractional functions. In Recent Advances in Global Optimization; Floudas, C. A., Pardalos, P. M., E&.; Princeton University Press: Princeton, NJ, 1992;pp 221-258. Floudas, C. A,; Ciric, A. R. Strategies for overcoming uncertainties in heat exchanger network synthesis. Comput. Chem. Eng. 1989, 13,1133-1152. Floudas, C. A.; Visweswaran, V. A global optimization algorithm (GOP) for certain classes of nonconvex NLPs-I Theory. Comput. Chem. Eng. 1990,14,1397-1417. Floudas, C. A.; Ciric, A. R.; Grossmann, I. E. Automatic synthesis of optimum heat exchanger network configurations. AZChE J. 1986, 32,276-290. Grossmann, I. E.; Floudas, C. A. Active constraint strategy for flexibility analysis in chemicalprocess. Comput. Chem. Eng. 1987, 11,675-693. Gundenen, T.; Naess, L. The synthesis of cost optimal heat exchanger network synthesis-an industrial review of the state of the art. Comput. Chem. Eng. 1988,12,503-530. Horst, R. Deterministic method in constrained global optimization: Some recent advancesand fields of application. Nav. Res. Logistics 1990,37,433-471. McCormick, G. P. Computability of global solutions to factorable nonconvex programs: Part I Convex underestimator problems. Math. Program. 1976,10,147-175. Quesada, I.; Grossmann, I. E. Global optimization algorithm for linear fractional and bilinear programs; Engineering Design Research Center, Report EDRC-0614392;Carnegie Mellon University: Pittsburgh, PA; November 1992. Sherali, H. D.; Alameddine, A. A new reformulation-linearization technique for bilinear programming problems. J. Global Optim. 1992,2,379-410. Swaney, R. E. Global solution of algebraic nonlinear programs. Presented at AIChE Meeting, Chicago, IL, 1990; Paper No. 22f. Westerberg, A. W.; Shah, J. V. Assuring a global optimum by the use of an upper bound on the lower (dual) bound. Comput. Chem. Eng. 1978,2,83-92. Yee, T. F.; Grossmann, I. E. Simultaneous optimization models for heat integration-11. Heat exchangernetwork synthesis. Comput. Chem. Eng. 1990,14,1165-1184. Yee, T. F.; Grossmann, I. E.; Kravanja, Z. Simultaneous optimization models for heat integration-I. Area and energy targeting and modeling of multi-stream exchangers. Comput. Chem. Eng. 1990, 14,1151-1164. Received for review June 1, 1992 Revised manuscript received November 20, 1992 Accepted December 4, 1992