The ascend-II system—a flowsheeting application of a successive quadratic programming methodology

Computers and Chemical Engineering Vol. 7, No. 5, pp. 615-630, 1983

0098-1354/83 $3.00+ .00 Pergamon Press Ltd.

Printed in Great Britain.

THE ASCEND-II SYSTEM--A FLOWSHEETING APPLICATION OF A SUCCESSIVE QUADRATIC PROGRAMMING METHODOLOGY MICHAEL H. LOCKEt and ARTHUR W. WESTERBERG* Department of Chemical Engineering, Carnegie Mellon University, Schenley Park, Pittsburgh, PA 15213, U.S.A. (Received 1 February 1983) Abstract--This paper describes ASCEND-II, an interactive equation based flowsheeting system capable of doing performance, design and optimization calculations with comparable ease. The optimization capability is a new, as yet unpublished, successive quadratic programming method operating in the reduced space which is the size of the degrees of freedom for the problem. Its performance is very encouraging. We spend considerable time describing the ASCEND-II system details because we believe the concepts underlying it make the system a powerful aid for design. It allows one to evolve toward the solution to a problem by allowing one to add and delete model complexity, interactively alter which variables are to be calculated and which fixed, single step through the convergence of a model, move easily from design to optimization and back, use any variable as an objective to be maximized or minimized, and to place an upper or lower bound on any variable desired. Scope Mathematical techniques have recently become available which make the computerized optimization of process flowsheets practical. These techniques find the optimal values of the flowsheet parameters while moving towards the solution to the set of algebraic equations defining a flowsheet. They are far more efficient than early techniques which employed a flowsheeting program as an "inner loop" with an optimizer choosing values of the parameters in an "outer loop"[1]. In these early methods, the optimizer would choose values for the unspecified parameters and feed them to the flowsheeting program which would then perform a complete flowsheeting simulation based on those parameters and return the value of the objective function to the optimizer which would readjust the parameters. Usually pattern search methods were used to find the optimal values of the unspecified parameters. It was not unusual for the flowsheeting program to go through as many as 1000 iterations. At a cost of several dollars per iteration, optimization was expensive. To make matters worse, the user might neglect to put a specification on the quantity of product and find that the lowest cost flow~heet was no flowsheet at all. Two points must be stressed here. The first is that optimization with an "inner loop" flowsheeting system and an "outer loop" optimizer is expensive. Newer methods which converge flowsheet equations as they move towards the optimum will save time and money. The second point is that engineers are human and they will make errors. Batch computing gives no feedback during the course of computation and allows the engineer to make large, expensive errors. Giving the user more control over the details of computation through interactive computing can help to minimize errors. An interactive flowsheeting program such as ASCEND-II, which will be described in this paper, may help to demonstrate the advantages of interactive equation-based computing. As mentioned previously, methods have been developed which allow simultaneous optimization and solution of process flowsheets. These methods remove the need for a flowsheeting program within an optimizer. All variables are updated each iteration. These new optimization methods seem to fall into two categories: those which interact with traditional sequential modular flowsheeting systems and those which are based on equation solving methods. Probably the most promising of the sequential modular based methods is the method recently developed by Biegler[2, 3]. His method, based on the Han-Powell algorithm [4], has been imbedded in the Flowtran flowsheeting system. An optimization code replaces the convergence acceleration code and chooses new values for the decision variables at the same time that new values of the torn variables are guessed. The Hessian matrix which he constructs is in the space of the decision and torn variables only, rather than in the space of all the variables as in Powell's original algorithm. The remainder of the flowsheet variables are eliminated by the flowsheet equations. Gradients are calculated numerically. *Author to whom correspondence should be addressed. $Present Address: Department of Chemical Engineering, University of Massachusetts, Amherst, MA 01003, U.S.A. 615

616

M. H. LOCKEand A. W. Wm1~t~RG Equation-based flowsheeting methods are rapidly gaining popularity, at least among educators. At the 1981 Annual Meeting of the AIChE in New Orleans, all papers in the "Computers in Process Design and Analysis" Symposium concerning the development of new flowsheeting systems described equation-based systems. For a recent review on the topic, see Shacham et al. [5]. Examples of such systems are QUASILIN (Gorcyznski & Hutchison[6]), SPEED-UP (Sargent & Perkins[7]), SEQUEL (Stadtberr & Hilton[8]) and FLOWSIM (Babcock[9] and Shacham et al. [5]). Equation based systems have many advantages over traditional flowsheeting systems. One advantage, especially for those based on Newton-Raphson or Quasi-Newton methods, is that they can be easily extended to perform dynamic simulation and optimization calculations (Sargent & Perkins [7] and Berna et al. [10]). A modification of the Han-Powell method (Locke, Edahl & Westerberg [l l ]) has been added to the ASCEND-II equation-based flowsheeting system. This modification sets up a quadratic programming problem (QPP) in only the space of the decision variables, rather than in the space of the entire set of variables. The theory behind the new algorithm will be briefly explained later. Other workers have also noticed the advantage of equation-based optimization procedures. In particular, Murtagh[12] has described a method which uses the MINOS-AUGMENTED (Murtagh & Saunders [13]) optimization package to find the optimal solution of flowsheeting problems. The method seems quite promising, especially if the problem is mostly linear. Conclusions and Significance--We have described many of the details behind a new flowsheeting system ASCEND-II. Computational evidence suggests it is a powerful design aid with a very efficient optimization capability. Typically convergence with this new quadratic programming algorithm occurs in a number of iteration comparable to three "design" calculations for problems having up to five degrees of freedom. Experience also shows that an experienced user can quickly reach solutions to difficult problems with the facilities provided, as well as use the tool to understand the correctness of the solutions found. 1. A NEW SQP O P T I M I Z A T I O N A L G O R I T H M

We have recently developed the theory (Locke et al.[11]) and implemented the resulting algorithm for a new successive quadratic programming (SQP) optimization algorithm. It is similar to Powell's implementation of an SQP algorithm (Powell[4]) except (and it is a very big "except") the storage space needed is drastically reduced. In an SQP algorithm one must build up a quadratic approximation of the objective function in terms of the problem variables. This entails either analytically evaluating or numerically estimating the "Hessian" matrix of second partial derivatives of the objective function. The drawback in implementing Powell's original algorithm is that one must have sufficient storage space for an estimated Hessian matrix the size of the number of variables in the problem. A problem with 10,000 variables would require one hundred million storage locations for the Hessian, far too many for even the most advanced machines. Berna et al. in 1980110] published an algorithm which is based on Powell's method and which decomposes the optimization problem so that a reduced quadratic program is set up in only the decision variables. The algorithm is rather complicated to implement. A new algorithm which we have just developed is very much like an incorrect algorithm published by Berna et al. in 1978114] and is an outgrowth of work done by Edahl[15]. A detailed description of the algorithm can be found in Locke et al.[11]. The major difference between this algorithm and Berna's (1980) algorithm is in the estimation of the reduced Hessian. Using Berna's algorithm the Hessian is first "estimated" in all the variables. This Hessian, which is the same as the one Powell would estimate, is then projected into the space of the

