Ind. Eng. Chem. Res. 2008, 47, 3771-3782
3771
PROCESS DESIGN AND CONTROL Global Capital/Total Annualized Cost Minimization of Homogeneous and Isothermal Reactor Networks Wen Zhou and Vasilios I. Manousiouthakis* Chemical and Biomolecular Engineering Department, UCLA, Los Angeles, California 90095
The minimum capital cost/total annualized cost (MCC/MTAC) problem is investigated, within the infinitedimensional state-space (IDEAS) framework, for homogeneous and isothermal reactor networks. The resulting mathematical formulation is an infinite-dimensional program with linear constraints and a separable concave objective function to be minimized. The global optimum of this infinite program is approximated through global solution of a series of finite-dimensional programs with concave separable objective functions and linear feasible regions. A branch-and-bound algorithm is used to solve globally each of these finite programs. The proposed methodology is illustrated with a reactor network synthesis case study aiming at capital cost minimization. 1. Introduction Reactor network synthesis (RNS) has been an active area of chemical engineering research for over half a century. Over this time period, it has become evident that RNS can be an instrumental tool in lowering capital investment and operating costs and in advancing environmentally conscious manufacturing and pollution prevention. Objectives commonly employed in RNS aim at efficiently using raw materials and resources and effectively reducing equipment size. Prominent among these objectives are maximum yield, maximum selectivity, minimum utility cost, minimum capital cost (MCC), and minimum total annualized cost (MTAC). Despite the fact that process synthesis has been developed into an active research topic, very few systematic methods have been proposed for the MCC/MTAC synthesis of reactor networks. One such method is superstructure-based optimization, which provides general process networks (flowsheets) that incorporate an a priori selected number and type of units and their connections.1-5 Solutions of MCC/MTAC synthesis problems, using the superstructure approach, were demonstrated mainly on reactor-separator-recycle systems, where the objectives were total annual cost or profit to be minimized or maximized.6-8 The capital cost of reactors was formulated as a linear function of reactor volume,6 a power function of reactor molar holdup,7 or a power function of reactor volume.8 Although the superstructure method leads to a systematic synthesis methodology, generally the resulting solution is merely guaranteed to be locally optimal and is only as rich as the proposed superstructure. Bedenik et al.8 tried to upgrade the capabilities of the superstructure method by applying a hierarchical multilevel mixed-integer nonlinear programming (MINLP) synthesis method, where a flexible superstructure was employed instead of a fixed superstructure. Profitable superstructures were found by analysis of the intermediate MINLP solutions in the attainable region (AR) and * To whom correspondence should be addressed. Tel.: +1-310-2060300. Fax: +1-310-206-4107. E-mail address:
[email protected].
used for the next iteration. Although it helps to abate the deficiency of the superstructure method, namely, that solutions are only as rich as the proposed superstructure, the method cannot claim that all alternatives are included during the search. In addition, efforts have been made to find the global optimum of certain special superstructure problems with compromise of the generality of the superstructure.9 The results verified that no general global solution procedure exists at present for a general MINLP program with convex/linear objective functions that result from the superstructure method. On the other hand, targeting-based methods for RNS are not capable of handling MCC/MTAC problems.10,11 The infinite-dimensional state-space (IDEAS) framework was first proposed by Wilson and Manousiouthakis12 and then applied to a variety of process network synthesis problems: (i) multiple-component mass exchange networks;12 (ii) distillation networks;13-15 (iii) reactor networks;16,17 (iv) reactive distillation networks;18 (v) power cycle networks;19 and (vi) separation networks.20 In those works, objectives are constructed as linear functions of the framework’s variables, such as maximum yield,16,17 minimum utility cost,12,13,15,19 minimum area,14,20 and minimum plate holdup.18 IDEAS gives rise to an optimization program whose feasible region is constrained by an infinite number of linear equations. When the objective is linear/convex, solutions to the resulting infinite-dimensional linear program (ILP) are guaranteed to be globally optimum. In this work, the IDEAS framework is applied, for the first time, to MCC/MTAC RNS problems, which are shown to give rise to concave, separable objective functions. For such problems, local optima are not guaranteed to be global. Nevertheless, the resulting infinite-dimensional concave linear program (ICLP) can be solved globally by employing a two-layer algorithm, consisting of finite approximation and branch-and-bound underestimation. In addition, techniques of dimensionality reduction are presented that can be used to abate the associated computational burden. The proposed algorithm and dimensionality-reduction techniques are illustrated on a capital cost minimization case study for a reactor network involving
10.1021/ie060653+ CCC: $40.75 © 2008 American Chemical Society Published on Web 04/05/2008
3772
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
reactor volume (m3); τ be the reactor residence time (s); and λ be a binary integer flag to identify the reactor type (i.e., λ ) 0 for PFR, λ ) 1 for CSTR). Then the reactor model can be written as follows:
{
}
if λ ) 0
∀k ) 1, ..., n}
if λ ) 1 (2)
dCk n ) ∀ k ) 1, ..., n ) Rk({Cj}j)1 dτj in ∀k ) 1, ..., n Ck|τh)0 ) Ck Ck|τh)τ ) Cout k
in out n {Cout k - Ck ) τ Rk({Cj }j)1)
(1)
V ) τF
Figure 1. IDEAS representation for homogeneous and isothermal reactor networks.
Trambouze reaction kinetics. The work is demonstrated on reactor networks where the constant-density assumption is valid. However, as shown at the end of the paper, the variable-density condition brings no conceptual limitations to the proposed method, but does bring additional complexity. The rest of this paper is organized as follows. In section 2, the IDEAS framework is reviewed, and the linearity of its feasible region is explained. In section 3, the optimization formulation of the MCC/MTAC RNS problems, within the IDEAS framework, is presented. In section 4, the two-layer solution procedure is described. In section 5, some properties of RNS problems are presented. In section 6, the case study is presented, and conclusions are drawn. 2. Linearity of the IDEAS Feasible Region IDEAS provides a general network representation that is able to address most process network synthesis problems. Different from traditional optimization approaches, in IDEAS stream extensive (quantity) properties (e.g., flowrates) and stream intensive (quality) properties (such as concentrations, molar enthalpies, molar fractions, etc.) are dealt with separately. The overall process network is decomposed into two parts: the distribution network (DN) and the process operator (OP) (Figure 1). The DN contains all possible mixing, splitting, recycling, and bypassing operations. The OP quantifies the action of the process unit operations under investigation. Indeed, chemical process models are represented through an infinite-dimensional linear OP that contains an infinite number of linear input-output maps. Each map, corresponds to an operating unit, has prior specified values of certain intensive properties, design parameters, and linearly related (i.e., linearly constrained) extensive variables. In addition, values of selected intensive properties remain unchanged when values of extensive variables change, as optimization is carried out. Under steady-state and constantdensity assumptions, homogeneous and isothermal reactor networks are considered, where reactor types of continuous stirred tank reactor (CSTR) and plug-flow reactor (PFR) are allowed. Next, it is shown that CSTR and PFR models can be represented by such an infinite-dimensional linear OP. Let F be the volumetric flowrate in and out of the reactor (Fin ) Fout ) F by constant-density assumption) (m3/s); Cin k (Cout k ) be the molar concentrations of the kth component at the reactor inlet (outlet), respectively (kmol/m3); Rk be the reaction generation rate of the kth component (kmol/(m3 s)); V be the
(3)
In this model, the volumetric flowrate and the volume of the reactor are extensive properties; the inlet/outlet concentrations are intensive properties, and the reactor residence time and the integer flag are design parameters. Let vectors u and y be defined as follows: T in out u ) [u1 u2]T; u1 ( [F]; u2 ( [Cin 1 ... Cn C1 λ] out T y ) [y1 y2]T; y1 ( [F V]T; y2 ( [Cout 2 ... Cn τ]
Then the following nonlinear input-output information map Φ: u f y can be established,
[]
[] [
u y Φ (u , u ) Φ:u ) u1 f y ) y1 ) 1 1 2 Φ2(u1, u2) 2 2
]
(4)
where Φ1 and Φ2 are viewed as components of Φ mapping u to y1 and y2, respectively. It is easy to verify that the map possesses two properties: • 1st property: y2 ) Φ2(u1, u2) ) Φ2(u2), i.e., there is no direct dependence between u1 and y2. • 2nd property: y1 ) Φ1(u1, u2) ) Φ1(u2)u1 where Φ1(*) is a linear operator for fixed u2. Indeed, (1),(2)
Φ2:u2 98 y2; Φ1(u2) ( [1 τ ]T Consequently, by fixing u2 at any possible condition, one can ∞ create an infinite sequence that is denoted as {u2(i)}i)1 . For each u2(i), y2(i) is calculated as Φ2(u2(i)), and a linear map Φ1(u2(i)), through which u1(i) maps to y1(i), is generated. Therefore, the action of quantifying process unit operations in the IDEAS OP is captured by a linear mapping whose domain and range are subsets of infinite-dimensional spaces. Given the linearity of the IDEAS OP, the linearity of the IDEAS feasible region can be ensured by establishing that the DN description is also linear. Indeed, the DN contains solely mixing and splitting operations. Since the intensive properties concerning any stream leaving or entering the DN are known, the resulting balance equations are linear. As a result, the linearity of the OP and the DN gives rise to an infinitedimensional linear feasible region for the IDEAS optimization formulation. It is worth pointing out that there are a number of ways to attain the linear OP (i.e., to define the vectors u and y); different ways of attaining linearity of the OP may result in different IDEAS optimization formulations with different computational characteristics. 3. IDEAS Formulation for the MCC/MTAC Problem Figure 1 illustrates the IDEAS representation of a reactor network, where the system has n components, NI network inlet
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3773
streams, and NO network outlet streams. Under steady-state, homogeneous, isothermal, and constant-density conditions, the infinite-dimensional linear feasible region for the corresponding mathematical formulation of the RNS problems can be written as
F ∑ i)1
F (j) )
OI
(i,j) +
F (i) )
F ∑ j)1
OI
(i,j) +
F (i) )
F ∑ j)1
(i,j)
∀i ) 1, ..., NO (6)
FOOˆ (i,j) + ∑ FˆIOˆ (i,j) ∑ i)1 i)1
∑ j)1
∀j ) 1, ..., NI (5)
∀i ) 1, ..., ∞
(7)
∞
NI
CkˆI (i)FˆI(i) )
OO ˆ
FˆII(i,j) + ∑ FˆIOˆ (i,j) ∑ j)1 j)1 NO
FOˆ (j) )
ˆII
∞
NI
ˆI
F (i,j) ∑ i)1 ∞
NI
O
CIk(j)FˆII(i,j) +
∀j ) 1, ..., ∞ (8)
∞
COkˆ (j)FˆIOˆ (i,j) ∑ j)1
∀ k ) 1, ..., n ∀i ) 1, ..., ∞ (9) ∀i ) 1, ..., NO (10)
(FO(i))low e FO(i) e (FO(i))upp ∞
NI
(COk (i))lowFO(i) e
∑ j)1
CIk(j)FOI(i,j) +
(COk (i))uppFO(i)
COkˆ (j)FOOˆ (i,j) e ∑ j)1
∀k ) 1, ..., n ∀i ) 1, ..., NO (11)
FOˆ (i) ) FˆI(i) V(i) ) τ(i)FˆI(i)
∀ i ) 1, ..., ∞ ∀ i ) 1, ..., ∞
(12) (13)
FI g 0; FO g 0; FˆI g 0; FOˆ g 0; FOI g 0; FˆII g 0; FOOˆ g 0; FˆIOˆ g 0; V g 0 (14) where sequences (CO)low ((CO)upp) and (FO)low ((FO)upp) are lower (upper) bounds of sequences CO and FO, respectively. Equations 5-8 correspond to total mass balances in the DN; eqs 9, 10, and 11 represent component balances in the DN; and eqs 12 and 13 represent the action of the OP, for constant-density reactors. In addition, COˆ and CˆI satisfy the PFR or CSTR models:
{(
)
∞
cν(i)V(i)e ∑ i)1 ∞
∞
NO
I
RNS problems can be expressed by eqs 17 and 18, respectively.
}
dCk(i) n ) Rk({Cj(i)}j)1 ) ∀ k ) 1, ..., n dτj ˆI ∀k ) 1, ..., n Ck(i)|τh)0 ) Ck(i) Ck(i)|τh)τ(i) ) COkˆ (i) if λ(i) ) 0 ∀i ) 1, ..., ∞ (15)
n {COkˆ (i) - CkˆI (i) ) τ(i) Rk({COjˆ (i)}j)1 ) ∀k ) 1, ..., n} if λ(i) ) 1 ∀i ) 1, ..., ∞ (16)
3.1. Formulation of Objective Function. The capital cost of a reactor is a function of its volume, which usually can be approximated by a power function: cνVeν, where cν is the capital cost coefficient (a constant) and eν ∈ (0, 1).21 On the other hand, the total annualized cost of a reactor is a function of its volume and flowrate, which can be approximated by cµVeµ + cˆ µF, where cµ is the annualized capital cost coefficient (a constant), cˆµ is the operating cost coefficient (a constant), and eµ ∈ (0,1). Accordingly, the objective functions for the MCC and MTAC
∑ i)1
ν
(17)
∞
cµ(i)V(i)eµ +
cˆ µ(i)FˆI(i) ∑ i)1
(18)
Combination of the feasible region constraints, shown in eqs 5-14, and the objective function in eqs 17 and 18 gives rise to the following infinite-dimensional concave linear programs (ICLPs) for the synthesis of reactor networks featuring MCC and MTAC, respectively.
MCC: ν ) inf eq 17 s.t. eqs 5-14
(19)
MTAC: µ ) inf eq 18 s.t. eqs 5-14
(20)
It is important to point out that objective functions 17 and 18 reflect economies of scale with respect to reactor volume. Indeed, eν ∈ (0, 1) implies that two reactors with volumes V1 and V2, respectively, cost more than one reactor with volume V1 + V2. This property of the objective functions 17 and 18 thus naturally leads problems 19 and 20 to possess globally optimal solutions typically featuring a small number of units. Of course, this desirable feature of problems 19 and 20 comes in tandem with the undesirable feature of an increased computational burden for the globally optimal solution. This is further elaborated in the next section. 3.2. Global Optimality Property. • Proposition 1. All self-recycle streams of CSTRs can be set to zero at the global minimum of the MCC/MTAC RNS problems. Proof. Consider an optimal reactor network containing a CSTR with self-recycle. Let Fr represent the self-recycle flowrate. Without loss of generality, we mix all streams entering the CSTR except the self-recycle stream and use F and C to represent the flowrate and concentration vector, respectively. The following model can be written as follows: in out n (F + Fr)(Cout k - Ck ) ) VRk({Cj }j)1) out (F + Fr)Cin k ) FCk + FrCk
∀k ) 1, ..., n
∀k ) 1, ..., n
Substituting the second equation into the first one gives rise to out n F(Cout k - Ck) ) VRk({Cj }j)1)
∀ k ) 1, ..., n
It implies that a variation from this optimal network, where the flowrate of the CSTR self-recycle stream is shut down, is feasible. In addition, the objective value after the variation remains the same for the MCC problem and even decreases for the MTAC problem. As a result, it can be concluded that selfrecycle streams of CSTRs can be set to zero at the optimum of the MCC/MTAC RNS problems. 9 4. Solution Procedure The aforementioned ICLPs may possess local optima that are different than the global optimum. The solution procedure of seeking the global optimum is, thus, carried out with a twolayer algorithm. 4.1. Outer Layer: Finite Approximation Strategy. In the outer layer, the solution of the ICLP is approximated by a series
3774
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
of finite-dimensional concave linear programs (FCLPs) of increasing size, whose sequence of solutions converges to that of the ICLP. For example, let L represent the number of OP inlet streams, outlet streams, and reactor units and let νL(µL) represent the optimum objective value of the corresponding MCC (MTAC) FCLP. Then the MCC (MTAC) FCLP is formulated as follows:
νL ) inf eq 17 s.t. eqs 5-14 with replacement of ∞ by L (21) µL ) inf eq 18 s.t. eqs 5-14 with replacement of ∞ by L (22) ∞ It can be established that {νL}L)1 is a nonincreasing sequence of upper bounds of ν and is guaranteed to converge to ν. ∞ Similarly, the nonincreasing sequence {µL}L)1 can be calculated and it converges to µ. The way chosen to obtain the linear OP is used to generate finite-dimensional programs. Indeed, instead of considering all possible conditions of u2 in eq 4, a finite number of possible conditions of u2, based on discretization, are utilized to generate a finite linear OP, consisting of a finite number of linear feasible units (i.e., reactors in this work). One possible way to generate these finitely many conditions of u2 in this work is to first identify lower/upper in out bounds on Cin 1 , ..., Cn , C1 . These bounds can be derived based on physical considerations (e.g., concentrations are non-negative) or from brief analysis of the kinetics. On the basis of these ranges, defined by the lower/upper bounds, one can then establish a (m + 1) point grid in each of the (n + 1) ranges for in out n+1 Cin 1 , ..., Cn , C1 . By doing so, one generates (m + 1) in in out T possible values of [C1 ... Cn C1 ] . As a result, by using eqs 1 and 2 to calculate y2 ) [Cout ... Cout τ]T for λ ) 0 and λ ) 2 n 1, respectively, one can generate a finite number of feasible CSTRs/PFRs and corresponding OP inlet/outlet streams up to the number of 2(m + 1)n+1 (some may be infeasible). It can be seen that the number of elements contained by u2 is critical in determining the size of the FCLP for a certain discretization level and, consequently, the computational burden for solving the resulting FCLP. Efforts to abate the size of the FCLP generated without compromising the global optimality for the same level of discretization will be discussed later in Section 5. 4.2. Inner Layer: Branch-and-Bound Underestimation. In solving a FCLP, solutions from employing standard nonlinear programming solvers cannot be guaranteed to be globally optimal, since the objective of a FCLP is not linear/convex. Falk and Soland22 proposed a branch-and-bound algorithm to solve separable nonconvex programming problems in the following general form:
m ˆ
min φ(x) ) x
φi(xi) ∑ i)1
s.t. x ∈ Ω xlow e x e xupp where x ) [x1, ‚‚‚, xmˆ ]T is a vector, and the term “separable” means that φi is only a function of xi. Assumptions for the algorithm include the following: (1) and xupp are finite. (3) Each The set Ω is closed. (2) xlow i i φi(xi) is lower semicontinuous on Ψi ) {xi|xlow e xi e i low e x e xupp} is xupp }. (4) Ω ∩ Ψ, where Ψ ) {x|x i nonempty. These assumptions are sufficient to ensure that φ(x) attains its minimum over the set Ω ∩ Ψ. The
algorithm solves a sequence of problems, in each of which the objective function is convex. More detailed description of the algorithm and its convergence theorems are given in Appendix A. Consider a FCLP (i.e., problem (νL) or (µL)), and let xi ) V(i), m ˆ ) L, and Ω be defined by the finite-dimensional linear feasible region of the FCLP; V(i)e is a separable, concave function of V(i) for e ∈ (0, 1). To apply Soland’s branch-andbound algorithm to the FCLP, finite lower/upper bounds for each V(i) are needed. An obvious lower bound is zero. However, it is not trivial to find a finite upper bound. One wants to set the upper bound high enough to contain all possible global optima and low enough to increase the convergence speed of the branch-and-bound algorithm. The upper bound could, but not necessarily, be identical for all V(i)s and be equal to a Vmax possibly provided by the reactor manufacturer in association with the cost functions 17 or 18. A possible tighter upper bound can be identified using the total reactor volume minimization problem: L
inf
V(i) s.t. eqs 5-14 with replacement of ∞ by L ∑ i)1
(23)
This problem formulation leads to a linear program, whose global solution can be determined through straightforward simplex procedures. Let Vtv(i) (FˆItv(j)) be the optimal value of V(i) (FˆI(j)) in the above problem (eq 23), and let Vcc(i) be the optimal value of V(i) in the problem (νL). Then the following holds (since cν > 0): L
cν(i)(Vcc(i))eν e
L
cν(j)(Vcc(j))e e∑ cν(j)(Vtv(j))e ∑ j)1 j)1 ν
ν
∀ i ) 1, ..., L
The upper bound of V(i) can then be set identical to ((1/cν(i)) L cν(j)(Vtv(j))eν)1/eν for the MCC problem. Similarly, Let Vtac∑j)1 (i) be the optimal value of V(i) in the problem (µL). Then the following holds: L
cµ(i)(Vtac(i))eµ e
∑ j)1
L
cµ(j)(Vtv(j))eµ +
cˆ µ(j)FˆI (j) ∑ j)1 tv
∀ i ) 1, ..., L
The upper bound of V(i) can be found for the MTAC problem. Given that the lower and upper bounds of the concave variables are finite, by Theorem A3, solution of a FCLP can be guaranteed in a finite number of steps using the weak refining rule. 5. Dimensionality Reduction for Single Inlet/Outlet Reactor Networks 5.1. All Components Are of Interest. It can be easily verified that FI(1) ) FO(1) for single inlet/outlet reactor networks (i.e., NI ) NO ) 1). Zhou and Manousiouthakis23 presented and proved a number of novel properties of RNS problems within the IDEAS framework for isothermal reaction systems; these properties were then used to develop dimensionality-reduction theorems. According to the theorems, it is possible to reduce the dimension of the concentration space, within which RNS (i.e., the AR/shrink-wrap methodology) is conducted, by investigating the reaction scheme. Although the primary interest is the dimensionality reduction of AR construction problems, as pointed out in that paper, the results are applicable to other
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3775
RNS methods. Results pertinent to this work are presented below, and proofs are referred to the original paper. • Theorem 1: If n
ηiRi ) 0 ∑ i)1
ηi ∈ R
(24)
holds for a reaction scheme, then for any feasible well-connected reactor network, containing a finite number (L) of active reactors, the following holds: n
ηk(COk (1) - CIk(1)) ) 0 ∑ k)1
(25)
n
∑ ηk(CkˆI (i) - CIk(1)) ) 0 k)1
∀i ) 1, ..., L
(26)
in out T u2 ( [Cin 1 ... CnzC1 λ]
n
∑ ηk(COkˆ (i) - CIk(1)) ) 0
∀i ) 1, ..., L
(27)
k)1
Thus, it can be concluded that any stream of the reactor network n holds the relation ∑k)1 ηk(Ck - CIk) ) 0, where Ck is the component concentration of the stream. The systematic means to identify conditions like eq 24 are provided in the original paper, including examining reaction stoichiometry. It helps to reduce the dimensionality of a FCLP in two ways: in the feasible reactor generation stage and in the optimization formulation stage. To demonstrate, let the following relations be valid for a reaction scheme, nb
Rj )
ηjlRl ∑ l)1
∀j ) nb + 1, ..., n
(28)
where ηjl is the linear coefficient. From Theorem 1, the following then holds,
Cj ) Cj
∀j ) 1, ..., nb
nb
(Cj -
C0j )
)
ηjl(Cl - C0l ) ∑ l)1
5.2. Part of Components Are of Interest. One of the important facts is that quite often only part of all components are of interest, and their generation rates do not depend on all component concentrations. This specific knowledge certainly can be used to reduce the dimensionality of a FCLP. Let Jx be a subset of Jn ) {j, j ) 1, ..., n} containing the indices of concentrations of interest, Jz be the least subset of Jn such that (1) it contains Jx and (2) Rj (∀ j ∈ Jz) depends solely on{Ci}i∈Jz. Thus, in the worst scenario, Jz ) Jn; and in the best scenario, Jz ) Jx. In addition, let Jy be a subset of Jz containing the indices of concentrations on which Rj (∀ j ∈ Jz) depend. Apparently, some indices of Jx may not belong to Jy, but the complement of Jx with respect to Jz has to belong to Jy. Without loss of generality, let Jz contain indices of the first nz elements. Thus, in the feasible reactor generation stage, u2 and y2 can be redefined as follows:
∀j ) nb + 1, ..., n (29)
where C0j is the concentration in the network feed (i.e., CIj ). Any feasible reactor in a FCLP, whose concentrations of inlet or outlet stream do not satisfy linear relations (eq 29), has to be inactive in the optimal network, i.e., the flowrate, reactor volume, and associated cross-flowrates are zero at optimum. In other words, these reactors could be eliminated from the FCLP in the feasible reactor generation stage. Considering that the way chosen to obtain the linear OP is also used to generate finite feasible reactors as demonstrated in section 4.1, the above statement is equivalent to redefining the u2 and y2 as follows: in out T u2 ( [Cin 1 ... CnbC1 λ] T out out y2 ( [Cninb+1 ... Cin n C2 ... Cn τ]
Furthermore, in the optimization formulation stage, only the first nb component balance equations are necessary since the remaining (n - nb) balance equations will be automatically satisfied according to linear relations (eq 29), when the first nb component balance equations are satisfied.
out T y2 ( [Cout 2 ... Cnz τ]
That is, some of the CSTR/PFR models are solved (i.e., model nz equations regarding {Cj}j)1 ) and other component concentran tions (i.e., {Cj}j)nz+1) remain unknown. In addition, in the optimization formulation stage, only the first nz component nz in eqs 9 and 11) are balance equations (i.e., regarding {Cj}j)1 used. This treatment is based on the reason that, since performance of reactors and reactor networks (residence time, objective function, concentration restrictions, etc.) does not rely on concentrations of these n - nz components (i.e., n {Cj}j)n ), the suboptimization problem excluding these conz+1 centration variables can be solved and gives rise to the same optimum solution. The values of these excluded concentration variables can be found later by solving the feasibility problem (i.e., running a simulation of the resulting network). Indeed, this method has been often used in RNS. For example, when the AR is constructed for the Trambouze reaction scheme (i.e., A f B, A f C, A f D), and the A and C components are of interest, the AR is constructed in the 2-dimensional (CA, CC) concentration space, since the generation rates of A and C components (i.e., RA, RC) depend solely on concentrations of A and C. It is worth pointing out that dimensionality-reduction techniques presented in section 5.1 can also be employed here as long as linear relations like eq 24 can be found within {Rj}j∈Jz. In the next two subsections, a common situation will be considered where, for components of interest that do not belong to Jy, only lower and/or upper bounds are enforced for these components’ concentrations at the network outlet, when the MCC/MTAC is sought. It can be shown that the dimensionality of a FCLP can possibly be reduced further. 5.2.1. Concentration Bounds on Components of Interest at Network Outlet. The following proposition can first be established. • Proposition 2: For single inlet/outlet reactor networks, the following holds: ∞
FO(1)COk (1) - FI(1)CIk(1) )
(FˆI(i)(COkˆ (i) - CkˆI (i))) ∑ i)1
∀k ) 1, ..., n (30)
i.e., the change of component molar flowrate from network inlet
3776
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
to network outlet is equal to the sum of the changes of component molar flowrate for all individual reactors in the OP. Proof. Let C be the concentration vector [C1...Cn]T. Since FI(1) ) FO(1),
F (1)C (1) - F (1)C (1) ) F (1)(C (1) - C (1)) O
O
I
I
O
O
FCLP can be reduced further. To demonstrate, let Jy (Jz) contain indices of the first ny (nz, nz > ny) elements. First in the feasible reactor generation stage, u2 and y2 can be redefined as follows: in out T u2 ( [Cin 1 ... CnyC1 λ]
I
(31)
out T y2 ( [Cout 2 ... Cny ∆Cny+1 ... ∆Cnz τ]
From eqs 6 and 11, the following holds: ∞
[FOOˆ (1,j)(COˆ (j) - CI(1))] ∑ j)1
FO(1)(CO(1) - CI(1)) )
(32)
In addition, in the optimization formulation stage, only the first ny ny component balance equations (i.e., regarding {Cj}j)1 in eqs 9 and 11) are used, and the following bounding constraints are used:
From eqs 7 and 9, the following holds: ˆI
∞
ˆI
F (i)(C (i) - C (1)) ) I
[F ∑ j)1
(COj (1))low e COj (1) e (COj (1))upp ˆIO ˆ
O ˆ
(i,j)(C (j) - C (1))] I
∀i ) 1, ..., ∞ (33)
(COj (1))low e
L
1 I
F (1)
[
(FˆI(i)∆Cj(i))] + CjI(1) e (COj (1))upp ∑ i)1
Summation of eqs 32 and 33 yields: ∞
FO(1)(CO(1) - CI(1)) + ∞
∑ j)1
∑ i)1
[(FOOˆ (1,j) +
[FˆI(i)(CˆI(i) - CI(1))] )
∞
(FˆIOˆ (i,j)))(COˆ (j) - CI(1))] ∑ i)1
(34)
Substitution of eq 8 gives rise to ∞
FO(1)(CO(1) - CI(1)) +
[FˆI(i)(CˆI(i) - CI(1))] ) ∑ i)1 ∞
[FOˆ (j)(COˆ (j) - CI(1))] ∑ j)1
(35)
Using FOˆ (i) ) FˆI(i) yields ∞
FO(1)(CO(1) - CI(1)) )
[FˆI(j)(COˆ (j) - CˆI(j))] ∑ j)1
9 Extension of Proposition 2 to a general reactor network is obvious and is given as follows without further proof. • Proposition 3: For a reactor network with NI network inlet streams and NO network outlet streams, the following holds NO
∑ i)1
NI
FO(i)CkO(i) -
FI(i)CkI(i) ) ∑ i)1
∞
[FˆI(i)(COk (i) - CIk(i))] ∑ i)1
∀k ) 1, ..., n (36)
i.e., the change of component molar flowrate from network inlet to network outlet is equal to sum of the changes of component molar flowrates for all individual reactors in the OP. From the constant-density CSTR/PFR models, it can be easily seen that, when a component’s generation rate (denote the index of the component as i1) is not a function of its own concentration, its concentration change through a reactor (denoted as ∆Ci1) can be calculated and can be viewed as a reactor parameter like residence time. As a result, balance equations of this component are not needed, if constraints involving this component (e.g., its bounds at network outlet) can be expressed in terms of its concentration change (i.e., ∆Ci1); thus, the dimensionality of a
∀j ) 1, ..., ny
∀ j ) ny + 1, ..., nz (37)
nz Values of {Cj}j)n can be found later by running a simulay+1 n . tion of the resulting network as well as {Cj}j)n z+1 5.2.2. Global Optimality Property. The variation-induced minimization (VIM) technique is used to derive global optimality properties that optimal reactor networks of this problem should satisfy based upon the mathematical formulation of a FCLP. VIM is a technique in which problem variables are varied through a series of feasible perturbations, so as to establish which ones could be set to zero at the global optimum, thus reducing the size of the FCLP under investigation. For those components of interest that do not belong to Jy, their lower/ upper bounds on network outlet concentrations can be of different types. Let index subset Jx1, Jx2, Jx3 and Jx4 be defined as follows: Jx1 ) {j:j ∈ Jx, j ∈ Jy}; Jx2 ) {j:j ∈ Jx, j ∉ Jy, (COj (1))low g 0, (CjO(1))upp ) ∞}; Jx3 ) {j:j ∈ Jx, j ∉ Jy, (COj (1))low ) 0, 0 < (COj (1))upp < ∞}; Jx4 ) {j:j ∈ Jx, j ∉ Jx1, Jx2, Jx3}. Considering the reactor network of the FCLP operated at an optimal steady state, the following proposition can be established. • Proposition 4: If there are two reactors (denoted as i1 and i2) in the reactor network of a FCLP for the MCC/MTAC RNS problems, such that (1) i1 and i2 have the same inlet/outlet concentration of component j, ∀ j ∈ Jx1 (i.e., CˆIj (i1) ) CˆIj (i2), COjˆ (i1) ) COjˆ (i2)); (2) i1 generates less or equal amounts of component j, ∀ j ∈ Jx2 (i.e., ∆Cj(i1) e ∆Cj(i2)); (3) i1 generates more or equal amounts of component j, ∀ j ∈ Jx3 (i.e., ∆Cj(i1) g ∆Cj(i2)); (4) i1 generates an equal amount of component j, ∀ j ∈ Jx4 (i.e., ∆Cj(i1) ) ∆Cj(i2)); and (5) cν(i1)τ(i1)eν g cν(i2)τ(i2)eν for the MCC problem; cµ(i1)τ(i1)eµ g cµ(i2)τ(i2)eµ and cˆ µ(i1) g cˆ µ(i2) for the MTAC problem, then the flowrate through reactor i1 can be set to zero at the global optimum of the FCLP for the corresponding problem. Proof. In proving this proposition, we compare an optimal state and its variation under which the reactor network is operated, using the variation variable δ. The optimal state variables satisfy constraints of the FCLP (eqs 70-79), which are given in Appendix B. Let us use [FˆI, FOI, FˆII, FOOˆ , FˆIOˆ ] and [F h ˆI, F h OI, F h ˆII, F h OOˆ , F h ˆIOˆ ] to represent the variables of the reactor network for the optimal state and its variation, respectively. We then proceed by contradiction. Let FˆI(i1) > 0. Consider the variation, F h ˆI(i1), of the flowrate ˆ I of the i1 reactor F (i1), such that
F h ˆI(i1) ) δFˆI(i1)
where δ ∈ [0, 1]
(38)
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3777
Let these perturbed variables be defined by the following equations:
F h )F OI
F h ˆI(i) ) FˆI(i)
OI
(39)
∀ i ) 1, ..., L; i * i1, i2
(40)
F h ˆI(i2) ) FˆI(i2) + (1 - δ)FˆI(i1)
(41)
F h ˆII(i, 1) ) FˆII(i, 1)
(42)
i * i 1 , i2
F h ˆII(i1, 1) ) δFˆII(i1, 1)
(43)
F h ˆII(i2, 1) ) FˆII(i2, 1) + (1 - δ)FˆII(i1, 1)
(44)
F h OOˆ (1, i) ) FOOˆ (1, i)
i * i1, i2
F h OOˆ (1, i2) ) FOOˆ (1, i1) + (1 - δ)FOOˆ (1, i1)
(47)
F h ˆIOˆ (i, j) ) FˆIOˆ (i, j) ∀ i ) 1, ..., L; i * i1, i2 ∀ j ) 1, ..., L; j * i1, i2 (48)
F h ˆIOˆ (i, i2) ) FˆIOˆ (i, i2) + (1 - δ)FˆIOˆ (i, i1) F h ˆIOˆ (i1, j) ) δFˆIOˆ (i1, j)
F h ˆIOˆ (i2, j) ) FˆIOˆ (i2, j) + (1 - δ)FˆIOˆ (i2, j)
(49)
∀i ) 1, ..., L (50)
∀j ) 1, ..., L
(51)
∀j ) 1, ..., L (52)
Straightforward algebraic manipulations establish that the variation, defined above through eqs 39-52, satisfies the constraints of the FCLP (eqs 70-79). For example, satisfaction of eq 76 is given in Appendix B. Let us first consider the MCC problem. The capital cost of the FCLP of the optimal state (νL) and its variation (νjL) are calculated as L
∑
ν)
cν(i)τ(i)eνFˆI(i)eν + cν(i2)τ(i2)eνFˆI(i1)eν +
i)1(i* i1,i2)
cν(i2)τ(i2)eνFˆI(i2)eν L
νjL )
∑
cν(i)τ(i)eνF h ˆI(i)eν + cν(i1)τ(i1)eνF h ˆI(i1)eν +
i)1(i*i1,i2)
cν(i2)τ(i2)eνF h ˆI(i2)eν L
)
∑
Since FˆI(i1) > 0, FˆI(i2) g 0, let z ) (FˆI(i1)/FˆI(i1) + FˆI(i2)), then z ∈ (0, 1]. As a result, when δ ) 0,
ν˜ L|δ)0 ) (FˆI(i1) + FˆI(i2))eν(cν(i1)τ(i1)eνzeν + cν(i2)τ(i2)eν(1 - z)eν - cν(i2)τ(i2)eν) Using the results in Appendix C and recalling that z ∈ (0, 1], the following can be established for the MCC problem:
if cν(i1)τ(i1)eν > cν(i2)τ(i2)eν, then ν˜ L|δ)0 > 0
(45) (46)
∀i ) 1, ..., L
cν(i1)τ(i1)eν(δFˆI(i1))eν - cν(i2)τ(i2)eν(FˆI(i2) + (1 - δ)FˆI(i1))eν
if cν(i1)τ(i1)eν ) cν(i2)τ(i2)eν, then
F h OOˆ (1, i1) ) δFOOˆ (1, i1)
F h ˆIOˆ (i, i1) ) δFˆIOˆ (i, i1)
ν˜ L ) cν(i1)τ(i1)eνFˆI(i1)eν + cν(i2)τ(i2)eνFˆI(i2)eν -
cν(i)τ(i)eνFˆI(i)eν + cν(i1)τ(i1)eν(δFˆI(i1))eν +
i)1(i*i1,i2)
{
ν˜ L|δ)0 ) 0 if FˆI(i2) ) 0 ν˜ L|δ)0 > 0 if FˆI(i2) > 0
It implies that, when cν(i1)τ(i1)eν > cν(i2)τ(i2)eν, the variation of δ ) 0 from FˆI(i1) > 0 (i.e., F h ˆI(i1) ) 0) is an improvement over the solution (i.e., νL > νjL). It is a contradiction to our assumption that FˆI(i1) is a variable of the optimal state. On the other hand, when cν(i1)τ(i1)eν ) cν(i2)τ(i2)eν, FˆI(i2) has to be zero at the optimal state to avoid the contradiction, if FˆI(i1) > 0 is assumed at the optimal state. However, it can be easily verified that the variation of δ ) 0 (i.e., F h ˆI(i2) ) FˆI(i1), F h ˆI(i1) ) FˆI(i2) ) 0) gives rise to the same objective value. Thus, FˆI(i1) can be set to zero at the optimum of the MCC problem. Similarly, for the MTAC problem, it can be established that
µ˜ L ) cµ(i1)τ(i1)eµFˆI(i1)eµ + cµ(i2)τ(i2)eµFˆI(i2)eµ + cˆ µ(i1)FˆI(i1) + cˆ µ(i2)FˆI(i2) - cµ(i1)τ(i1)eµ(δFˆI(i1))eµ cµ(i2)τ(i2)eµ(FˆI(i2) + (1 - δ)FˆI(i1))eµ - cˆ µ(i1)(δFˆI(i1)) cˆ µ(i2)(FˆI(i2) + (1 - δ)FˆI(i1)) and when δ ) 0,
µ˜ L|δ)0 ) (FˆI(i1) + FˆI(i2))eµ(cµ(i1)τ(i1)eµzeµ + cµ(i2)τ(i2)eµ(1 - z)eµ - cµ(i2)τ(i2)eµ) + FˆI(i1)(cˆ µ(i1) - cˆ µ(i2)) When cˆ µ(i1) g cˆ µ(i2), the second part of µ˜ L|δ)0 [i.e., FˆI(i1)(cˆ µ(i1) - cˆ µ(i2))] is always non-negative. Since the proof and arguments made for ν˜ L|δ)0 are valid for the first part of µ˜ L|δ)0, it can be concluded that FˆI(i1) can be set to zero at the optimum of the MTAC problem. 9 6. Case Study Consider the following Trambouze reaction scheme taking place in a homogeneous, isothermal reactor network, with a feed of pure reactant A with volumetric flowrate of 1.0 m3/s and CA0 ) 1.0 kmol/m3 (thus, NI ) NO ) 1). 0th order, k1 ) 0.025
A 98 B 1st order, k2 ) 0.2
cν(i2)τ(i2)eν(FˆI(i2) + (1 - δ)FˆI(i1))eν
A 98 C
where the first terms in νL and νjL represent capital costs associated with reactors excluding reactor i1 and i2, and they are not affected by perturbation and are identical in both equations. By defining ν˜ L as the difference between νL and νjL, then
A 98 D
2nd order, k3 ) 0.4
The reaction rates take the form
RB ) k 1
3778
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
RA ) -k1 - k2CA - k3CA2
{
τ)
RC ) k2CA RD ) k3CA
2
We look for the MCC of the reactor network when the conversion of A is greater than or equal to 30% and the component concentration of C is required to be less than or equal to 0.1 kmol/m3. The capital cost coefficient is c ) $104/ m1.8, and the economies of scale exponent eν ) 0.6. Before any optimization effort is made to solve the MCC problem, the AR is constructed to examine the feasibility of this problem. If there is no intersection between the AR and the specifications [i.e., (COA(1))upp ) 0.7, (COA(1))low ) 0, (COC (1))upp ) 0.1, and (COC (1))low ) 0], the MCC problem is infeasible. Figure 2 shows the AR in (CA, CC) concentration space, the CSTR/PFR trajectories, and the specified region which the network outlet concentrations have to belong to. Indeed, the AR contains the whole specified outlet region; thus, the MCC problem is feasible. By examining the reaction scheme, it can be seen that U, the set of reaction stoichiometric vectors, is a linearly independent set and S, the set of reaction rates, is also a linearly independent set. Thus, m ˜ b ) mb ) ms ) m ) 3, and nb ) 3. The following holds,
RA + RB + RC + RD ) 0
(
in -4k3(Cout A - CA )
in (2k3Cout A + k2) (2k3CA + k2) in ∆CC ) Cout C - CC ) (k2/2k3) 2 + k2τ + 2Cin A k3τ -k2τ + 2 ln 2
(
(
{
τ)
))
(
Thus, if all four components are of interest, one is able to conduct the optimization problem in a three-dimensional concentration space. In this particular problem, however, only components A and C are of interest, and Jz ) Jx ) {A,C}. Thus, u and y can be defined as T T T in out out u1 ( [F], u2 ( [Cin A CC CA λ] , y1 ( [FV] , y2 ( [CC τ]
{
The reactor model then becomes
( )
dCA ) -k1 - k2CA - k3CA2; CA|τj)0 ) dτj out Cin A , CA|τj)τ ) CA dCC out ) k2CA; CC|τj)0 ) Cin C , CC|τj)τ ) CC dτj
( )
{
out out 2 in Cout A - CA ) τ(-k1 - k2CA - k3(CA ) )
Cout C
-
Cin C
)τ
k2Cout A
}
}
)}
in Cout A - CA
(58)
if λ ) 1 (59)
V ) τF
(60)
The corresponding mathematical formulation of the reactor network synthesis problem featuring MCC can be rewritten as follows for a FCLP, L
νL ) inf(FˆI,FOI,FˆII,FOOˆ ,FˆIOˆ ) 104
(τ(i)0.6FˆI(i)0.6) ∑ i)1
(61)
L
1 ) FOI(1,1) +
FˆII(i, 1) ∑ i)1
(62)
L
1 ) FOI(1,1) +
which then implies
(54)
}
if λ ) 0
out 2 -k1 - k2Cout A - k3(CA ) in out ∆CC ) Cout C - CC ) τk2CA
(53)
CA + CB + CC + CD ) CA0
)
ˆI
FOOˆ (1, j) ∑ j)1
(63)
L
ˆII
F (i) ) F (i, 1) +
FˆIOˆ (i, j) ∑ j)1
FˆI(j) ) FOOˆ (1, j) +
FˆII(i, 1) ) CAˆI(i)FˆI(i) -
∀i ) 1, ..., L
(64)
L
FˆIOˆ (i, j) ∑ i)1
∀j ) 1, ..., L (65)
L
CAOˆ (j)FˆIOˆ (i, j) ∑ j)1
∀i ) 1, ..., L (66)
L
0 e FOI(1,1) + if λ ) 0
(55)
CAOˆ (j)FOOˆ (1, j) e 0.7 ∑ j)1
(67)
L
0e
if λ ) 1
V ) τF
(56) (57)
In addition, Jy ) {A} (i.e., generation rates of components A and C depend only on CA), and only lower/upper bounds at the network outlet are enforced for CC. As a result, the dimensionality can be reduced further, and the vectors u2 and y2 can be redefined as out T T u2 ( [Cin A CA λ] , y2 ( [∆CCτ]
Indeed, the aforementioned model can be solved analytically and gives rise to the following results:
∆CC(j)FˆI(j) e 0.1 ∑ j)1
FˆI g 0; FOI g 0; FˆII g 0; FOOˆ g 0; FˆIOˆ g0
(68) (69)
where L is the number of feasible units in the OP, and FI(1) ) FO(1) ) 1, CAI(1) ) 1, (CAO(1))upp ) 0.7, (CAO(1))low ) 0, CCI(1) ) 0, (CCO(1))upp ) 0.1, (CCO(1))low ) 0, V(i) ) τ(i)FˆI(i), and FOˆ (i) ) FˆI(i) have been substituted. The global optimality property presented in Section 5.2.2 is valid and employed in this study. Indeed, Jx1 ) {A}; Jx2 ) 0; Jx3 ) {C}; Jx4 ) 0. Thus, if there are two reactors (i1 and i2) in the reactor network of a FNLP, such that they have the same inlet/outlet component concentrations of A (i.e., CAˆI(i1) ) CAˆI(i2), CAOˆ (i1) ) CAOˆ (i2)), and i1 generates more or equal amounts of component C (i.e., ∆CC(i1) g ∆CC(i2)) and requires more or equal residence times (i.e., τ(i1) g τ(i2)), then the flowrate through reactor i1 can be set to zero at the optimum.
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3779
Figure 2. Feasibility check for the MCC problem of Trambouze reaction scheme.
Figure 3. IDEAS convergence for the MCC problem of Trambouze reaction scheme.
To demonstrate how the aforementioned dimensionalityreduction techniques can be used and their effectiveness, consider a discretization level of 16 subintervals. Hereafter, instead of the number of feasible reactors, the number of subintervals after discretization is used to represent the size of a FCLP (denote its objective value as νI). First, consider the following definition: in out out T T T u1 ( [F], u2 ( [Cin A CC CA λ] , y1 ( [FV] , y2 ( [CC τ]
Discretization is carried out such that the range of CA (i.e., [0.0, 1.0]) and the range of CC (i.e., [0.0, 1.0]) are divided uniformly to create 16 subintervals. It generates 17 × 17 ) 289 grid points, and after the feasible reactor generation stage, 4112 feasible CSTR/PFRs are generated. The resulting optimization formulation contains 4385 linear constraints and 932 655 flow variables. The FCLP is solved using the aforementioned branch-and-bound algorithm and gives rise to an objective of $1.5379 × 104.
Then consider the definition: T T T out u1 ( [F], u2 ( [Cin A CA λ] , y1 ( [FV] , y2 ( [∆CC τ]
Discretization is carried out such that the range of CA (i.e., [0.0, 1.0]) is divided uniformly to create 16 subintervals. It generates 272 feasible CSTRs/PFRs (i.e., reduced by 3840, 93%). The resulting optimization formulation contains 55 linear constraints and 561 flow variables. The solution gives rise to the same objective value. If Proposition 4 is used, the number of feasible reactors becomes 180 (i.e., reduced by 92, 34% more), and the number of linear constraints and variables in the resulting FCLP are 55 and 469, respectively. The same optimum objective is obtained again. From the comparison, it can be seen that the aforementioned dimensionality-reduction techniques are very effective. The numerical approximation result (νI) for the MCC problem is presented in Figure 3. The values of the capital cost presented
3780
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
RNS problems. Future work will focus on MCC/MTAC synthesis of nonisothermal reactor networks. Acknowledgment Figure 4. Optimal reactor network for the MCC problem of Trambouze reaction scheme.
in Figure 3 are the solutions to a series of FCLPs, each involving an increasing level of discretization, which in turn increases the number of available feasible reactors and, consequently, the size of the FCLP that is optimized. All runs are carried out by using the LP function in MINOS6.0 on a Linux machine with an AMD 2400 + MP CPU (running at 2 GHz clock speed). The termination condition employed is related to the percentage change in two consecutive values. The value used in this work is 5 × 10-3. The percentage changes between the last three FCLPs (νI)56 ) $1.4662 × 104, νI)84 ) $1.4569 × 104, and νI)114 ) $1.4546 × 104) are 6.34 × 10-3 and 1.58 × 10-3, respectively. The times to run the first νI)4 FCLP and the last νI)114 FCLP are 0.6 s and 14 days, respectively. It is worth pointing out that the branch-and-bound algorithm is naturally parallelizable. Every remaining branch can be checked simultaneously, and an upgraded lower bound can be applied immediately. Thus, if the program can be carried out on a 100node cluster, the time would be 0}j)1 intervals. Thus, after the branching step, the strong refining rule may divide Ψk into as many as 2m rectangular subsets, while the weak refining rule yields two subsets. In this work, the weak refining rule is employed. In general, the algorithm does not converge in a finite number of steps. The theorems about the convergence of the algorithm under different requirements on the problem functions are provided as follows, the proofs are referred to the original paper. Theorem A1. If φ is lower semicontinuous, Ω is closed, and the strong refining rule is used, then any limit point of {xk} is a solution of problem P. Theorem A2. If φ is continuous, Ω is closed, and the weak refining rule is used, then any limit point of {xk} is a solution of problem P. Theorem A3. If φ is concave, Ω is a linear polyhedron, and the weak refining rule is used, then the sequence {xk} is finite and terminates in a solution of problem P. Appendix B The feasible region of a FCLP with different types of concentration bounds is defined as follows:
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008 3781
Thus,
L
FI(1) ) FOI(1,1) +
FˆII(i, 1) ∑ i)1
(70) L
L
FO(1) ) FOI(1,1) + FˆI(i) ) FˆII(i, 1) +
FOOˆ (1, j) ∑ j)1
∀i ) 1, ..., L (72)
FˆIOˆ (i, j) ∑ i)1
∀j ) 1, ..., L (73)
L
COkˆ (j)FˆIOˆ (i, j) ∑ j)1
L
1) +
e
CIk(1)FOI(1,
e
(COk (1))uppFO(1)
L
e
∆Ck(j)F (j) + ∑ j)1 ˆI
∑ j)1
COAˆ (j)FOOˆ (1,
j)
∀ k ∈ J x1
(75)
CIk(1)FI(1)
∀ k ∈ J x4
FˆI g0; FOI g 0; FˆII g 0; FOOˆ g 0; FˆIOˆ g 0
where FI(1) ) FO(1), and substitutions of V(i) ) τ(i)FˆI(i) and FOˆ (i) ) FˆI(i) are used. Satisfaction of eq 76 for perturbed variables is as follows:
∑ j)1
L
∑
∆Ck(j)F h ˆI(j) )
(∆Ck(j)F h ˆI(j)) +
j)1(j*i1,i2)
∆Ck(i1)F h ˆI(i1) + ∆Ck(i2)F h ˆI(i2) L
∑
)
(∆Ck(j)FˆI(j)) + ∆Ck(i1) δFˆI(i1) +
j)1(j*i1,i2)
∆Ck(i2)(FˆI(i2) + (1 - δ)FˆI(i1)) L
)
(∆Ck(j)FˆI(j)) + FˆI(i1)(1 - δ)(∆Ck(i2) - ∆Ck(i1)) ∑ j)1 L
g
(∆Ck(j)FˆI(j)) ∑ j)1
∀k ∈ Jx2
where e is a constant and e ∈ (0, 1) (i.e., e ) eν for the MCC problem, e ) eµ for the MTAC problem). It can be easily carried out and yields the following:
{
a2, and jz * ) 0 if a1 > a2 min (zj ∈[0,1]) f(zj) ) f(zj *) ) a2, and jz * ) 0, 1 if a1 ) a2 a1, and jz * ) 1 if a1 < a2 Nomenclature
IDEAS Variables
L
∆Ck(j)FˆI(j) + CIk(1)FI(1) ∑ j)1
e (CkO(1))uppFO(1)
min(zj ∈[0,1]) f(zj) ) a1jze + a2(1 - jz)e
(79)
(76)
∀ k ∈ Jx3 (77)
L
Appendix C
(78)
∀ k ∈ J x2
∆Ck(j)FˆII(j) + CIk(1)FI(1) e (COk (1))uppFO(1) ∑ j)1
(COk (1))lowFO(1) e
∀ k ∈ J x2
Φ ) input-output information map representing a process model u ) input vector of an information map y ) output vector of an information map rj ) jth reaction rate (kmol/(m3 s)) Ri ) ith component’s generation rate F ) volumetric flowrate (m3/s) V ) reactor volume (m3) τ ) reactor residence time (s) Ci ) ith component’s molar concentration (kmol/m3) Cin i ) ith component’s inlet molar concentration for a unit Cout i ) ith component’s outlet molar concentration for a unit C0i ) ith component’s initial molar concentration entering a reactor network k ) reaction rate constants
L
0e
CIk(1)FI(1) g (COk (1))lowFO(1)
Consider the following optimization problem,
∀i ) 1, ..., L ∀ k ∈ Jx1 (74)
(COk (1))lowFO(1)
(∆Ck(j)FˆI(j)) + ∑ j)1
L
CkˆI (i)FˆI(i) ) CIk(1)FˆII(i, 1) +
(COk (1))lowFO(1)
L
L
FˆIOˆ (i, j) ∑ j)1
FˆI(j) ) FOOˆ (1, j) +
(71)
∑ j)1
∆Ck(j)F h ˆI(j) + CIk(1)FI(1) g
NI ) number of network inlet streams NO ) number of network outlet streams CIk(j) ) kth component concentration in the jth network inlet k ) 1, n; j ) 1, NI COk (i) ) kth component concentration in the ith network outlet k ) 1, n; i ) 1, NO CˆIk(i) ) kth component concentration in the ith OP inlet k ) 1, n; i ) 1, ∞ COkˆ (j) ) kth component concentration in the jth OP outlet k ) 1, n; j ) 1, ∞ FI(j) ) jth network inlet flowrate j ) 1, NI FO(i) ) ith network outlet flowrate i ) 1, NO FˆI(i) ) ith OP inlet flowrate i ) 1, ∞ FOˆ (j) ) jth OP outlet flowrate i ) 1, ∞ FOI(ij) ) flowrate from the jth network inlet to the ith network outlet i ) 1, NO; j ) 1, NI FˆII(ij) ) flowrate from the jth network inlet to the ith OP inlet i ) 1, ∞; j ) 1, NI FOOˆ (ij) ) flowrate from the jth OP outlet to the ith network outlet i ) 1, NO; j ) 1, ∞ FˆIOˆ (ij) ) flowrate from the jth OP outlet to the ith OP inlet i ) 1, ∞; j ) 1, ∞ V(i) ) volume of the ith OP unit (reactor) i ) 1, ∞ τ(i) ) residence time of the ith OP unit (reactor) i ) 1, ∞
3782
Ind. Eng. Chem. Res., Vol. 47, No. 10, 2008
λ(i) ) technology flag of the ith OP unit (reactor) (i.e., 0 for PFR, 1 for CSTR), i ) 1, ∞ u(i) ) input of the ith OP unit (reactor) information map i ) 1, ∞ y(i) ) output of the ith OP unit (reactor) information map i ) 1, ∞ Index Subsets Jn ) whole component index set { j, j ) 1, ..., n} Jx ) subset of Jn containing the indices of components of interest Jz ) least subset: it contains Jx; Rj(∀ j ∈ Jz) depend solely on {Ci}i∈Jz Jy ) subset of Jz containing the indices of concentrations on which Rj (∀ j ∈ Jz) depend Jx1 ) {j:j ∈ Jx, j ∈ Jy} Jx2 ) {j:j ∈ Jx, j ∉ Jy, (COj (1))low g 0, (COj (1))upp ) ∞} Jx3 ) {j:j ∈ Jx, j ∉ Jy, (COj (1))low ) 0, 0 < (COj (1))upp < ∞} Jx4 ) {j:j ∈ Jx, j ∉ Jx1, Jx2, Jx3} Note Added after ASAP Publication. Because of production errors, the version of this paper that was published on the Web 4/5/2008 had errors, primarily involving some superscripts and subscripts throughout the paper. The paper was corrected and reposted to the Web 4/28/2008, and the final corrected version of this paper was reposted to the Web 5/14/2008. Literature Cited (1) Achenie, L.; Biegler, L. T. A Superstructure Based Approach to Chemical Reactor Network Synthesis. Comput. Chem. Eng. 1990, 14, 2340. (2) Kokossis, A. C.; Floudas, C. A. Optimization of Complex Reactor Networks:. Isothermal Operation. Chem. Eng. Sci. 1990, 45, 595-614. (3) Kokossis, A. C.; Floudas, C. A. Optimization of Complex Reactor Networks:. Nonisothermal Operation. Chem. Eng. Sci. 1994, 49, 10371051. (4) Schweiger, C. A.; Floudas, C. A. Optimization Framework for the Synthesis of Chemical Reactor Networks. Ind. Eng. Chem. Res. 1999, 38, 744-766. (5) Mehta, V. L.; Kokossis, A. C. Nonisothermal Synthesis of Homogeneous and Multiphase Reactor Networks. AIChE J. 2000, 46, 22562273. (6) Kokossis, A. C.; Floudas, C. A. Synthesis of Isothermal ReactorSeperator-Recycle Systems. Chem. Eng. Sci. 1991, 46, 1361-1383. (7) Luyben, M. L.; Luyben, W. L. Design and Control of a Complex Process Involving Two Reaction Steps, Three Distillation Columns, and Two Recycle Streams. Ind. Eng. Chem. Res. 1995, 34, 3885-3898. (8) Bedenik, N. I.; Pahor, B.; Kravanja, Z. An Integrated Strategy for the Hierarchical Multilevel MINLP Synthesis of Overall Process Flowsheets
Using the Combined Synthesis/Analysis Approach. Comput. Chem. Eng. 2004, 28, 693-706. (9) Esposito, W. R.; Floudas, C. A. Deterministic Global Optimization in Isothermal Reactor Network Synthesis. J. Global Optimization 2002, 22, 59-95. (10) Glasser, D.; Hildebrandt, D.; Crowe, C. A Geometric Approach to Steady Flow Reactors: The Attainable Region and Optimization in Concentration Space. Ind. Eng. Chem. Res. 1987, 26, 18031810. (11) Hildebrandt, D.; Glasser, D.; Crowe, C. Geometry of the Attainable Region Generated by Reaction and Mixing: With and Without Constraints. Ind. Eng. Chem. Res. 1990, 29, 49-58. (12) Wilson, S.; Manousiouthakis, V. IDEAS Approach to Process Network Synthesis: Application to Multicomponent MEN. AIChE J. 2000, 46, 2408-2416. (13) Drake, J. E.; Manousiouthakis, V. I. IDEAS Approach to Process Network Synthesis: Minimum Utility Cost for Complex Distillation Networks. Chem. Eng. Sci. 2002, 57, 3095-3106. (14) Drake, J. E.; Manousiouthakis, V. I. IDEAS Approach to Process Network Synthesis: Minimum Plate Area for Complex Distillation Networks with Fixed Utility Cost. Ind. Eng. Chem. Res. 2002, 41, 4984-4992. (15) Holiastos, K.; Manousiouthakis, V. I. Infinite-dimensional statespace (IDEAS) Approach to Globally Optimal Design of Distillation Networks Featuring Heat and Power Integration. Ind. Eng. Chem. Res. 2004, 43, 7826-7842. (16) Burri, J. F.; Wilson, S. D.; Manousiouthakis, V. I. Infinite DimEnsionAl State-space (IDEAS) Approach to Reactor Network Synthesis: Application to Attainable Region Construction. Comput. Chem. Eng. 2002, 26, 849-862. (17) Zhou, W.; Manousiouthakis, V. I. Variable Density Fluid Reactor Network SynthesissConstruction of the Attainable Region through the IDEAS Approach. Chem. Eng. J. 2007, 129, 91-103. (18) Burri, J. F.; Manousiouthakis, V. I. Global Optimization of Reactive Distillation Networks Using IDEAS. Comput. Chem. Eng. 2004, 28, 25092521. (19) Justanieah, A. M.; Manousiouthakis, V. I. IDEAS Approach to the Synthesis of Globally Optimal Separation Networks: Application to Chromium Recovery from Wastewater. AdV. EnViron. Res. 2007, 549562. (20) Martin, L.; Manousiouthakis, V. I. Globally Optimal Power Cycle Synthesis via the Infinite-DimEnsionAl State-space (IDEAS) Approach Featuring Minimum Area with Fixed Utility. Chem. Eng. Sci. 2003, 58, 4291-4305. (21) Brown, T. R. Capital Cost Estimating. Hydrocarbon Process. 2000, 79, 93-100. (22) Falk, J. E.; Soland, R. M. An Algorithm for Separable Nonconvex Programming Problems. Manage. Sci. 1969, 15, 550-569. (23) Zhou, W.; Manousiouthakis, V. I. On Dimensionality of Attainable Region Construction for Isothermal Reactor Networks. Comput. Chem. Eng. 2008, 32 (3), 439-450.
ReceiVed for reView May 23, 2006 ReVised manuscript receiVed January 2, 2008 Accepted January 3, 2008 IE060653+
