Ozawa, Y., Bischoff, K. B., ibid., 7, 72 (1968). Tan, C. H., Fuller, 0. M., Can. J . Chem. Eng., 48, 174 (1970). Weekman, V. W., Jr., I d . Eng. Chem., Process Des. Develop., 7 (I), 90 (1968). Weekman, V. W., Jr., Nace, D. M., AIChE J., 16 (3), 397
SUBSCRIPTS u = unknown catalyst s = standard catalyst Literature Cited
(1970). \-".
Blanding, F. H., I d . Eng. Chem., 45 (6), 1186 (1953). Froment, G. F., Bischoff, K . B., Chem. Eng. Sci., 16, 189 (1961). Gustafson, W. R., Amer. Chem. SOL Preprints, 14 (31, €356 (1969). Henderson, D. S., Ciapetta, F. G., Oil Gas J., 65 (42), 88 (1967). Nace, D. AI., I d . Eng. Chem., Process Des. Develop., 8 (l), 24 (1969). Nace, D. M., Voltz, S. E., Weekman, V. W., Jr., ibid., 10 (4), 530 (1971).
Whitakir, A. C., Kinzer, A. D., Ind. Eng. Chem., 47 (lo), 2153 (1955 '. Wojciechowski, €3. w., Can. J . Chem. Eng., 46, 48 (1968).
RECEIVED for review December 6, 1971 ACCEPTED July 28, 1972
Optimization of a Simulation Model of a Chemical Plant Paul Friedman1and Kenneth 1. Pirider2 Department of Chemical Engineering, University of British Columbia, Vancouver 8 , B.C., Canada
A model of a chemical process plant for the polymerization of olefins in the production of polymer gasoline was simulated on a digital computer using the Chemical Engineering Simulation System (CHESS). This model was then optimized by three automatic optimization algorithms which were modified for this purpose: Deflected Gradient-Created Response Surface, Pattern Search, Complex Method. The results of the plant model simulation were within the range of reported plant data. The constrained Pattern Search Method did not handle the constraints very well, while the Deflected Gradient-Created Response Surface Technique was extremely sensitive to the values chosen for the search parameters and thus quite difficult to use. The Modified Complex Method seemed to handle the constraints quite well, was simple to use, and did not require an excessive amount of computer time.
S e v e r a l studies have been reported in the literature where the design of chemical process equipment and chemical plants have been optimized (Hwa, 1965; Hens, 1966; Singer, 1962; Bracken and McCorniick, 1968; Beauduin, 1967; Umeda, 1969; Umeda and Ichikawa, 1971; Komatsu, 1968; Barneson e t al., 1970; Hughes e t al., 1963). Concurrently with t h e above optimization studies, systems have been developed for the steady-state simulation of chemical process plants (Evans e t al., 1968). The use of some of these chemical process plant simulation systems has also been reported (Hughes e t al., 1963; James et al., 1966; Shannon et al., 1966; Johnson et al., 1968; Crowe et al., 1971; Menzies and Johnson, 1971; Johnson e t al., 1971). With t h e exception of a preliminary study by Johnson et al. (1971), these two parallel developments have not been combined. The optimization studies reported generally do not utilize full-scale detailed process models. I n most cases the processes or equipment units whose design is to be optimized ilre represented by either a series of nonlinear algebraic equations or by relatively simple FORTRAN models and are quite specific in nature. In the case of detailed chemical process plant simulation, the optimization studies reported have been
of a manual nature; case studies have been run for several different values of independent variables. Automatic optimization methods have not been utilized because of the large amount of computer calculation time required for the steadystate simulation of a large chemical process plant with recycle. Purpose and Scope
The objective of the present work mas to study the utilization of automatic methods of optimization on a chemical plant which was modeled using a chemical process plant simulation system. The process plant chosen for the study was the gasoline polymerization plant a t the Shellburn Refinery, Shell Canada Ltd., in Burnaby, B.C., Canada. The Chemical Engineering Simulation System (CHESS), developed by N o t a r d e t al. (1968) a t the University of Houston, was utilized for the process plant simulation. Three optimization methods were tested: Deflected Gradient-Created Response Surface Rlethod, Constrained Pattern Search, and Modified Complex Method. Extensive modifications were made t o the optimization algorithms so t h a t they could be successfully used for t h e chemical plant model optimization. Chemical Engineering Simulation System (CHESS)
Present address, Facultad de Tecnologia, La Universidad de Habana, Habana, Cuba. 2 To whom correspondence should be addressed. 1
512 Ind.
Eng. Chem. Process Des. Develop., Vol. 1 1 , No.
4, 1972
CHESS is composed of several parts, the first being t h e executive system of programs which transmit information
POLY MERIZA11011 REACTORS
FRESH FEED
OEPROPANIZER
a
DEBUTANIILR
a
FROM POL"
FEED THEATE
Table 1. Range of Feed Compositions Mol
$&ER
POL" FEED
EXCHANGE EFFUIEN
Ethane Propane Propene Isobutane %-Butane 1-butene Isobutene Trans-2-butene Cis-2-butene Isopentane %-Pentane 1-Pentene 2-Methyl-1-butene Trans-2-pentene 2-Methyl-2-butene Total olefins
+
Figure 1. Flow diagram of gasoline polymerization plant
through t h e process streams and store t h e calculated results. The executive will read the input data, call the unit computations in accordance with some calculation sequence, utilize convergence routines for material and heat balancing, and print out the results. T h e second part of CHESS consists of the unit computations or process modules which are mathematical models of the plant equipment. A third segment of CHESS is t h e physical properties evaluation system. The above-mentioned parts of the system are of general applicability. For a given plant, a set of data must be supplied which describes the initial stream variables and t h e equipment parameters to be used in the unit cornputations. CHESS uses the concept of modularity which allows changes to be made in the ordering of units independently of what t h e units are. This is possible since information is passed between the units in 21 standard format. I n the physical properties prediction package, basic physical constants are supplied for 62 common substances, most of which are hydrocarbons. Physical properties can also be computed for components whose basic properties are not stored, by supplying the system with these constants. CHESS utilizes t h e technique of successive substitution for the solution of t h e material and heat balances when recycle exists. The authors of CHESS utilize a n algorithm developed b y Kegstein (1958) for convergence forcing. It should be noted t h a t the CHESS system is quite large. It contains abaut 6000 source language cards and occupies approximately 40,000 words of I B M 360 core storage. However, i t was possible to run the complete program as a single load module on the I B M 360/67 system at t h e University of British Columbia. Process Description
The catalytic polymerization plant at the Shellburn refinery, Shell Oil Co. of Canada Ltd., is typical of the units based on the Universal Oil Products Co. process. A flow diagram of the unit is shown in Figure l . The poly feed (see Table I for a typical range of feed compositions) enters the unit a t about 300 psia and 90°F. It is pumped a t a pressure of about 550 psia through t h e poly feedleffluent heat exchanger and the poly feed heater and enters t h e polymerization reactors a t a temperature and pressure of 360°F and 515 psia. The temperature of t h e reacting mixture which is undergoing dimerization of propene, butenes, and pentenes is controlled b y the introduction of liquid propane from the depropanizer, which serves as a quench between the catalyst beds. If t h e concentration of olefins in t h e feed is too high, the propane stream can also be used as a recycle t o dilute the feed. The reactor effluent then passes through the feed/effluent
?&
0 1-0 9 21-24 14-20 23-26 10 5-20 6-7 2 4 2-5 8 2 8-4 3 3 3-7 2 0 1-0 2 0 1-0 2 0 2-0 3 0 1-0 2 0 1-0 2 28-36
exchanger and a pressure reduction valve, and enters the depropanizer where propane is removed as t,he top product. The propane is condensed and a part of it is recycled to the reactor while the rest is dried and sent to storage. The depropanizer bottom product is sent to the debutanizer where butane is removed as the top product. The butane is condensed, cooled, and sent to storage. The bottom product of the debutanizer is polymer gasoline which is cooled and sent to storage. The most important process variables are those which apply to t h e polymerization reactors. They are: catalyst bed temperature, feed composition, water content of feed, space velocity of feed, and catalyst age. The water content of t h e feed is maintained within prescribed limits so t h a t the phosphoric acid on the catalyst will not volatize causing a loss of catalyst activity nor will the cat'alyst be overhydrated resulting in low reaction rates. The catalyst loses its activity with time because of the formation of high polymer products on its surface. For any given feed composition a t a specified feed flow rate and catalyst age, the process variable of major interest is the catalyst bed temperature. This will, within certain limitations, determine the olefin conversion and the production of polymer gasoline. The product splits and purities are directly influenced by t h e reboiler loads and reflux ratios of the debutanizer and depropanizer. Strategy for the Development of the Simulation Model
Perhaps the most important question t h a t must be answered before t h e simulation of a chemical plant is to be attempted is: Why should this plant be simulated? The answer to this question should determine the steps to be taken in the development and use of the simulated model. The answer for t h e present study is: The plant should be simulated so t h a t t h e simulation model can then be used to optimize the plant operation, which means finding the set of important independent process variables which will maximize (or minimize) some predetermined objective function (such as conversion, cost of production of a product or products, or profit). With the above-mentioned goal in mind several criteria for t h e simulation can be set: (1) The important process variables must be located. These are t h e independent process variables t h a t might appreciably affect t h e value of the objective function. I n our case where the objective function is t h e olefin conversion in the reactor, the amount of propane quench would be a n important variaInd. Eng. Chem. Process Des. Develop., Vol. 1 1 , No. 4, 1972
513
IW REICTOR INLET TEMPERATURE * 815 .R
Table II. Fit of Reactor Model to Plant Data Typical plant data a t start-up
095
Figure 2. Maximum overall conversion of olefins vs. time .%
1
080-
-i 0 7 5 0415
ON
STREPM
ble while the amount of cooling water used in the butane cooler is only a minor significance. (2) A second factor to be considered is the equipment limitations. Because of these limitations, the methodology used in the simulation of an operating plant is different from t h a t used in the design of a new plant. The equipment design characteristics (such as heat exchanger area and number of plates in a distillation column) are fixed in a n existing plant and thus the possible operating ranges are much more limited than in the design case. Since the equipment limitations define the bounds of some of t h e process variables, this should be considered in the development of the process equipment models. For example, if it is known t h a t a certain heat exchanger cannot heat a given stream to more than some maximum temperature because of equipment and process limitations, then perhaps a simple model fixing this temperature could be used rather than a complex overall heat exchanger model. (3) A very important factor t h a t must be considered in the development of a plant model for optimization is computer calculation time. Since the simulated model may be called on hundreds (and perhaps thousands) of times by the optimization program, the computer calculation time for each plant simulation must be kept as small as possible. On the other hand the simulated model must represent the operating plant in a n adequate manner if the results of the optimization are t o be of any value. Polymerization Reactor Model
The polymerization reactors are the most important process units in the polymerization plant and thus the unit computation developed for t h e reactor must be more comprehensive than those used for the other process units. The basic plug flow reactor equation is
The reaction rate can be represented b y the expression developed by Friedman and Pinder (1971) and multiplied by the effectiveness factor, E, and the flow efficiency factor, F,, to obtain, dzoL 2.87 X 105 e--6840/T A,F,E (1 - z ) ~ __ dz Q (1 XI'
+
Equation 2 represents the change in total olefin conversion with catalyst bed depth. =i second relationship can be derived from a n enthalpy balance
The pressure drop through the bed was modeled by a n equation derived by Ergun (1952). 514 Ind.
Eng. Chem. Process Des. Develop., Vol. 1 1 , No. 4, 1972
Temp,
OF
In
out
1 2 3 4
20 350-360 390-400 38 360-370 400-4 10 60 400-420 430-450 64 430-440 440-450 Overall conversion = 0.85-0.90 Model data at start-up flow efficiency factor F , = 0.60 Bed No.
In
Temp,
out
Cumulative conversion
1
330 360 410 455
395 445 475 480
0.20 0.35 0.55 0.88
2 3 4
dz =
OF
-[e
(PVO2)(1
+ 1.75][
-4
144 DPgce3
]
(4)
where
A fourth-order Runge-Kutta method was used to solve this system of ordinary differential equations. Plant data were utilized to develop relationships for the change in the effectiveness factor and the bed porosity with time. Since the first bed acts as a filter and is deactivated much more rapidly than t h e remaining beds, a different relationship for the change in effectivity with time was used: Bed 1:
E
=
1.0 - 0.10 D ,
(5)
Beds 2, 3, 4:
E
= 1.0
- 0.0015 D,
(6)
where D , = time on-stream, days. A value was found for t h e flow efficiency factor by fitting the plug flow calculation to plant data. A comparison of the range of plant data and t h e model results are shown in Table 11. The comparison is better than that indicated in t h e table since the inlet and outlet temperatures in the plant are obtained from thermocouples which are 4-6 in. in from each end of the catalyst beds. il parametric study was made of the polymerization reactor by manually searching for the best combination of propane quench to t h e four catalyst beds so t h a t a maximum conversion would be obtained a t different times during the IO
(2)
Bed Depth, in.
No.
REACTOR
INLET TEMPERbTURE
a
815.R
i TOTAL RECYCLE RECYCLE BED 3
Figure
3'
Recycle
Figure 4. Information flow diagram for gasoline polymerization plant
lifetime of the catalyst. The results are shown in Figures 2 and 3. Plant Model
The information flow diagram for the polymerization plant model is shown in Figure 4. Some differences exist between i t and t h e plant flow diagram (Figure 1). The two polymerization reactors are combined into one unit to lower the computer calculation time. Since the two reactors normally run under the same conditions, this is quite reasonable. The reactor has been divided into four mixing units and four reactor units (one for each catalyst bed). The poly feed/ effluent heat exchanger has been divided into two simple heat exchangers, and the propane driers have not been considered. The unit computations which were utilized were all supplied by t h e CHESS system with the exception of t h e ADD3 reactor module which was previously described. A brief description of these unit computations is shown in Table 111. System Convergence Error
CHESS utilizes the method of successive substitution for the solution of the material and heat balances when recycle exists. This means t h a t the stream variables (component flows, temperature, pressure, enthalpy) of the recycle streams Table 111. CHESS Unit Computations Used in Plant Model Program
Function
DVDR
Divides one input stream into a specified number of output streams A simple linear split b y components is made of the input stream to the top and bottom output streams Mixes a specified number of input streams to produce one output stream Calculates temperature and vapor fraction of output stream a t enthalpy of the inlet stream and the specified output stream pressure Simple mode: Vapor fraction and enthalpy of output stream calculated a t the specified output stream temperature Work necessary to pump a liquid or compress a fluid to a desired pressure is calculated
DIST
MIXR VXLV
HXER
PUMP
are initially set to zero. The material and energy balance calculations converge when t h e stream variables of all the streams show a fractional change between successive loop calculations of less t h a n a specified convergence tolerance. When t h e convergence tolerance is large, the number of recycle loops required for convergence will be smaller (as will the computer calculation time) but the final values for the stream variables will not yet have reached their correct value. (As the convergence tolerance is decreased, the stream variables asymptotically approach their limiting values.) Thus there is a trade-off between computer calculation time and the accuracy of the results of the simulation syst,em. This problem is not too important when the simulation model is used for case studies where low convergence tolerances can be used because t h e differences in computer calculation times are relatively unimportant. However, it is extremely important when the simulation model is t o be used in optimization studies. Another difficulty arises in t h e use of the simulated model as a black box calculator within a n optimization program. This is, t h a t the values calculated for the objective function, for the same set of independent variables, depend on the starting values of the process stream variables. I n a case study, the recycle stream variables are init'ially set to zero, but when used in a n optimization scheme, the initial values for the stream variables for a given calculation are t h e final values from the previous calculation. This effect is shown in Figure 5 for a series of runs in which the values of the independent variables were changed (Table IV) .
m]
1
0 001,o 0001
0.001, 0 0001 00001 0001 0 01
4ool 4 1 8
0
0
O V
B
395
5 1
Figure 5. Objective function values for repeated calculations
8
a
3901
I
38 5
,
Ind. Eng. Chem. Process Des. Develop., Vol. 1 1 , No. 4, 1972
5 15
Table IV. Values of Independent Variables for Figure 5 Independent variable
Run no.
1
2
3
4
5
6
7
8
0.60
0.62
0.60
0.70
0.60
0.61
0.60
0.60
0.20 0.30 0.30
0.05 0.31 0.41
0.20 0.30 0.30
0.10 0.30 0.40
0.20 0.30 0.30
0.21 0.31 0.31
0 20 0.30 0.30
0.20 0.30 0.30
DVDR 24 split stream 17 DVDR 26 split stream 19 20 21
Table V. Total Computer Calculation Time for Test Runs Case
No. R-k
Approx CPU, sec/cycie
steps
3 5 10 50
1 2 3 4
Conv tol.
3
Total cycles, runs 1-8
0.05 0.01 0.001 0.0001
3.7 ... 11
Total CPU, sec (approx)
44
81 163
71 107
1177
27
...
the number of steps used in a given bed length. Several conh s i o n s can be derived from Figure 5 . The first is that the xlculated results are very similar for the number of R-K 2 5 ind convergence tolerance 5 0.01. The second conclusion is that even though the calculated result is dependent on the initial values of the process stream variables, in the area near the optimum where only small changes are made in the values of the independent variables b e h e e n simulation runs, the objective function calculations are more exact because of the repeated calculations made a t very similar values of the independent variables. The importance of establishing the maximum acceptable tolerances is shown in Table V. The values of R-K steps and convergence tolerance shown in case 2 were utilized for this work.
R u n 1 starts n i t h all process streams (except the poly feed input) a t zero. I n the following runs (2-8) the final values of the process stream variables from the preceding run were utilized as the initial values in the new recycle calculations. The independent variables were changed for each run to simulate the way a n optimization main program would operate when most of the independent variable changes between runs were quite small. A large change was inserted between runs 3 and 4 to test the system's response. In the polymerization reactor model, -4DD3, the set of ordinary differential equations (2-4) and Equations 5 and 6 is solved by the Runge-Kutta method which operates on continuous functions in a stepwise manner. Thus the accuracy of the method and the computer calculation time are related to
MINIMUM POINT ALONG DIRECTIONSEARCH IS DETERMINED, F T N CALLED FOR POINT CALCULATIONS
MAKE STEP,
GRADIENT DEL0
b =TURN
li 1945 - T I IF C I I l C 0 0 ,
Clli
READ FP INPUT DATA
IOOCIII
0 6
SET PENAL'W FACTORS AND TOLERANCES
SUBSET
,
ENERGY 8ALANCE CALCULATIONS UTILIZES OTHER KCTOL
Figure
516 Ind.
REACHED
6. Computer
program linkage between optimization and simulation subroutines
Eng. Chem. Process Des. Oevelop., Vol. 11, No. 4 , 1972
0 FUNC. CONPAT
1
1
YO
\
FUN
I
NEW 8bSE POINT
n
..
ZIII Will WII) Xlll YO * YY
MADE "El
YY
.
PbTTERN
YO
XI11
m
XU1 t DIII XI11
.
2 xi11
- Zlil
CONSTRAINT SATISFIED
($
O
N
CONSTRAINT S I 0T I S F I E 0
0
W
Figure 7. Computer flow diagram of Constrained Pattern Search Method
Optimization Problem
Several objective functions were optimized in the original study (Friedman, 1971) : plant profit, total dimer from plant, and olefin conversion. This report will discuss only the olefin conversion optimization since it was found to be t h e controlling objective function. I n the profit calculation when only small changes are made in the process variables, the labor, maintenance, and overhead costs can be assumed C O W stant. The utility costs for the specific case of the polymerization plant do not shox much variation for such small changes and so for any given raw material (poly feed), the objective function can be considered as the total amount of products sold. This led to the study of the total dimer from the plant as the objective function. Since the four distillation column variables do riot affect the product'ioii of dimer, but rather teiid to slow u p and interfere with the search, the major optimization study was made using olefin conversion as t,he objective function. The independent variables considered in the model (Figure 4) were the recycle f l o ~splits; one independent split in DVDR24 and three independent splits in DVDR26. The explicit constraints are t h a t the sum of the splits in each DVDR must equal unity. Limits were also set 011 the temperatures of the streams leaving the ADD3 units thus accouiitirig for four implicit constraints. T o solve the nonlinear constrained optimization problem, three methods were utilized. X brief discussion of each is given below. Deflected Gradient-Created Response Surface Technique. T h e Deflected Gradient Method developed b y Fletcher a n d Powell (1963) and based on t h e earlier work of Davidoii (1959) chooses t h e search direction by updating a n approximation t o t h e inverse matrix of second derivatives a t t h e minimum in such a way as to assure t h a t these matrices are positive definite. I t is used for unconstrained nonlinear optimization problems. T h e
above method can be uqed for constrained optimization problems b y utilizing a technique developed by Carroll (1961) who proposed t h a t a created response surface be defined a s :
The minimum of this surface is then found and used as a starting point for the minimization of another response surface corresponding to a reduced value of K . The method is repeated for a sequence of respoiise surfaces corresponding to successively smaller values of K and finally for a terminal value of K = 0. Thus a single coiistraiiied minimization has been replaced by a sequence of unconstrained minimizations. Figure 6 shows the program linkage between t,he main calling program, opt'imizatiori subroutine, and CHESS system subroutines. Constrained Pattern Search. Hooke aiid Jeeves (1961) developed an accelerated climbing technique with ridge following propert'ies, sometimes referred to as t h e P a t t e r n Search Method. T h e method is based o n a local exploration in t h e directions parallel t o t h e coordinat,e axes at points a given distance from t h e base point. A larger move in a direction determined by these tests is then attempted. If this move is successful, aiid if a new local exploration indicates t h a t t h e direction is a good one, a further move is made with a much larger step size. So even t'hough the technique starts cautiously with short steps, the step size increases with repeated successes since it is assumed t h a t the previous moves are worth repeating. If t'lie objective function does not improve along a given direction, smaller steps are used again and if this is not successful, a new stage and a new pattern is started. The above-described method has been combined Tvith added logic to handle the constraints (Figure 7 ) . Each new point is Ind. Eng. Chem. Process Des. Develop., Vol. 11, No.
4, 1972 517
_ _ _ _ _ ~
~
~~~~
Table VI. Comparison of Optimization Methods Tested
Method
X(1)
X(2)
DGCR COXPAT MOCOMP
0.700 0.700 0.700 0.650 0.670 0,680 0.675 0.675 0.618 0.618 0.622 0.617 0.613 0,620
0.100 0.200 0.150 0.100 0.120 0.175 0.175 0.046 0.046 0.035 0.051 0.055 0.042
DGCR COXPAT DGCR COSPAT MOCOMP
X(3)
X(4)
Initial function volue
0.300 0.300
0.400 0.400 0 100 0.250 0.200 0.200 0.187 0.187 0.400 0.400 0.407 0.403 0,404 0.412
33.098 33.098 35,032 38.245 36,646 35.298 36,470 36.470 40.367 40.367 40.376 40.310 40.317 40.340
Initial point
-~
0.100
0.100 0.350 0.300 0.150 0.225 0.225 0.300 0.300 0.298 0.299 0.298 0.304
I
initially checked for violation of the explicit constraints. If this occurs, the independent variable which violated the constraint is placed just within the constraint boundary. This is followed by a check of the abortive implicit constraints which when violated (as in the case of a negative flow split) may cause the simulation program to abort. (Convergence test fails in the case of negative flows.) If this type of constraint is violated, the point is immediately rejected. If these constraints are satisfied, the simulation program is called on, and the objective function and the implicit constraints are calculated. If the implicit constraints are violated, the point is rejected. Modified Complex Method. T h e Complex Method developed b y Box (1965) handles constraints b y t h e use of a flexible figure of more t h a n n 1 vertices (where n is t h e number of independent variables) which can expand or contract in all directions. T h e objective function is evaluated a t each vertex a n d t h e vertex W yielding t h e poorest function value is rejected a n d replaced by a point which is located a distance a, times as far from t h e centroid C of t h e remaining points as t h e distance from t h e centroid t o t h e reiected vertex in a direction defined b y a vector pointing from t h e rejected point to t h e centroid.
+
Z1,N
=
zi,c
+
ae(Zi,c - Zi,tO)
(8)
where
X(1)
X(2)
X(3)
X(4)
Final function value
0,599 0.640 0.619
0.029 0.035 0.021
0.316 0.250 0.324
0.505 0.350 0.406
40.493 39.899 40.343
1090 812 912
0,632 0.654 0.611 0.611 0.602
0.178 0.145 0 047 0 046 0 017
0 229 0 200 0 ,302 0 305 0 ,333
0.193 0.167 0.403 0.401 0.464
39.350 38.762 40.410 40.407 40.565
494 766 307 447 441
Final point
some instances where it had previously stopped. The Modified Complex Method is shown in Figure 8. Comparison of Optimization Methods Tested
Table VI lists the results obtained with each of the three optimization methods tested on the polymerization plant model. For the specific set of equipment parameters and process steams defined in the polymerization plant model, the optimum seems to be: DVDR
Stream
24 26 26 26
17 19 20 21
xi,N = '/Z(X%,N
+
(9)
21,~)
A number of modifications were made (Friedman, 1971) to the original Comples Method because of the large amount of computer calculation time necessary for objective function calculation. ii new logic step was added to reflect the centroid through the best point thus allowing the search to continue in 5 18
Ind. Eng. Chem. Process Des. Develop., Vol. 1 1 ,
No. 4, 1972
Flow fraction
0 0 0 0
602 016 333 464
The constrained temperatures which must be less than 945'R are: Stream no.
Temp, OR
6 8 10 12
891 944 944 945
Xiid the objective function value, the number of moles of dimerized product:
f(z)
If this new point N is a feasible point, the function is evaluated there and the above-described expansion step is repeated. If, however, the function value of the new point is worse than t h a t of the previously determined worst point or if a n implicit constraint is violated, the vertex is moved halfway in toward the centroid until a valid point is obtained.
CPU time, sec
=
40.565
I n the first case, in which the starting point (or complex) is quite distant from the optimum, the three methods showed improvement in t h e value of the objective function but the constrained Pattern Search Xethod (COSP-IT) is definitely inferior to the other two. The Deflected Gradient-Created Response Surface Technique (DGCR) achieves good results but this is quite fortuitous as only a prior knowledge of the response surface a l l o w one to choose values for the search tolerances and penalty factors so that the search can proceed. Since this method utilizes the gradient, the system convergence error, previously discussed, has a large effect on the search direction. The Alodified Complex Method (lIOCO1\IP) gave good results and was still progressing toward the opt,imum rvlien the search was terminated. The computer central processing unit (CPV) calculation time required is of the same order of magnitude for each of the three methods.
u
Nomenclature
Q
MOCOMP
MOVE BACK
TO
1/2(5,
",I
+
A, A,
C A L C U L A T E CENTROID O F REMAINING POINTS
I
x
I
f
ii$L{YE: OONSTRAINTS ATISFIED
11 1 -
-1
CONSTRAINTS
-1
REFLECT CENTROID G FEH :;; ,BEST I,
.
I
EVALUATE
+
0,
(i-
c
,
I
I
REPEAT
CALCULATE CENTROID O F A L L POINTS
I YES
expansion factor in Complex Method cross-sectional area of catalyst bed, ft2 exteriial surface area of polymerization reactor, ft2 C(s) = constraint function D, = catalyst particle diameter, ft D, = number of days on-stream E = effectiveness factor f(z) = objective function I; = molar flow rate, mol/hr F , = flow efficiency factor gc = conrerqion factor, mass to force units h , = radiation heat transfer coefficient, Btu 'hr-ft*-OF H = enthalpy, B t u K , = number of vertices in Complex Xethod K = penalty factor Ji = number of coiistraints n = iiumber of indepeiiden1 variables P = pressure, psia P ( z , K ) = created response surface function Q = volunietric flon rate, ft3, hr I' = rate of reaction, mol olefin/ hr cc catalyst Re, = particle Reynolds number T = temperature, O R T , = ambient temperature, OR V , = linear gas velocity, ft/lir T T = weighting factor z = independent variable $01. = fraction of olefins converted into product ~ O L = mole fraction of olefinr z = length of catalyst bed, ft AH,,, = average heat of reactioii, Btu/'mol a,
WAY CENTROID
.
POINT I N COMPLEX
1/2
I
Figure 8. Computer flow diagram of Modified Complex Method
(CPU time was lowest for COXPAT because the search failed earlier than with the others.) I n the second case, the starting point is still fairly far from the optimum. The D G C R technique is shown to be much superior to COXPAT. However, the D G C R method still stopped quite a distance away from the optimum. The initial starting point (or complex) in the third case is very close to the optimum and lies on one constraint. Both the D G C R technique and CONPAT did not progress very far before they stopped and only MOCOMP seemed to work effectively near the constraints. I n fact, the optimum point calculated by MOCO?*lP lies on the intersection of the three active constraints. Conclusions
It has been shown t h a t it is practicable to utilize a n automatic optimization method for the optimization of a chemical plant model with recycle streams which has been simulated by means of a chemical plant simulation system. T o be able to utilize a n automatic optimization method without a n excessive requirement of computer calculation time a strategy of model building must be used t h a t simplifies the calculations wherever possible but a t t h e same time supplies reasonable estimates for t h e major decision variables and constraints. The Constrained Pattern Search method did not handle the constraints very well and did not give good results in t h e present study. The Deflected Gradient-Created Response Surface technique achieved better results but was extremely sensitive to the values chosen for the search parameters and penalty factors and is thus quite difficult to use. The Modified Complex X e t h o d worked quite well, even in the region close t o the constraints. It is very simple to use and does not require excessive amounts of computer time.
= = =
GREEKLETTERS E
=
p
=
p
=
void fraction density, lb 'ft? viscosity
literature Cited
Barneson, R . A , , Brannock, X. F., Moore, J. G., lIorri. G., Sprague, C. I{., ibirl.) 64 (4) 39-46 i\-"--,. lOAR\
Fletcher, It.,Powell, 11. J. D., Computer J . , 6 , 163-8 (1963). Friedman, P., PhD thesis, Cniversity of Britiyh Columbia, 1971. Friedman, P., Pinder, K. L., Incl. Eng. C h c m Process Des. Decelow.. 10 14). 548-51 11971). Hens, A,' AI., AIChE-ICE Symp. Series, S o . 4, 110-18, London, England, 1965. Hooke, R., Jeeveq, T. -2,J . Ass Comp. J l a c h , 8 (2), 212-29 f\ -l-O- -A ,. l \
Hughes, R. R., Singer, E., Souders, JI.> TT-orlrl Petrol. Congr. Proc. 6th (Frankfurt),Sect. TTI, Paper 17, 1963;Hwa, C. S.,AIChE-ICE Symp. Sene-, S o . 4, (2-83, London, England, 196.5. James, J. L., Gardner, N. F.,Ileiiihart, L. It.,Helleiiack, L. J., ICE Symp. Series, S o . 18, 11-18, London, England, 1966. Johnson, A . E., Aizawa, JI.,Petryschok, W. F., Brit. Chon. Eng., 13 ( l o ) , 1432-8 (1068). Johnson, A. E., Cheng, 15'. K., Karalenas, J., Leavoy, W., Leslie, E., C u d l i ~ D. , A , , Brit. Chen,. Eng. Process Tcchnol., 16 (lo), 923-9 (f971). Komatsu, S.,Ind. Eng. Chem., 60 (2) 36-43 (1968). AIenzies, 11. .4,,Johnson, A. I., Can. J . Chem. Eny., 49 (31, 407-11 119i1). Ind. Eng. Chem. Process Des. Develop., Vol. 1 1 , No. 4, 1972
519
nrotard, R. L., Lee, H. If.,Barkley, R., W., Ingels, D. At., "CHESS, Chemical Engineering Simulation System, Systems Guide," Tech. Publ. Co., Houston, Tex., 1968. Shannon, P. T., Johnson, A. I., Crowe, C. Jf., Hoffman, T. W., Hamielec, A. E., Woods, D. R., Chem. Eng. Progr., 62 (6), 49-59 (1966).
Singer, E., Chem.Eng. Progr. Symp.Series, Vol. 58 (37), 62-74 11962).
Umeda, T., Ind. Eng. Chem.Process Des. Deoelop., 8 (3), 308-17 (1969).
Umeda, T., Ichikawa, A., i b i d . , 10 ( 2 ) , 229-36 (1971). Wegstein, J. H., Comm. Ass. Comp.Mach., 6, 524-44 (1958).
RECEIVED for review December 20, 1971 ACCEPTF~D June 1 , 1972 Financial support was provided by the National Research Council of Canada and the University of British Columbia. Computer programs for optimization routines used in this paper will appear following these pages in the microfilm edition of this volume of the Journal. Single copies may be obtained from the Business Operations Office, Books and Journals Division, American Chemical Society, 1156 Sixteenth St., N.W., Washington, D.C. 20036. Refer to the following code number: PROC7 2 4 1 2 . Remit by check or money order $4.00 for photocopy or $2.00 for microfiche.
Studies on Tubular Flow Reactor with Motionless Mixing Elements Varadarad Jagadeesh and Mandavilli Satyanarayana' Department of Chemical Engineering, Indian Znstitute of Technology, Jdadras-36, India
A novel design of an annular tubular reactor incorporating motionless mixing elements of helical shape i s studied. The results of pressure drop experiments and residence time distribution are presented which show the equipment to be a better plug-flow reactor than an empty tube. The conversions obtained at low flow rates for the saponification reaction of ethyl acetate and sodium hydroxide are given, and the implications of experimental results discussed.
In-line mixing devices are of great interest to chemical engineers inasmuch as they have certain inherent advantages in mixing, heat transfer, and reactor applications. Pattisori (1969) describes an interesting piece of equipment, the 110tionless Inline Mixer, containing right- and left-handed helical mising elements, fixed alternately in a tube. 3Iising takes place owing to a unique flow pattern resulting from the geometry has been suggested by Jagadeesh and Rajendran (1970) arid t,he present paper deals with a preliminary investigation on this modified version from a chemical reactor viewpoint. The equipment is described as well as experiments t h a t were conducted on (a) pressure drop characteristics, (b) residence time distribution, and (e) conversion studies for the saponification reaction between ethyl acetat'e and sodium hydroxide. Equipment
The equipment consists of an outer pyrex tube 1210 mm in length and 35 mm i d . , provided with two inlets a t one end, opposite to each other and inclined a t 45" with the axis of the tube. The inlets are 15 mm in diam. Four intermediate outlets are provided 212 mm apart along the length of the t'ube. These exits serve as manometer tappings for pressure drop experiments and could be used to run out react,or effluent a t intermediate lengths of the tube when conversion experiments are done. dlternatively, they could be used to insert thermometers, n-ith slight modifications when heat transfer experiments are conducted. Those exits not in use are normally
To whom correspondence should be addressed. 520 Ind.
Eng. Chem. Process Des. Develop., Vol. 1 1 ,
No. 4, 1972
kept closed by spring-type pinch cocks on flexible rubber tubing fitted onto the side exits. Concentric with the outer pyreu tube is an inner brass tube 22.4 mm o.d., 18.5 mm i.d. il brass flange brazed to the tubing is bolted to a corresponding Perspex flange fastened to the pyres tube with an adhesive (draldite). Figure 1 s h o w in detail the liquid-tight gland a t the opposite end, as \\-ell as the details of the reactor. T o the brass tube are brazed the helical mising elements. Figure 2 shows clearly the geometry of the elements. Each element consists of two 1/a2-in. brass sheets, 6.2 mm wide, carefully cut and twisted to form part of a helix 60 mm pitch. The elements are half-a-pitch in lengthLe., 30 mm. These are brazed to the outer circumference of the tube 180" apart; alternate elements are right- and lefthanded and each one is a t right angles to the adjacent ones. I n all, the reactor contains 35 elements, 17 left-handed and 18 right-handed. The action of the mixing elements is somewhat similar to that of the Static Inline Mixer. The two streams entering the reactor get split into two a t the first element. Each qplit stream encountering the helix is twisted clockwise (if the element is right-handed) : a t the next element it is twisted counterclockwise (since the element is left-handed). The stream splitting and recombination contribute to effective mixing between fluid elements as they travel along the reactor. The split layers are also constantly moved back and forth with respect to each other since the helical path length traversed by the fluid element depends on the radial distance from the center. This is referred to as layer-layer displacement. Apart from this is the mixing action due to the swirling motion induced by the helices.