decision variables. The problem of overloading available core is circumvented by never actually storing the Hessian, but rather by saving and using the vector updates which would be used to estimate the Hessian matrix during previous iterations. This new algorithm estimates the reduced Hessian directly. We never deal with the full Hessian matrix. Gabay[16] recently and independently published an algorithm which is very similar to ours. He presents his algorithm for equality constrained problems only, but clearly his intention was the more general problem which includes inequality constraints that our version explicitly addresses. His paper also considers at length the convergence characteristics of the algorithm he describes. The optimization problem can be stated as follows: rain O(z) s.t. g ( z ) = 0 Zmi n ~_~ Z ~ Zma x

(PI)

Z E R s+r

g : R~+'.-,R

~

(I): Rn+r---*R I.

For a typical flowsheet calculation of the type we are considering, n may be on the order of 10,000 with r on the order of 10. This formulation is quite general as general inequality constrains can be converted to equality constraints through the use of bounded slack variables. By partitioning the variables into two sets, the independent or decision variables, u, and the dependent or pivoted variables, x, we now show how a QPP can be set up in the decision variables only. Let xk and uk be the values of the dependent and independent variables respectively at the present iter-

The ASCEND-II system--A flowsheeting application

617

ation. Linearizing the equality constraints about this point gives:

no change in the values of the decision variables and the equality constraints are satisfied.

g(xk +,, uk + ~) ", g(xk, uk) + (Og /Ox r),,axk

