Synthesis of Optimal Dynamic Process Systems by a Gradient Method Naonori Nishida*’ Science University of Tokyo, Shinluku-ku, Tokyo, Japan
Atsunobu lchikawa Tokyo Institute o f Technology, Meguro-ku, Tokyo, Japan
T h e structure parameter approach for the synthesis of dynamic process systems where t h e behavior of a processing unit is described by a set of ordinary differential equations is presented. By means of ordinary differential calculus, the necessary condition for t h e optimal system is derived which requires that the structure parameters and t h e design vectors satisfy the weak maximum condition of an integral of Hamiltonian-type scalar function and that t h e control vectors satisfy t h e maximum condition for t h e Hamiltonian-type scalar function. T h e gradient method in function space for attaining the optimal structure is developed on the basis of the necessary condition. To illustrate the use of the gradient method, synthesis of t h e optimal semibatch reaction system is conducted.
Introduction Process system synthesis is an act of determining the optimal interconnection of processing units as well as the optimal type and design of the units within a process system. The interconnection of processing units is called the structure of the process system. The following nature of process systems is of particular interest in connection with the system synthesis problem. When a performance of processing units and a structure of the system are specified, the performance of the system would be uniquely characterized, but not vice versa; that is, when the performance of the system is specified a performance of the processing units and the structure of the system are not to be determined uniquely. Out of this large number of alternatives, an optimal system is to be selected depending upon suitable criteria of the system. Several works on the consideration of the optimal synthesis of the chemical processing systems have been done by, e.g., Rudd and his coworkers (1969, 1970, 1971, 1972), Westerberg et al. (1972, 1973), King et al. (1972), and the authors‘ group (1971, 1972, 1973). These papers, however, concern only the optimal synthesis of steady-state process systems. In this paper, the synthesis of optimal dynamic process systems in which the state of the system changes with time is discussed. A wide variety of important industrial problems fall into the dynamic process system, such as the start-up of equipment, batch or semibatch operation of processing units, the change from one set of operating conditions to another, and the perturbations which develop as process conditions fluctuate. Furthermore, Douglas and Rippin (1966) demonstrated that certain types of processing units may give better performance by oscillatory rather than steady-state operation. It is therefore of interest to take into account of the effect of the dynamic behavior of a process system at the stage of process system synthesis rather than merely synthesizing the system assuming steady-state operation and then compensating the effect of dynamics by means of process control. The purpose of this paper is to develop a method of synthesizing a process system where the dynamic operation is inevitable or preferable. The structure parameter approach, which has been already introduced by the authors (Ichikawa et al., 19691, is employed to develop a
’ Address
Chemical 36830.
correspondence to this author c / o Professor Y . A. Liu. Engineering Department, Auburn University, Auburn, Ala.
236 Ind. Eng. Chem., Process Des. Dev., Vol. 14, No. 3, 1975
synthesis method of process system where behavior of the processing units is described by a set of ordinary differential equations. The problem will be formulated by introducing the structure parameter. The necessary condition for the optimal dynamic system is given on the basis of ordinary differential calculus. An iterative gradient method in function space is developed to reduce the difficulty in numerical computation that may arise in the solution of a two-point boundary-value problem. Statement of the Problem Consider a processing unit whose dynamic behavior is represented by a set of ordinary differential equations
y,(O) = yio (i = 1, 2,
. . ., 1V)
where x , t R n and y l t R n are, respectively, the input stream state vector and the output stream state vector expressed in terms of extensive variables, d, the design vector, u, the control vector, Q L C Rs the admissible region of the design, TL.C R’ the admissible region of the control, g , t R n the vector valued function, t c [ O , r ] the time variable, T the final time of process operation, and N the total number of subsystems. We shall limit our consideration to the oneinput one-output system where each subsystem as well as the system has only one input stream and one output stream. Extension to the multi-input multi-output system such as the separation units will be discussed in Appendix. The design vector is assumed to remain constant during the period of process operations. When the design vector is within the preassigned admissible region Q,, it is called the admissible design. The control varies with time so that optimal control is attained. When the control is piecewise continuous during the period of operation [0, r] and is within the preassigned region Tl,it is called admissible control. It is convenient to introduce an imaginary input subsystem and an imaginary output subsystem to account for the system-input and the system-output, respectively. The input subsystem is to receive the system-input and the output subsystem is to deliver the system-output. No transformations of stream vectors occur in the imaginary input and output sybsystems, that is
..
yi(t) = x , ( t ) (i = 1, 2 , . , 1v) (2 1 Associated with the subsystem we have a specified objective function which is a scalar function of the input and
output stream state vectors, the design vector, and the control vector. Ji
= L T v i i . ( x i , y i , d i , U i+ )
f i ” ( y i , ~ . ~ ) } d (it = 1, 2,
.. ., N )
(3)
where f,’ accounts for operating and investment cost of subsystem i; f L ” accounts for that of the output streams going other subsystems from subsystem i. The objective function of the imaginary input or output subsystem may be used to represent the objective function associated with the input or output stream of the system, respectively, such as the cost of the raw material stream or the sale of the product stream. The complete structure of the process flow scheme can be organized by joining the outputs of the unit to the inputs of other units in order to form an interconnected structure. By virtue of the structure parameters introduced in the previous papers (Ichikawa et al., 1972, 1973, 1974), the state of the input stream of a subsystem is expressed as follows N
cuijYj(t) (i = 1, 2,
x,(t) = j=i
. .., 147)
5
( 2 , j = 1, 2,
Qij
culj =
o
( j = 1, 2,
..., IV) . . ., iV)
N
T
{ f i ’ ( x i , y i , d i , u i )+ fi”(yi, a . i ) } d t (6) On the basis of the above definition of the system, the problem of the optimal synthesis of a dynamic process system can be stated as follows: “Find the optimal admissible structure, design. and control, so that the system objective function eq 6 is minimum under the constraints eq 1,4, and 5.” i.1
. . . , N)
The admissible structure and the admissible design are said to satisfy the weak maximum condition, if for given values of x,,y L ,A , , and p r , the functions
and
=;
J =g J i = i=l
X,(T) = 0 (i = 1, 2,
(5)
0 ( i = 1, 2, ..., iV) In essence eq 5 specifies the admissible region of the structural parameters. The structure parameters that satisfy eq 5 are called the admissible structure. It is assumed that the admissible structure is determined at the stage of process system synthesis and remains constant during the process operation. Notice that the interconnection given by eq 4 is conceptual and does not necessarily correspond to any real part or process unit of the system such as a mixing vessel or a pipe line. If it corresponds to an actual part or process unit, it should be regarded as a subsystem. There is no need, therefore, to consider the interconnection having dynamic behavior. For convenience, let us define all of a L ,(i, j = 1, 2, . . . , N ) by a matrix a = {a,Iand LY.,= (cyp, az,, . . . ,a.v,). The objective function of system is to be described as algebraic addition of the objection functions associated with the subsystems. CYiN
The associated adjoint equations are
(4)
where cy,, (i, j = 1, 2, , . . , N) are scalar-value variables called structure parameters and satisfy the constraints
0
vectors are uniquely determined when the input stream state vectors, the design vectors, and the control vectors are given. The control vectors are piecewise continuous in the interval [0, r ] . To obtain a necessary condition for the optimal structure of a system, the scalar functions analogous to the Hamiltonian in the theory of the maximum principle are defined
0
Necessary Condition for the Optimal Dynamic System Knowledge of the following assumptions of the problem are indispensable for dealing with the problem by means of variational calculus. The functions f,‘, f,”, and g, (i = 1, 2 , . . . , N ) are twice continuously differentiable and satisfy the Lipschitz conditions. The output stream state
are either stationary with respect to the structure and the design that lie in the interior of the respective admissible regions or attain their maximum when the admissible structure and/or the admissible design lie on the boundaries of the admissible regions. The admissible control is said to satisfy the maximum condition if at any instant [0, TI the function H , attains its maximum value for given values of x,, y,, A[, and p [. Then the necessary condition for the optimal dynamic system is given by the following theorem. Theorem 1. If a dynamic system is optimal, the admissible designs satisfy the weak maximum condition and the admissible control satisfies the maximum condition. Theorem 1 and the succeeding theorem 2 are direct extensions of analogous theorems of a previous paper on the synthesis of steady-state systems by Ichikawa and Fan (1973). The proofs for these theorems are straightforward and need not be detailed, but are discussed elsewhere (Nishida, 1972). Let us here investigate a special case where the objective functions associated with output streams of subsystems, f,” ( i = l L 2 , . . . , N),equal zero. Then, the weak maximization of SLwith respect to N .,, that is
1, 0
5 CYki
5
l}
(i = 1, 2,
. . . , NJ
(12)
is of particular importance in the determination of an optimal system structure. wmaxlA; B, C, . . .) stands for the weak maximization of A under the constraints of B, C, . . . The solution of eq 12 is given as Ind. Eng. Chem.. Process Des. Dev., Vol. 14, No. 3, 1975
237
N
0
5
aki 5 1 a n d x a k i = 1 if k E K i
(13)
k=l
aki = O if k @ K i (i = 1, 2,
where e-, is the positive-valued step-size factor. Step 6. Return to Step 2 and repeat the procedure until
. . . , N)
where
J‘R+”
Mki =
and
lT
pkTyidt (i = 1, 2,
. . .,
(14)
Ki = { k , Mki = max{Mfi}} 1
This leads to the following theorem. Theorem 2. If a stream, say y,, is fed to the input stream XJ in the optimal structure, then the value of M,, is the maximum among all MkL ( k = 1, 2, . . , , N). If split of a stream y, occurs, then M R , ( k E K ) have the same maximum value, where K , is a set of the subsystem numbers which are incident to the stream y,. Direct use of the necessary condition for the optimal dynamic system is, however, not very promising. Both the weak maximum condition and the maximum condition require optimal values of the adjoint vectors which are the solutions of the adjoint equations. The adjoint equations which are to be in turn involve the optimal values of CY,, determined. This problem is a two-point boundary-value problem formulated by the system equations (l),(2), ( 4 ) , associated by the adjoint equation (9) and the maximum condition. The solution to the problem will be very difficult to find if we note the high dimensionality and nonlinearity of chemical processing systems. Therefore, an efficient algorithm to find the optimal value of as well as the design and control variables has to be developed in order to reduce the computational difficulty.
Gradient Method for Attaining the Optimal Dynamic System The iterative procedure considered emphasizes the maximization of the functions (i = 1, 2 , . . . , N). The gradient direction of S, with respect to the structure parameters can be easily estimated as
s,
s,
where V, ,S, is the gradient of in the structure parameter space. Now, let us outline the overall iterative procedure. First, it is noted that eq 15 specifies that a change in CY., in the next iteration is determined based on the gradient of evaluated at the previous set of structure parameters. Thus, an iterative procedure of the dynamic version of the gradient method can be suggested as follows. Step 1. Assume initial estimates, a t R ) ,k = 0, of the admissible structure parameters which satisfy the constraint, eq 5. Step 2. Carry out optimization over the restricted structure of a t R ) and find the optimal values of designs, d,, controls, u, and state vectors, x, and y,, (i = 1, 2, . . . ,
s,,
J‘R’.
Step 2 of the above procedures may be carried out by appropriate hill climbing methods, such as the dynamic version of the gradient methods (Bryson et al., 1962, Kelley et al., 1963) and the conjugate gradient method (Lasdon et al., 1967). The proposed algorithm based on the gradient method will obviously yield only a local solution to the problem, if it exists. In order to obtain a global optimum, the iterative optimizations have to be tried after changing the initial point of hill climbing.
Synthesis of Semi-Batch Reaction System To illustrate the use of the gradient method for the optimal system, synthesis of a simple reaction system is conducted as an example. Two types of the reaction are considered here. One is an irreversible isothermal reaction K 1
A-B
and the other a reversible nonisothermal one K1
A
e
B
KZ
where A is the raw material and B the desired product. On account of space consideration, the problem will be formulated for the reversible system hereafter since such formulation will cover the irreversible reaction by taking the reaction rate K2 = 0. Operation of both the systems is to be carried out in the following manner. Initially the entire reaction system is filled with an inert fluid, C, and at time t = 0 operation of the system is initiated by continuously supplying reactant A. The operation is continued for a certain time period, [0, r ] . The performance of the system is measured by an integral of revenue from the product minus cost of the operation for this time period. Processing units available for the system are continuous stirred tank reactors (CSTR) and ideal separators (SEP). Among possible combinations of the available process units, we choose a restricted structure which is used to find the optimal structure parameters, the design and/or control variables, as shown in Figure 1. Although the other more complicated combination than that may be considered, the illustrated combination will generate all the typical cases such as the sequential type of system, the parallel type, and so on. The behavior of two CSTR’s is described by the following set of equations.
N).
Step 3. Evaluate the objective functionJ‘R),eq 6. Step 4. Compute the adjoint vectors A , and b L (i = 1, 2, . . . , N) through the adjoint equations using the above optimal design and control variables of the restricted structure. Step 5. Modify the ( k + 1)th trial of the structure parameters, a t R l’, by a formula +
238 Ind. Eng. Chem., Process Des. Dev., Vol. 14, No. 3, 1975
where V is the reactor volume, C the concentration of a material, F the volumetric flow rate of the reaction mix-
ture, K 1 and K z the reaction rate constants and u:! and u3 the temperatures of the CSTR's. The superscript i stands for the input stream, the subscripts A, B, and C of components in the reaction mixture, and the overall flow rates of the reaction mixtures are expressed in terms of molar flow rates of the components as -
F =
+
(XiA
xiBs
+
+
+
( X ~ A
X ~ B
xic)p
=
Y4l
YiB
=
1
-
-(XiA ' i
1 -(xiB ' i
-
+
YiB)
-
YjA)
Xic)p
(KlYiA
(K1YiA
= -(xjc
-
=
=
'4A2
YjC)
(19)
b ( X 4 ~+ X ~ Bf
XdB1
=
Vi/(XiA
XiB
+
xic)
(i = 2, 3)
is the mean residence time of the reactor. The objective function associated with a CSTR is given in the form of operating cost as
( i = 2, 3) (20) f i = avi = a ' i P ( x i A f XiB + XiC) where a is the cost factor of the CSTR. The adjoint equations for the CSTR, defined by eq 9, are written below for the convenience of carrying out the synthesis.
p3Aa34
2
+
p2Acy242
(25)
Y N B ON4
+
13Ccy342
+
p2Ccy242
We shall define the imaginary input and output subsystems to account for the system input, system output, and their associated objective functions. Thus = x1;h = 0
Yi
=
(24)
x4c)P
1
=
p(NtI)Cff(N+1)42
where ' i
+
p\Nt1)Aa(Nt1)42
- K~Y~B (i) = 2, 3)
= Y ~ B ' y, i c ( 0 ) = y i c o
3.'iAo,3'iB(o)
(23)
$4
where b is the cost factor of the separator. The adjoint equations for the separator are
K z y i ~ )
' i
YIA(0)
=
;]x4
f4
'4C2
1
YiC
-
y:
where y4l and y42 are, respectively, output stream 1 which contains only B and output stream 2 which contains A and C. The cost of the separator is given by
where x and y are the molar flow rates of the components, the molar volume which is assumed to be the same for all the components and does not change during reaction. Substituting eq 18 into eq 17 =
1 0 0
[; ;
p
YiA
[; ;
0 0 0
Xi
I - -
cis
The ideal separator removes material B completely from the rest of reaction mixture. Dynamic lag of the separator is considered negligible when compared to that of the CSTR's. Thus, the system equations for the ideal separator are
(26)
yN = xN l* fN =
-VTyN
(27)
=
-VNtlTYNtl
(28)
=
(29)
YN+l
INT
xN+l; f N + l
=
VNT;
=
IN+IT
vN+lT
where v is the price vector associated with the product stream. Now the synthesis problem is to find the system structure as well as the design and control of processing units so as to minimize the total objective function T N+l
J =
(30)
xfidt i=l
0
= o
The following numerical values are used throughout the example: a = 1.0, b = 0.04, p = 1.0, V N = (0.0, 1.0, O.O)T, v,y+I = (0.0, 0.0, O.O)T, T = 1.0, XI ( t L 0) = (1.0,0.0, O.O)T, x1 i t < 0 ) = (0.0, 0.0, l . O ) T , yLo = (0.0, 0.0, l.O)T (i = 2,3 ) . The following scalar functions defined by eq 8 will be used to find the gradient directions with respect to the structure parameters.
=o =o = o
3
st
+ pkB(YklYIB)
(PkAOk 1YiA k=2
s2
=
5
s3
=
p4Aff43Y3A
4
= o '3,
-
1 +
Q3C
P4Cff43
= 0
(i = 2, 3) 1 P2s
=
[ A 2s ( xZ A +
Y ~ A )+ % B ( x ? B
p3s
=
1
XZB
-
v3
xzc) +
Y ~ B )+
- [ & 3 ~ ( ~ 3 A + x3B
Y3A)
+
+
%A(x~A
I
-
Xzc(x2c - ~ 2 c ) l
x3C)
+
'3A(X3A
+ ' 3 B ( X 3 B - 3'38) + ' 3 C ( x 3 C
-
(22)
3
s4' =
+
(pkAffk2Y2A
+ 2
p4Bcy43Y3B
(pkAak2 Y4A k=2
2
+ pkCcyk2y?C)
pkBC'k2Y2B
+ p4ca43Y3c
+ p k C cyk42Y4C2) + 2
p ( N t 1)AaINt1)4 Y4A2
+
2
2
p ( N t l ) C a ( N + 1 ) 4 Y4C
Two types of dynamic synthesis problem are examined here. Problem A. The first problem is to find the optimal structure parameters together with the optimal design variables of two CSTR's, that is, V:! and Va, so as to minimize the objective function J, eq 30. The reaction system taken here is an irreversible isothermal one
- Y3C)l
Ind. Eng. Chem., Process Des. Dev., Vol. 14, No. 3, 1975
239
0.8
1
a54
Figure 1. System structure for the example. 0.b 2.0
NUMBER OF TRIALS
1.8
Figure 4 . Convergence rates of structure parameters of Problem A.
1.6 wl.4 + w
61.2 = 21.0 LL q0.8 2 g0.6
580
7
570
0.u
560
0.2
0.0 00
02
16
O L
550
C8
t TIWE
Y
3 z
Figure 2. Optimal profiles of state vectors of Problem A.
540 a r w c
530
n
520
-0.44
510 500
0.0
0.2
0.4
0,6
0.8
1.0
t,TlME
Figure 5. Optimal temperature controls of Problem B. -3 52
5
13 15 20 25 NUYBEQ OF TI;IALS
30
Figure 3. Convergence rate of the objective function of Problem A.
where the value of reaction rate constant K1 used in the computation is 10.0 and K2 = 0.0 in eq 17. Problem B. The second problem is to find the optimal system such that the structure parameters and the temperatures of two C S T R s minimize the objective function J, eq 30. The reaction system considered here is reversible one K1
A--B K2
where the numerical values of reaction rate constant used in the computation are K1 = 0.15 X lo6 exp(-lO,OOO.O/Rui), K2 = 0.10 x 1010 exp (-2O,OOO.O/Rui) (i = 2, 3, and R = 2.0). Both the reactor volumes of CSTRs, V 2 and V3, are 0.1 and they are not the decision variables in this problem. The following constraint is imposed on the control variables u i ( t ) 5 573°K f o r all 240
tc[O, T ] (i = 2, 3)
Ind. Eng. Chem., Process Des. Dev.. Vol. 14, No. 3, 1975
Computed Results Problem A. The computed optimal structure parameters are as follows: a 2 1 = 1.0, ( ~ 3 1= 0.0, a 3 2 = 1.0, a 4 2 = 0.0, (Y2d2 = 0.9305, ( ~ 3 = 4 ~0.0, ( Y , , Y + I = ~ ~0.0695. ~ The optimal values of the objective function and the design variables are, respectively, J = -0.5026, VZ = 0.0732, and V3 = 0.0736. The corresponding optimal values of state vectors are shown in Figure 2. To illustrate the rate of convergence, the convergence rates in the objective function and the structure parameters for (Y242 and (Y(,v+1)4*are shown in Figures 3 and 4, respectively. The slow convergence after 18 trials is due to the singularity of the structure parameters for (Y242 and 0 ( ( ~ V + l i 4 2 .In Step 5 of the algorithm, the step size factor used is constant at each iteration (c- = 0.15). Problem B. The computed optimal values of the objective function and structure parameters aie J = -0.4852, 0 2 1 = 1.0, (Y31 = 0.0, e 3 2 = 1.0, Cy42 = 0.0, 0(242 = 1.0, (Y 342 = 0.0, and ( Y ( N + 1 ) 4 ' = 0.0. The time variations in the optimal temperature controls, u2(t) and u 3 ( t ) , and state vectors are, respectively, shown in Figures 5 and 6 . The convergence rate of the objective function and the structure parameters are shown in Figures 7 and 8. As can be seen from the figure, eighteen trials have been allowed before all the structure parameters converge to zero or one. Figure 9 summarizes the optimal structure of Problems A and B schematically.
2
3
C S T R
C S T R
(0.9305))(O.Oh951 I
Ca) O P T I I A L STRUCTURE OF PROBLEM A .
3 -
I
06
c~
08
13
t,iIVE
Figure 6. Optimal profiles of state vectors of Problem B
7
O’I
,
1
b) O P T I M A L
STRUCTURE OF PROBLEM
B
Figure 9. Optimal structures of illustrative examples. -0.1
-0.2
Appendix: Multi-Input and Multi-Output System We shall consider that the behavior of multi-input and multi-output subsystems is described by a set of differential equations
‘ 1
-0.3
-0.4
2
yi = &Q(Xil,Xi, . . . ,
-0.5
0
2
4
E
8 10 12 14 16 18 20 NUYEER OF T R I A L S
Figure 7 . Convergence rate of the objective function of Problem B.
. . . ,yiQi,di,ui)Y,”(o) = Y t O , (i = I, 2, . . . , N ; = 1, 2, . . . , Qi)
xipi,yii,yi2,
(AI)
and the objective function of subsystem i is given by
where the superscripts indicate the stream number, and
Pi and Qi are the number of input streams and that of output streams, respectively. The structure of the system is expressed by 9
2
4
6
8 10 12 14 16 18 20 OF T R I A L S
YlJMBER
Figure 8. Convergence rates of structure parameters of Problem B.
Conclusion The structure parameter approach is presented for finding the optimal structure of dynamic process systems where the behavior of a processing unit is described by a set of ordinary differential equations. By means of variational calculus, the necessary condition for the optimal dynamic system is derived to develop an algorithm to attain the optimal structure. The application of this algorithm is in the design and control of a semi-batch reaction system. Numerical results on this application show that the algorithm gives rise to reasonable structure of a dynamic process system with stable convergence.
where cq,pq is the structure parameter defined as a portion of the qth output stream of subsystem j which is fed to the p t h input stream of subsystem i. The scalar functions are defined as
( i = 1, 2,
. . ., N ; q
= 1, 2,
. . .,
Q i ) (A4)
where the adjoint vectors are defined by the adjoint equations Ind. Eng. Chem., Process Des. Dev., Vol. 14, No. 3, 1975
241
0 = holding time of reactor X,p = adjoint vectors defined by eq 9 T = admissible region of control vector D = admissible region of design vector
( i = 1, 2,
T
. . ., N ;q
. . ., N; p 2, . . ., N ;q
pip = 0 (i = 1, 2,
X,"(T)
= 0 (i = I ,
= 1, 2 ,
. . . , Qi)
. . ., Pi) 2, . ., Qi)
= 1, 2, = 1,
(-45)
Nomenclature A = substanceA a = cost factor of reactor B = substanceB C = substanceC c = concentration of a material d = design vector E L E2 = activationenergy F = volumetric flow rate f', f " = integrand of objective function g = nght-hand side of subsystem equation H, = scalar function defined by eq 7 H, = integrand of scalar function defined by eq 10 J = objective function defined by eq 6 K , = set of subsystems incident to y1 K1, K Z = reaction rate constants Klo, K20 = preexponential factors M = scalar function defined by eq 14 N = number of subsystems P, = number of output streams at subsystem i Q, = number of input streams at subsystem i R = gasconstant Rn = n-dimensional Euclid space S, = scalar function defined by eq 8 3, = integrand of scalar function defined by eq 11 T = final time of process operation u = controlvector V = volume of reactor v = pricevector x = input stream state vector y = output stream state vector
Superscripts T = transpose of matrix on vector ( k ) = the kth iteration 0 = initialvalue p , q = streamnumber Subscripts A = substanceA B = substanceB C = substanceC i, j , k, . . . = number of subsystems N = total number of subsystems or output subsystems 1 = inputsubsystem Symbol A = stands for definition
Literature Cited Bryson, A. E., Denham. W. F., Trans. ASME, Ser. E, 29, 247 (1962). Douglas, J. M., Rippin, D. W . T., Chem. Eng. Sci., 21, 305 (1966). Ichikawa, A., Nishida. N., Umeda, T.. paper presented to the 34th Annual Meeting, The Society of Chemical Engineers, Japan, 1969. Ichikawa, A., Fan, L. T., Chem. Eng. Sci.. 28, 357 (1973). Kobayashi. S.. Umeda, T., ichikawa, A,, Chem. Eng. Sci., 26, 1367 (1971). Kelley, H . J., Kopp, R . E., Moyer, H. G., AlAA Astro-dynamics Conference, Yale University. New Haven, Conn.. 1963. King, C. J:. Gantz, D. W., Barnes, F. J . . Ind. Eng. Chem.. Process Des. Dev.. 11, 271 (1972). Lasdon. L. S., Mitter, S. K.. Waren. A. D., IEEEAC-72. No. 2 (1967). Lee, K. F., Masso, A. H . , Rudd, D. F.. lnd. Eng. Chem., fundam.. 9, 48 (1970). Masso, A. H.. Rudd, D. F . , A / C h E J . . 15, 10 (1969). McGalliard. R. L., Westerberg, A. W.. 71th National Meeting of the AIChE, Dallas, Texas, 1972. Nishida. N., Kobayashi. S., Ichikawa, A.. Chem. Eng. Sci.. 26, 1841 (1971). Nishida, N.. Ph.D. Thesis, Tokyo institute of Technology, Tokyo, Japan, 1972. Nishida. N., Ichikawa, A., Tazaki. E., ind. Eng. Chem.. Process Des. Dev.. 13, 209 (1974). Powers. G. J.. Rudd, D. F., 71th National Meeting of the AIChE, Dallas, Texas, 1972. Siirola, J . J . , Powers, G. J . , Rudd, D . F.. AlChEJ., 17, 677 (1971). Umeda, T.. Hirai, A., lchikawa. A.. Chem. f n g . Sci., 27, 795 (1972). Westerberg, A. W.. "Decomposition of Large-scale Problems,' D. M. Himmelblau, Ed., North-Holland, 1973.
Greek Letters a = structure parameter e = step size factor defined by eq 16
242
Ind. Eng. Chem., Process Des. Dev., Vol. 14, No. 3, 1975
Received for review J u n e 12, 1974 Accepted J a n u a r y 27, 1975