1.1 Sparse matrix considerations Note that we must compute the jacobian matrix [dg/dx rk at each iteration of the SQP algorithm. This is the same matrix which must be computed when performing a design or simulation calculation without optimization. It is generally a sparse matrix, containing about three to five nonzeros per row. Because of it structure, we obviously use sparse matrix methods to store and manipulate it. We never actually calculate the inverse of this jacobian. Instead we calculate the sparse L/U factors using the MA28 series of routines from the Harwell library[17]. See Westerberg [18] for a discussion of many of these ideas.

+ (Og/~u~)kauk. Setting g(x~÷t, Uk+O equal to O, and requiring that (dg/gxr)k be nonsingular yields: Ax~ = -- ( ag /Oxr)~- t {g(xk, Uk) + ( ~g /~Ur)kAUt,}.

For any choice of Auk, a value for Axk can be calculated. The problem can then be stated as one to calculate Au~', the optimal change in the decision variables: min F(Auk) ~ at + bkrAuk + l/2AukrCAuk Auk s . t . Umin - - Uk <__A U k <__ Umax - - Uk

Xmi. -- xk <-- -- (Og/tgXr)kl {g(uk, Xk)

+ (Og/Ou~kAuk} <_ Xmax - xk

(I'2)

Auk E E"

xk ~ E ~ g : E"+r--,E ~

where

F(Auk) = @(xk + Axk, uk + Auk). F(Auk) is then expanded in a Taylor series about the point Aut = 0 giving rise to the QPP associated with (P2). The linear term in the QPP is the constrained derivative bk = 6016u = (Oc~ldu) - { (~g/~x 7) - ' (Og/Ou q } T(~O/Ox). The Hessian matrix C is an approximation to the Hessian of the Lagrange Function and is initialized as the identity matrix. The QPP is in the space of the decision variables. Lower and upper bounds on the pivoted or x variables are converted to linear inequality constraints on the decision variables. Once the QPP is solved for the optimal changes in the u variables, the changes in the x variables can be found from: ax~ = - (~g /ax

q~-' {g(x~, u,) + (~g /~u')kaUk}.

Usually the full step predicted, a {Ax, Au } with a = 1, is taken. However, it has been found that in the early iterations especially the step must be scaled down. To find the scaling factor, a, we use the ~F function suggested by Powell. We can optionally turn offusing this test and always take a full step, a move which more frequently than not has proved to be a good idea. The tp function test often slows down convergence to the right answer even when close to it. After the variables are updated, the constrained derivative is calculated at the new point, and the Hessian matrix C is updated. The new QPP is set up and solved and the process is repeated until there is

2. THE ASCEND-IIFLOWSHEETINGSYSTEM We shall next use a large part of this paper to describe ASCEND-II, an equation based flowsheeting/optimization system. We believe that this system offers features necessary for process optimization calculations. ASCEND-II has a basic design philosophy which is to allow a design engineer to evolve interactively to his problem definition and its solution. We shall describe this aspect of ASCENDII first. The system also contains many features to allow for convenient problem definition by the engineer and to initialize variable values. Finally it has an interesting internal data structure which we believe is a significant conceptual contribution; we shall describe this aspect, too.

2.1 Evolving to a problem solution Most large scale flowsheeting systems run in batch mode, allowing the user no freedom to control the calculations once they have begun. In this mode, the user prepares his input as a disc file or a stack of cards which the program then reads. After some period of time, from a few minutes to several hours, he receives a solution to his problem. Batch processing encourages a user to attempt to solve a large problem in one attempt since turnaround time is large. These large attempts fail more often than not. A more efficient approach is to solve the large problem in stages, beginning with one or a few pieces of equipment and working up to the complete flowsheet. When it comes to specifying the type of physical property calculation necessary for flowsheet simulation, most users automatically choose the most rigorous type of calculation available. They know that this choice will be the most accurate, even though it will also be the most expensive, and may not be necessary. The problem here is that users generally do not have a convenient way to move back and forth between simple and rigorous physical property calculations to see which is suitable. In designing ASCEND-II, we envisioned a flowsheeting package that would allow the user to evolve towards the solution of his problem along three independent axes, as shown in Fig. 1: a m o d e l complexity axis, a calculation type axis, and a computation detail axis. In the direction of model complexity, the user can begin with simple models of one or two process units.

618

M. H. LOCKEand A. W. Wma~XBERG

Control of Computation

Modeling

B

Typecalculaflfon .~ Fig. I. Modeling, calculation and computation axes.

This model will (hopefully) converge quickly and can be used as a basis for further calculations. The user can then probe this model, for example by adding more complex physical property calculations. The converged values of the simplified model are used as a starting point for this more rigorous calculation. If the more complex physical property calculation is necessary, it can be left in. If not, the user can readily return to the simple calculation. Another way of increasing model complexity is by adding more process units to a flowsheet. The user can add process units one by one, each time probing the complexity necessary in the physical properties, until the entire flowsheet is modeled. The converged values obtained from an earlier flowsheet would in each instance be used to initialize the current flowsheet model. Proceeding in this manner, by making small moves, probing, initializing with a previous solution, helps to avoid creating models which are not legitimate and to insure convergence at each stage. This capability also gives the user confidence in the numbers which the machine generates since he is able to re-examine the progress of the model at each stage. The second axis in Fig. 1 represents user control of the type of calculation. The simplest type of calculation performed by ASCEND-II is a simulation calculation. This is the type of calculation most easily performed by sequential modular simulators. It is often a "rating" calculation. All flowsheet inputs and unit parameters are specified by the user. The system calculates unit outputs and intermediate stream flows. Because the user can more easily visualize that a unit must operate if its input streams, sizes, and operating levels are given, a simulation calculation is the type of calculation for which one's intuition is strongest. A person will more likely handle the degrees of freedom correctly in this mode. A more difficult calculation is a design calculation. In this calculation, the user specifies an output variable value and asks the system to calculate a unit parameter or flowsheet input. It is more difficult in two ways: 1. The user may not specify easily a legitimate set of variables, leaving a singular set of equations. 2. It is possible to put a specification on the outlet stream which is impossible to meet, in which case the equations do not converge or they converge to a nonsensical solution, and the user is left wondering what happened. A simple example of this second case is a distillation column. In the simulation calculation, the user specifies the feed, feed tray location, reflux rate,

and number of trays. The simulator calculates distillate and bottoms stream compositions and internal stream compositions. Given the column inputs and unit parameters, the user is confident that outputs can be calculated. In a design calculation for a column, the user may specify the feed, feed tray location, number of trays, and a mole fraction in the distillate. The program would calculate reflux rate and all unspecified stream flows and compositions. The problem here is that the number of trays specified may be below the minimum number needed to perform the required splR. No value for reflux rate can be calculated under these conditions, and the equations will not converge. The correct way to solve the above problem using ASCEND-II is to start out with a simulation calculation; i.e. specify reflux rate, number of trays, etc. to get a converged solution to use as a starting point. Then, one can adjust the reflux rate and number of trays until the target distillate mole fraciton is close to its desired value. Now one can switch to a design calculation, specifying the mole fraction and releasing the specification on the reflux rate. The equations should converge quickly to a meaningful solution. The current version of ASCEND-II is able to perform simulation, design, and optimization calculations. Kuru[19] has laid the groundwork for adding dynamic simulation. The idea of progressing from a simple to a more difficult calculation extends to optimization as well. A user should not rush in and try to perform an optimization immediately. He should first become familiar with his flowsheet by performing a simulation. He should move the specified variables around to see the effect on the objective function. This way he will have some idea of what should occur when he finally runs the optimization. The third axis in Fig. 1 represents user control of the computation sequence. At the lowest level the user can be totally unaware of the details. We have provided a disc file containing the commands necessary to solve a standard problem. The user can tell ASCEND-II to read commands from this file to solve his problem. An experienced user would exercise more control over the details of the computation. Rather than telling ASCEND-II to read commands from a file, he can give the commands himself interactively. He can input his flowsheet interactively and give the commands to convert it to executable form. Once this is done, he can interactively input needed initial values and give the command to initialize the rest of the flowsheet variables based on this input, some of which may be guessed. Now he is ready to solve the equations. He may wish to perform an iteration of the solution algorithm and then rescale the variables and equations. If he is having trouble converging, he may isolate the troublesome unit model and work only on it for a few iterations. He gives these commands interactively, in effect controlling the details of the computation. An analogy might be the driving of a car where the user learns about his problem and how to maneuver around the obstacles. It can be a task done well without ones knowing details on how the engine and transmission work. We shall now spend considerable time now describing many of the essential details of the

The ASCEND-II system--A flowsheeting application ASCEND-II system. We believe that the details of such a system can make or break it as a system useful for doing design and optimization calculations. Unfortunately we can only argue intuitively that the features to be described "make" ASCEND-II a useful system. 3. STRATEGY FOR SOLVING

3.1 Fixed and calculated variables In any chemical engineering computer simulation there are more variables in the problem than there are equations. The user is required to fix a correct set of as many variables as are necessary to produce a "square" set of equations, i.e. one with as many variables to be calculated as there are equations. This task is far from simvle for a user. To aid in the mechanics of doing it, each variable has a calculate/fix flag, with defaulting used to free the user from having to specify this status for every variable. The calculate/fix flag for physical properties such as molecular weight, for instance, default to fixed. Stream flowrates default to calculated. These flags are easily altered by the user. Any variable in the flowsheet can be set either to calculated or fixed. This is one place where ASCEND-II differs from traditional simulation packages. In sequential modular simulators, the unit models are written to calculate outputs from fixed inputs. The program decides ahead of time which variables are fixed and which are calculated. ASCEND-II is an equation solver which treats all variables the same. This fact has allowed us (in fact has almost required us) to develop an evolutionary strategy for solving flowsheets. Changing the degrees of freedom for a problem within ASCEND-II can be done without modifying the flowsheet. With the approach taken, a design calculation is solved as easily as a simulation calculation. 3.2 Inputting variables The user of ASCEND-II can interactively change variable values at any time (using the RUNVAL command). This feature allows the user to change the value of one or several variables in a problem and do a rerun in just a few seconds time. 3.3 Initialization At the start of the solution procedure every variable in a problem must have an initial value. It would be quite difficult for the user to initialize every variable with a value close enough to the solution for the equations to converge. F o r this reason initialization subroutines have been written for each unit model. These routines require that inputs to a unit are initialized either by the user or by another unit initialization routine. The routines calculate initial

619

values for outlet variables based on the inputs. They work in much the same was as a sequential modular flowsheeting package. The initialization routines can be as simple or as complex as the writer of the model wishes. They need not solve the equations of a unit. They must only get an initial guess close enough so that the Newton-Raphson procedure can solve the equations. To illustrate how the unit initialization routines work, consider the splitter of Fig. 2. Given values (which may only be guesses) for the variables of stream S1 and for the split fraction, a, the split initialization routine gives initial guesses to the variables of streams $2 and $3. The flow of $2 is initialized to a times the flow of Sl. The flow of $3 is the flow of S1 minus the flow of $2. Mole fractions of $2 and $3 are initialized to the mole fractions of S1. The initialization of the flowsheet in Fig. 3 would work as follows: 1. The user gives initial guesses for inlet stream Sl and recycle stream $4. 2. Given S1 and $4, the initialization routine of the mixer calculates $2. 3. Given $2, the reactor initialization routine calculates $3 and any internal variables that can be calculated. The user may be required to guess some internal variables of the reactor. 4. Given $3, the column initialization routine calculates initial guesses for all the internal flows in the column and stream $5. It may override the user's initial guess of stream $4 if he does not set a flag forbidding this. Once again the user may be required to set some internal variables of the column. To simplify the user's input, a feature called Macro Expansion was devised. Macro Expansion allows complex structures to be created from simple input. It is described in Section 7.2.3. A user could, for example, creat a macro called REACRECY to input the entire flowsheet of Fig. 3. One unit, called REACRECY, would "expand" into the mixer, reactor, and distillation column. When the user initialized a flowsheet with a REACRECY unit in it, he would have two options. He could specify an initialization sequence which named the individual units created by expansion ,,

$3

Fig. 2. Simple splitter.

S4

S2

-| -I

REACT

S2

! I

$3

Fig. 3. An example flowshoet.

1

$5

620

M. H. LOCKE and A. W. WESTm~aERG

• The user is familiar with a flowsheet and is able to give initial guesses which are close to the solution. Initialization flags are used to tell the initialization routines not to initialize a variable. If set to C, the variable will be initialized. If set to F, it will not. The default for every variable is to be initialized by the initialization routine. The user may override the default flags using the appropriate commands. The default initialization flags are valid until initialization is completed. At that time, the flags are interpreted as calcuiate/fix flags for the computation. If the user re-initializes a flowsheet after starting a calculation to solve the equations, the calculate/fix flags act as initialization flags. We have found that using the same flags is actually very convenient for the user. Variables to be fixed for a computation should not be changed if reinitialization is being used to aid convergence.

(mixer, reactor, and column) in the order that he wanted them initialized, not mentioning the REACRECY unit at all. The user's second option is to name the REACRECY unit when he gives the initialization sequence and not name any of the units it expands into, clearly a simpler alternative. The creator of the REACRECY macro may designate the system initialization routine for macros, and it simply will replace the RECREACY unit name on the initialization list by a list of unit names into which it expands. The order reflects the order of the units given in the macro definition. The two approaches are essentially equivalent in this case. Alternatively the creator of the REACRECY macro may choose to prepare a special initialization routine for this unit. This special initialization routine can take advantage of the known nature of the unit to estimate initial values for it. In this latter case the two options would give quite different answers. Initialization routines have proven to be very useful. In most cases they provide an initial guess from which the Newton-Raphson equations have no trouble converging. The initialization routines can be called at any time after the user initializes the appropriate variables-even midway through a calculation. For the flowsheet of Fig.3 the user may solve the original problem then, perhaps, set stream S1 at a different value. He could re-initialize the flowsheet based on this value by calling the initialization routines, which may provide a better initial guess than using the values at the previous solution. Often when a flowsheet is being initialized the user does not want some of the variable values to be set by the initialization routines. This situation can occur for several reasons: eThe user is initializing most of a flowsheet with stored values from a previous run. The retrieved values are closer to the solution than the values which the initialization routine would calculate. eThe user has specified some of the degrees of freedom in the flowsheet and does not want their values changed. He could be initializing a new flowsheet or re-initializing a solved flowsheet in which one or more of the fixed variables have been changed slightly.

4. DATA CENTERED

One could look at ASCEND-II as a flowsheeting system which allows the user to perform several operations on a file of steady state data, as shown in Fig. 4. This program design philosophy has turned out to be a simple and convenient way to permit us to create the ASCEND-II system and to make it truly flexible to use. How does one design and create complex software in a way that gives the user freedom to control the details of input, solution, and output? Most systems require that a set procedure be followed when solving a problem. ASCEND-II is different in that the user is able to control and mix the order of the phases of input, the solving of a problem, and the generating of output. ASCEND-II is actually several independent programs which operate on a file of variable values and pointers set up to describe a flowsheet. By setting up ASCEND-II in this way, one has great freedom to add and delete new commands. The commands do not interfere with one another, and for the most part, are not dependent upon each other. 5. DISCIPLINE INDEPENDENCE--THE CONTEXT SETI'ING FILE OF INPUT

ASCEND-II is an equation solver. The executive system does not know what kinds of equations it is

mass mum.

]

restofs

I

,.,,,a,,.

.o.,..,..

L

\

/

1

,n,", --r

.°,o

[ constrained I derivatives J

Ll,c. o I

'o'J

i ,qu.

rlmldual~

I Nut

L

~..~

Fig. 4. Operations on a steady-state file.

~l~t l

The ASCEND-II system--A flowsheeting application solving. ASCEND-II could just as well solve algebraic equations modeling Civil or Mechanical En#oneering systems. Two things must be done to set up the environment for the calculations of a particular discipline: The Fortran subroutines describing the model equations must be added to the program, and an input file which contains definitions for those models must be assembled. The file containing definitions is called a

context setting file. The casual user of ASCEND-II would never see this file. It would be assembled by those in charge of maintaining the system. The file contains several types of information: definitions of variable packets, definitions of equation packets, types of variable and equation packets associated with a unit, system supplied Macro Expansion templates, and system supplied Default Cards. In Section 7 we shall discuss the details of this input further. We discuss this input in more general terms here. 5.1 Variable packets A problem which had to be overcome in designing ASCEND-II was the problem of accessing variables in a meaningful way. A large flowsheet may contain 10,000 variables. We wanted to include in ASCENDII a type of user interaction which allowed him to access each and every individual variable easily. How could these variables be grouped so that both the user and the system could easily keep track of them? It would be absurd to ask the user to think up individual names for each variable. Instead, the user is allowed to give names to variable packets, which are associated groups of variables. The grouping is arbitrary and is established by the (system's) programmers who write the unit routines to be included in the system. While a flowsheet may contain many variable packets, there are not a great many variable packet types in a flowsheet. With the units we have currently included, all conventional streams, for instance, are variable packets of type STREAM. If a flowsheet were to contain 5000 variables, it might only contain 100 to 200 variable packets of perhaps some 10-20 different types. Definitions of variable packet types reside in the context setting file. Packets may contain (at present) either all double precision variables or all logical variables. Several pieces of information are stored for each variable in a variable packet definition. First, the name of the variable is #oven. The user refers to this name when accessing it. Also stored for each variable is a word telling whether the variable is dimensioned. Mole fractions, for instance, are dimensioned to the number of components in a stream. Default calculate/fix flags for each variable are also specified here. These flags free the user from having to specify the flag for every variable. He only must specify the exceptions to the defaults. Lastly, values are given for default upper and lower bounds of the variable. Once again, the user can override these values. The information in the variable packet definitions is stored once for each variable packet type, not once for each variable packet. A point is set from each variable packet to the vector holding the names of the variables in it. The names of the variables are stored during an entire run, so the user can print the variable

621

values along with their names at any time. The convenience of this total access cannot be overemphasized. Playing "detective" with a problem not behaving itself makes this access indepensible. 5.2 Equation packets Information for equation packets is established in the "context setting file" which is analogous to the information established for variable packet types. Associated with each equation is a name. The value associated with the equation is the residual of that equation. The grouping of equations is also done by the system's programmer to form equation packets and is similar to the grouping of variables. A user can print out equation residuals by naming an equation packet. Since each equation is (usually) associated with a unit, it is easy to tell which units of a flowsheet have converged and which have not by printing out equation packets, allowing the user to pinpoint problems. Again, we did not fully appreciate the value of this access when designing ASCEND-II, but we do now after having solved several flowsheeting problems. Besides containing equation names, equation packet type definition statements specify which equations are dimensioned. As with variables, equations can have more than one dimension. 5.3 Variable and equation packets associated with a

unit Part of the context setting file input defines what types of equation and variable packets are associated with a unit and determines the order in which they should appear in the user's input. Its use to discover user input errors is obvious. Also it is needed so the system knows where to generate default names for packets not explicitly named by the user. 5.4 Macro expansion templates This context setting file also contains entities called "Macro" definitions created by the system's programmer which allow the user to create complex models from simple input. The user can add his own macro expansion templates if he chooses. These would be part of his input which follows the input from the context setting file. 5.5 Default cards The final section of the context setting file is a description of the standard default cards included with the system. User input, as we have seen, is in the form of naming things (units, variable packets, etc), and defaulting significantly reduces much of this input that the user might otherwise need to enter to define a problem. As with macro expansion, the user can add his own defaults to the system at the start of his input file. 5.6 Advantages of a context setting file Adding this typing information in the context setting file adds flexibility to the system while saving time and space. The user can include only those definitions relevant to his problem. If he has no input of a particular type, he would not be required to include its definition, saving time in reading and processing the input, and freeing space in the computer.

622

M.H. LOCKEand A. W. WESTERnF.~O

The method also allows for changes in these definitions without changing computer code, removing the necessity of recompiling and linking the program when small changes are made. This approach has proven to be especially useful and time saving when debugging new macro expansion templates. Changes can be made in the templates without changing computer code.

as those to start a new problem, stop running ASCEND-II, etc. The final group of commands allow the user to perform optimization calculations. They deal with inputting, saving and setting upper a n d lower bounds, running the optimization, and setting the objective variable.

6. OVERVIEW OF THE COMMANDS

Developing a Newton-Raphson based flowsheeting package seems a simple enough idea. One of the reasons that people have been skeptical about such an approach to flowsheeting is the significant problem of internally representing a flowsheet in a manner which is also convenient to the user. A flowsheet may contain several hundred to several thousand variables. How can they be stored in a manner such that they are easily accessed by the model subroutines, while also being grouped conveniently for the user? We have found that the very simple but very flexible GEV (Generators, Equation Packets, Variable Packets) representation is a solution to this problem (Berna[20]). This structure made implementing ASCEND-II very straight forward. Indeed we feel it is a major contribution in concept and worth describing in some detail here.

Tables 1-7 briefly describe the commands now available in ASCEND-II. The commands are divided into groups based on their functions. The first group of commands appear in the context setting file. These commands define the unit models, variable packets, and equation packets. The second group of commands contains those necessary to set up a problem. These commands perform functions such as processing the user's input and checking that each variable has been initialized. The third group of commands includes those used to input a problem. These commands allow the user to input his flowsheet and variable values. Once a user has input his flowsheet and initialized variables, he is ready to converge the equations. He does this with the commands in the fourth group. The fifth group of commands is used for report generation. Variable values and flowsheet descriptions can be displayed to the user at the terminal or saved in a file. The sixth group of commands includes miscellaneous commands, such

7. THE GEV REPRESENTATION FOR A FLOWSHEETING PROBLEM

7.1 Three representations of a flowsheet Figure 5 shows three representations of a flowsheet. The top diagram is the user's flowsheet.

Table 1. Commands in the context setting file DEFAUL

A l l o w s default cards to be input

EQTYP

Equation packet definitions f o l l o w

EXFORM

Macro expansion templates f o l l o w

WPEPK

Types of variable and equation packets associated with a unit are entered here

VARTYP

Variable packet definitions f o l l o w

Table 2. Commands to set up a problem EXPAND

Perform macro expansion, set up internal GEV structure

PREINT

Get ready to start initializing variables

AFTINT

Check that all variables have been initialized

PRESOL

Prepare necessary information for solving or optimizing

The ASCEND-II system--A flowsheeting application

623

Table 3. Commands for user input Read input/send output to a file

RESET

Integer seeds f o r system-created

SEEDS

variable packet, equation packet, and unit names f o l l o w

FSHEET

Flowsheet description f o l l o w s

VPCONX

Pointer system linking variable packets to other variable packets f o l l o w s Initial variable values and/or

INVAL

calculate/fix flags f o l l o w

•FETCH

Fetch variable values f r o m a file

GETCF

Retrieve calculate/fix flags f r o m a file

SEQUEN

Initialization sequence is given here

RUNVAL

Variable values and/or calculate/fix flags f o l l o w

This flowsheet represents a flash unit with a feed stream named SF, a vapor stream named SV, and a liquid stream named SL. The middle diagram is the internal representation of the flowsheet. This is called a Generator, Equation packet, Variable packet (GEV) Diagram [20]. Generators are subroutines. They calculate partial derivatives and equation residuals necessary for the Newton-Raphson based solution and optimization scheme. Equation packets (represented by vertical lines---think "rows of the Jacobian matrix") are groups of equations. They are grouped together because they may represent the equations of the

Jacobian matrix reserved for a single generator. Variable packets (horizontal lines--think "columns of the Jacobian matrix") are groups of associated variables, grouped together for convenience. One type of variable packet in our current implementation is a stream variable packet. It contains the flowrate, temperature, pressure, enthalpy, entropy, and component mole fractions of the stream. The person who writes the generators, a system's programmer for instance, defines the equation and variable packet when developing the unit models for the system. The lower diagram in Fig. 5 is the sparse Jacobian representation of the flowsheet in the top figure. The

Table 4. Commands to solve a problem

UNITIN SETWTS

Call the initialization routines Set variable and equation scaling factors

RHMAG

Calculate the norm o f the equation residual vector

SOLVE

Start the Newton-Raphson solution procedure

CONDER

Calculate sensitivity o f one variable with respect to another

M. H. LOCKEand A. W. W m ~ . B ~ G

624

Table 5. Commands for report generation PRINT

Print variable values

RHSDS

Print equation residuals

SAEXPD

Save fully expanded version o f the input in a file or display to terminal

SAINTR

Save or display intermediate version o f input

Save or display original

SAORIG

version o f input

Save calculate/fix flags in a format

SAVECF

that can be read by the GETCF command

Save variable values in a format that

SAVEVP

can be read by the FETCH command

equation packet FLEP contains equations MB (material balance), EQ (equilibrium) and ENTH (enthalpy balance). These equation correspond to particular rows of the Jacobian. The columns of the Jacobian are related to the variable packets. Tying the equation packets to the Jacobian tells the generator for which rows of the Jacobian it generates partial derivatives and residuals. The variable packets tied to the generator tell it for which columns it is generating partial derivatives.

A major "commandment" for ASCEND-II is that every number is a variable. There are no "parameters" as distinct from "variables" as distinct from "constants" in ASCEND-II. 7.2 Flowsheet input A flowsheet can be characterized as a set of generators with associated variable packets and equation packets. The input reflects this observation. It is a

Table 6. Miscellaneous commands

ALLOW

A l l o w s the user to give commands neglecting any necessary order

EOF

Last command in an input file, returns control to the user at the terminal

HELP

Prints list o f currently valid commands

LOOKUP

Tells variable name and name o f packet in which it resides for specified variable integer position

REWIND START

A l l o w s user to rewind an input file

Re-initializes necessary values to run a new f lo w sh e e t

STOP

Stops the ASCEND-II system

The ASCEND-II system--A flowsheeting application

625

• Table 7. Commands for optimization Perform optimization according to

OPTMIZ

modified Han-Powell algorithm

GETLO

Retrieve lower bounds from a file

GETUP

Retrieve upper bounds from a file

SAVELO

Save lower bounds in a file

SAVEUP

Save upper bounds in a file

OBJECT

Tell the system what variable is to be minimized or maximized

LOWBND

Input lower bounds at the terminal

UPBND

Input upper bounds at the terminal

VIOLAT

Report which upper or lower bounds have been violated

• SV

SL

SF MB

X

EQ ENTH

X

SV

SL

X

X

X

X

X

X

Fig. 5. Three representations of a flowsheet.

/ery general and easily comprehended format for nput. Figure 6 shows a flowsheet consisting of two ;tage units. Below the flowshect is the GEV diagram vhich represents it. The ASCEND-II statement sets needed to model he flowsheet of Fig. 6 are shown below.

C* C* C* FS V G E C* C* C* C* FS V G E C*

First, the lower STAGE unit:

Next unit is STGI, a stage unit. Variable packets are $5, $6, S I and $2. The subroutine to model the unit is "STAGE". The equation packet for the unit is called EPSTG 1. Now, the upper STAGE unit:

STGI STAGE $5 $6 SI $2 STAGE EPSTG1

STG2 STAGE $4 $5 $2 $3 STAGE EPSTG2

Note that the input exactly reflects the GEV representation. Not illustrated are " N - • e s " which permit single variables rather than whole variable packets to be referenced, "I-lines" which permit integer information to be entered (for example, the number of stages in the top section of a distillation column), and "D-lines" which allow one to pass along default card information when a unit macro expands. 7.2.1 System generated names. Often a user does not wish to give names to untis, variable packets, or equation packets. This situation occurs when a user feels he may never access an equation or variable packet, and most frequently, when units are created by Macro Expansion. If a word is left blank, ASCEND-II will fill in the blank with a systemcreated name. A "$" as the first character indicates the name was generated by the system rather than suppfied by the user. 7.2.2 Hierarchy in naming words. There is a definite hierarchy in the way that words are filled in

626

M.H. LOCKEand A. W. WESa'E~m~G

$4

$5 $8

S6

SI I

VLV

¢,<1

$7

~~FLASH

oj

T s9

PMP

Fig. 6. Absorber with flash unit. on the FS statement set. Words entered by the user have the highest priority. If a word is left blank by the user, the system attempts to fill it in with the first then second then third default statement set named (words 4--6 on the user's FS line). If the word is still blank, the system will generate a name for that word. 7.2.3 Macro expansion. Macro Expansion is used to create complex structures from simple input. It is part of the input language to ASCEND-II as a convenience for the user. Suppose, for example, that in Fig. 6 two stages were not enough to accomplish the required split. The user may find that five stages are necessary. He may replace the two stages with five. His input for the five-stages (with equation packet names defaulted) might look like:

FS ST5 STAGE V S16 SI $6 S17 G STAGE Inputting these five stages separately is tedious work, and it is easy to make a mistake. Changing the number of stages requires making a number of changes to the flowsheeting input. If stage ST3 is later removed, the user must be sure that the streams connecting ST2 to ST4 match up correctly. Instead of inputting the stages one at a time, the user can input a unit of type SECTN. SECTN is recognized by the current version of ASCEND-II as a unit type which "expands" through Macro Expansion, and it would then be expanded into the individual STAGE units, creating names for the stages, variable packets, and equation packets, while also insuring that the streams connecting the stages are in the appropriate positions. Instead of inputting five individual STAGE units as in the previous example, the user would input the SECTN card:

FS V G

STI STAGE S 3 S 1 0 S l l S4 STAGE

FS V G

ST2 STAGE S10 S12 S13 SII STAGE

FS V

SECTI SECTN $3 SI $6 $4

I

5

FS V G

ST3 STAGE S12 S14 S15 S13 STAGE

FS V G

ST4 STAGE S14 S16 S17 S15 STAGE

7.3 Macro expansion of a distillation column One of the most commonly used pieces of equipment in chemical plants is the distillation column. ASCEND-II models a distillation column as a collection of mixers, flash units, stage models, and a splitter, as shown in Fig. 9. Since this model is used so often, it was worth the time to write a macro to

$2

~

$5

s6

['q'~.,

,

Fig. 7. GEV representation of absorber with flash.

The ASCEND-II systern--A flowsheeting application

s3 l ONOO, ] $ VPO01~'-~

VPO02

~

VPO04

$ vPoo3

I ,o,oo, 1

$

VPO05 "-'~~VPO06

E0o,oo, 1

$ VPO07~ - ~ "

vPO08

F,o,ooo 1 Sl ~-'~'~S6

Fig. 8. Absorber created by macro expansion. model it. The user input to the column model is far simpler than entering statement sets for each of the units making up the column individually. The user enters a COLUM card, specifying variable packet names for tray efficiency and number of stage models for the top and bottom column sections, and naming the inlet and outlet streams of the column, plus a variable packet which contains the reflux ratio. He has the option of naming some of the intermediate streams in the column. If he doesn't name them, the system will generate those names. The user's COLUM input expands into three flash units, two mixers, a splitter, and two column sections.

COLUMN

SECTION

=~FLASH~ MIX "[ 3M'X COLUMNI SECTION

FLASH Fig. 9. Units grouped to form a distillation column. CACE 7:5-E

627

Once this expansion is done, a second stage of Macro Expansion is required because the units into which the column first expands are not basic units and must themselves be macro expanded. The flashes, mixers, and splitter further expand into the generators to do the material and equilibrium calculations, plus the generators to force the sums of the mole fractions of each outlet stream to be 1.0. Also at this time, each column section expands into the number of individual stage units specified on the I line of the column input. In the last (third) stage of expansion, each of the stage units created through Macro Expansion of the column sections expands into stage generators (material and equilibrium equations) plus sum of mole fractions generators for the stage outlet streams. Had complex physical property models been a part of this model, further levels of expansion would occur, i.e. until the units are all basic units. 7.4 Macro expansion with default cards In the previous example (Macro Expansion of a COLUM statement set) the types of the units created by expansion were set in the expansion template. The column always expands into three F L U N statement sets, two SECTN statement sets, one SPLUN statement set, and two MIXUN statement sets. There is no other choice in the expansion of a COLUM if the Macro Expansion template is written as it is. Using the D-line to pass default cards permits the creation of OR branches in a Macro Expansion template because a unit is able to pick up its type from a default card name. 7.5 Routing physical properties An important issue in the design of a chemical engineering simulation program is allowing the user to choose easily among different physical property options. The temperature, pressure, and composition of a stream determine which physical property calculations should be done. Older simulation programs, such as FLOWTRAN[21] require that the same physical property option be done on every stream. Often it is useful to do different physical property calculations on different process streams. The internal data structure of ASCEND-II, represented by the GEV diagram, is flexible enough to allow different calculations on different streams. A routing tree (Seider, et al.[22]) can be used to represent the physical property options which can be selected from among the options available. In this example, we will not go into detail describing the routing process, but will show the type of tree structure which can be created by ASCEND-II to route physical properties. In this example, suppose the user has the choice of two equations of state; the Virial equation and the Prausnitz-Chueh version of the Redlich-Kwong equation (PCRK). (Only the PCRK equation is now in ASCEND-II.) The equation of state is used in the calculation of vapor phase fugacity coefficient, enthalpy, and entropy. Suppose also that he has a choice of the Wilson equation or the Van Laar equation for liquid phase activity coefficients. (Again only the Wilson equation is actually in ASCEND-II.) A third choice is to assume ideal conditions and use no equation of state and/or not calculate liquid phase activity coefficients.

628

M . H . LOCKE and A. W. W n s ~ o FLUN I

AND KVAL

OR

!

PCRKEX

I

I

I

I

fD'ES

AND

AND

I

I

I

i

STRM

PCRKV ENTHPC ENTPC i AND

S I"~RM

PCRKV ENTHPC ENTPC

I

VAF~ID

I

AND

I I

I

I

LIQID

AND

1

,&c

VII~LEX

FuC'OPCLOPC.K VPCRK I t ,1

I

WlLSN VANLR

AND

I

ANTON

OR

FLASH

1

I

I

,

STRM ENTHID

ENTID

I

AND !

STRM I FUCOVI LQVIRL VPV(RL ,

AND r

STRM

I

I

!

ENTHID ENTIO

I

AND

V~RL ENTHVI ENTVI STRM

V~RV ENTHVI ENTVI

Fig. 10. Macro expansion tree for a FLUN unit.

An OR branch of the routing tree can be represented by Macro Expansion used with a Default Card, where the unit type is picked up from the Default Card. An AND branch is created by Macro Expansion without a Default Card. Figure 10 shows the type of routing tree which can be created by ASCEND-II. Any branch preceded by the word AND is produced purely by Macro Expansion. Any OR branch is produced by Macro Expansion with Default Cards to set the unit type. 8. SAMPLE PROBLEMS

8.1 Optimization of a distillation column The first test problem involves the eventual optimization of a distillation column. The word eventual is used because the evolutionary strategy was used to get to the solution. The distillation column structure is shown in Fig. 9. This particular column has 5 stages in each of the

two column sections. The feed is a 3-component mixture with relative volatifities for components A, B, and C of 3, 2., and 1., respectively. The feed is set at 100moles/hr. with a composition of 30% A, 30% B, and 40% C. The objective is to minimize the composition of the heaviest component, component C, in the distillate, stream S2. For the entire series of calculations, the flow of stream S12 is specified at 0.0 (all liquid feed), and the flow of $4 is also set to 0.0 (total condenser). Constraints other than the fact that all flowrates, mole fractions, and split fractions must be non-negative, will be discussed later. The column is modeled by ASCEND-II with 131 equations, most of them nonlinear. The way to find the optimal operating conditions of this column is to evolve towards the solution by first solving the column as a simulation problem, then moving to the optimization calculation. In the simulation calculation, the feed is fixed, the reflux rate is

Table 8. Summary of distillation optimization example

Tray Efficiency

Flow of $2

Objective

Comment

1. (F)

1. (F)

50. (F}

.0627458

Base Case

1. (F)

1. (F)

25.9987 (C)

.0517580

Optimum

1. (F)

1. (F)

26.2 (F)

.0517583

Verification

1. (F)

1. (F)

25.8 (F)

.0517563

Verification

1. (F)

1.62385 (C)

21.5233 (C)

.0368205

Optimum

9. (C)

1. (F)

23.9651 (C)

.00104832

Optimum

g. (F)

1. (F)

24.1 (F)

.00104837

Verification

9. (F)

1. (F)

23.9 (F)

.00104834

Verification

Reflux Rate

The ASCEND-II system--A flowsheeting application set to 1.0, the number of trays top and bottom are fixed, and the flowrate of the distillate stream (stream $2) is fixed at 50.0 moles/hr. With these specifications ASCEND-II has always found a solution to the resulting set of equations. This solution can then be used as a base case for further calculations. At these specifications, the composition of component C in stream $2 is found to be 0.0627458 (excuse the accuracy shown). From thxs base point, the user may run some experiments. He can show himself that raising the Murphree stage efficiency (has an effect similar to increasing the number of stages (Brosilow, et aL [23]) or the reflux rate lowers the objective mole fraction. Lowering the flow of stream $2 seems also to lower that mole fraction. In a few minutes the user can get a feel of the way things should operate. When he does turn on the optimizer, he will be able to tell whether or not the solution makes sense. While running the experiments off the base case, the user could perform a design calculation. He could, for instance, release the specification on the flow of $2 or on the reflux rate and specify the composition of component C in $2. This switch can be done interactively in a few seconds without changing the flowsheet in any way. He can lower the mole fraction from the base case and see the effect on the other variables. Eventually a point will be reached where he has specified a mole fraction too low: in this case the equations will not converge or a solution with negative flows will be found. Once the user has learned how the variables behave, he is ready to turn on the optimizer. Several optimization calculations were performed from the base case of this flowsheet. In the first run, the reflux rate is fixed at 1.0, the tray efficiency is fixed at 1.0, and the specification of the flowrate of stream $2 is released. ASCEND-II finds an internal optimum with respect to the flows of $2 at a value of 25.9987 moles/hr. The mole fraction of component C in $2 is 0.0517580. The user may now ask himself the question, "How do I know that this is really the optimal flow of $2?". Using ASCEND-II he can demonstrate to himself the validity of this solution. Interactively specifying the flowrate of $2 at 26.2 gives a mole fraction of 0.0517583, higher than the optimum found by ASCEND-II. Lowering $2 to 25.8 also gives a mole fraction of.0517583. These runs can be done in 2 or 3 min--while the question being asked is still fresh in the user's mind. At this point, the user may wish to release another degree of freedom, for instance the tray efficiency of the top column section. Now there are two degrees of freedom, tray efficiency and flowrate of $2. With these two variables to work with, the optimizer seems not to converge. Rather tray efficiency increases without bound. The user can conjecture that more trays will always improve the objective, which is intuitively correct. The user may now wish to free the specification on the reflux rate while fixing the tray efficiency. In the next run, he does that, but fixes the tray efficiency at 1.0 and puts an upper bound of 9.0 on the reflux rate. ASCEND-II pushes the reflux rate to its upper bound, and finds an internal optimum with respect to the flowrate of $2 at 23.9651 moles/hr. The target mole fraction is lowered to 0.00104832, the lowest so

629

far. Verification as before proves that this is the optimum for this problem, an optimum that seems counterintuitive. Raising $2 to 24.1 gives a mole fraction of 0.00104837. Lowering it to 23.9 gives 0.00104834. Clearly both are slightly higher than the optimum. The internal minimum on distillate product flow can be explained. We leave it to the reader to 1) decide if he believes it and 2) explain it. We have optimized many flowsheets of modest size (fewer than 1000 equations), each having from one to five degrees of freedom. Typically 10-20 total iterations are required for convergence. When not optimizing, the number of iterations is typically 4-6, so it appears the optimization algorithm used takes only about three times the number of iterations as needed for a design calculation--a remarkably low number. When convergence for these simple problems has not been achieved, we have always found we messed up the problem definition--an encouraging result. The examples have not used complex physical property models. Incidently, we have taken care with respect to the scaling of variables and equations and the avoiding of division as an arithmetic operator either in evaluation of errors or derivatives (see Westerberg Is for details on these ideas). Also all calculations are done in double precision (DEC-20 computer). While seemingly mundane ideas, they pay handsomely we believe in aiding to get problem convergence. This paper is based on the PhD thesis of Locke [24] and on an addenda to that thesis produced during a year of postdoctoral work. The reader is referred to these two documents for a detailed description of ASCEND-II. Acknowledgements--This work was funded by NSF Grant CPE 780189 and by the Exxon Corporation. REFERENCES

1. P. Friedman & K. Pinder, Optimization of a simulation model of a chemical plant. Ind. Engng Chem. Proc. Des. Dev. 11, 512 (1972). 2. L. T. Biegler & R. R. Hughes, Approximation programming of chemical processes with Q/LAP. Chemical Engng Progress (April 1981). 3. L. T. Biegler & R. R. Hughes, Infeasible path optimization with sequential modular simulators. AIChE J. 1982 (Accepted for Publication). 4. M. J. D. Powell, A fast algorithm for nonlinearly constrained optimization calculations. Presented at 1977 Dundee Conf. on Numerical Analysis (1977). 5. M. Shacham, S. Macchietto, L. F. Stutzman & P. Babcock, Equation oriented approach to process flowsheeting. Comput. Chem. Engng 6(2), 79-95 0982). 6. E. W. Gorczynski & H. P. Hutchison, Towards a quasilinear process simulator--I. Fundamental ideas. Comput. Chem. Engng 2, 159-196 (1978). 7. R.W.H. Sargent & J. D. Perkins, A computer program for steady-state and dynamic simulation and design of chemical processes. Paper 49a, Annual Meeting, AIChE, New Orleans, Louisiana (1981). 8. M. A. Stadtherr & C. M. Hilton, Development of a new equation-based process flowsheeting system. Paper 49b, Annual Meeting, AIChE, New Orleans, Louisiana (1981). 9. P. D. Babcock, Equation-oriented flowsheet simulation: Opening a new world of process simulation. Paper 49c, Annual Meeting, AIChE, New Orleans, Louisiana (1981).

630

M . H . LOCKEand A. W. WESa'F~BEgG

10. T. J. Begna, M. H. Locke, & A. W. Westberg, A new approach to optimization of chemical processes. AIChE J. 26(0, 37--44 (Jan. 1980). 11. M. H. Locke, R. H. Edahl, & A. W. Westerberg, An improved successive quadratic programming optimization algorithm for engineering design problems. Accepted for Publication, 1982. 12. B. A. Murtagh, On the simultaneous solution and optimization of large-scale engineering systems. Cornput. Chem. Engng 6(1), 1-5 (1982). 13. B. A. Murtagh & M. A. Saunders, MINOS/ AUGMENTED User's Manual. SOL 80-14, Systems Optimization Laboratory, Stanford University (1980). 14. T. J. Berna, M. H. Locke & A. W. Westerberg, Optimization of large chemical process systems. Large Engng Systems (Edited by G. J. Savage and P. H. Roe), Vol. 2, pp. 157-161. Sandford Educational Press (May 1978). 15. R. H. Edahl, Sensitivity Analysis in Nonlinear Programming. Ph.D. Dissertation. Carnegie-Mellon University, Pittsburgh, Pennsylvania (1982). 16. D. Gabay, Reduced quasi-Newton Methods with feasibility improvement for nonlinearly constrained optimization. Math. Prog. Study 16, 18-44 (1982).

17. I. S. Duff& J. K. Reid, Some design features of a sparse matrix code. Trans. Math. Software 5, 18 (1979). 18. A. W. Westerberg, Close encounters of a sparse kind. Chem. Engng Education Spring, 72-77 (1980). 19. S. Kuru, Dynamic Simulation with an Equation Based Flowsheeting System. Ph.D. Dissertation, CarnegieMellon University, Pittsburgh, Pennsylvania (1981). 20. T. J. Berna, A Newton-Raphson Based Procedure for Optimization of Large Chemical Processes. Ph.D. Dissertation, Carnegie-MellonUniversity, Pittsburgh, Pennsylvania (1979). 21. J. D. Seader, W. D. Seider and A. C. Pauls, FLOWTRAN Simulation--An Introduction. CACHE Corporation (1974). 22. W. S. Seider, U B. Evans, B. Joseph & E. Wong, Routing of calculation in process simulation. I&EC Process Design Develop. 18, 292-297 (1979). 23. C. Brosilow, R. Tanner & H. Tureff, Optimization of staged counter current processes. 63rd National Meet ing, AIChE, St. Louis, Missouri (1968). 24. M. H. Locke, A CAD Tool which Accommodates an Evolutionary Strategy in Engineering Design Calculations. Ph.D. Dissertation, Carnegie-Mellon University (Sept. 1981).