Oredo2Swlpl) Dwlrqri&Khplfdo (Qjlqhhulqj 3Odqwve/Phdqvri (Yroxwlrqdu/$Ojrulwkpv
Oredo2Swlpl) Dwlrqri&Khplfdo (Qjlqhhulqj 3Odqwve/Phdqvri (Yroxwlrqdu/$Ojrulwkpv
Oredo2Swlpl) Dwlrqri&Khplfdo (Qjlqhhulqj 3Odqwve/Phdqvri (Yroxwlrqdu/$Ojrulwkpv
3ODQWVE\PHDQVRI(YROXWLRQDU\$OJRULWKPV
Michael Emmerich*, Bernd Gro+, Frank Henrich+, Peter Roosen+, Martin Schtz*
*
ABSTRACT:
The availability of high speed parallel computers and process simulators capable of rigorous
flowsheet calculation has given way to new methodologies for the optimization of chemical processes.
Evolutionary algorithms (EAs), derived from biological adaptation paradigms, have proven to be a
powerful and robust optimizing technique even for structural optimization problems.
We have adopted this method for constrained mixed-integer and structural process optimizations and
therefore developed EA-Modules with an interface to the ASPEN PLUS simulation system, thereby
combining the evolutionary optimization with the simulators detailed process models and algorithms
for realistic cost target function calculations.
The method has been applied successfully to typical test problems, such as heat exchanger network
and distillation sequence synthesis and even the optimization of whole processes.
Introduction
The optimization of chemical engineering plants is still a challenging task. Target functions
are often nonlinear. Discontinuities, multi-modalities and a strong interaction between
flowsheet parameters can often be found. Furthermore the optimizer has to deal with explicit
and implicit restrictions and different types of variables (i.e. real-valued, integer or discrete).
The overall process optimization comprises the search for the optimal process structure as
well as the adjustment of the unit operation parameters.
There are many papers published that describe the application of mathematical programming
techniques on simplified models for chemical engineering plants, or their subsystems (e.g.
heat exchanger networks or separation sequences). An overview can be found in Grossmann
1996.
In many cases the optimization of subsystems, like separation sequences and heat exchanger
networks, is considered on its own. Simplifying assumptions like fixed feed conditions or
fixed overall pressure make it difficult to combine these approaches for an overall plant
optimization, because there is often a correlation between parameters of different subsystems.
Most of the mathematical programming methods are based on gradient evaluations and
therefore do not work on unsteady topologies. Therefore they estimate capital costs as
continuous functions of equipment capacities such as exponential functions. Realistic system
designs must select components offered in discrete sizes as Maia et al. (1995) states.
The simulator ASPEN PLUS enables the detailed economical evaluation of a process
flowsheet. The inclusion of the process structure in an optimization step is not possible
without modifying the structure by hand. Even the parameter
optimization is a difficult task, e.g. if economical target functions are chosen with discrete
item sizes resulting in an unsteady target function.
In this paper a combination of a rigorous simulation model and a robust optimization
algorithm, which enables the optimization of difficult target functions, is presented. A selfdeveloped external optimizer based on Evolutionary Algorithms (EAs) is coupled to the
process simulator by a modular and extendible interface.
After a general introduction of EAs for numerical optimization, two applications are
introduced that illustrate the use of EAs in chemical plant design.
The first example deals with the optimization of heat exchanger networks with EAs and the
second example with the overall optimization of a chemical engineering plant, applying
ASPEN PLUS as a flowsheet simulator. Finally we give an overview on further applications
of EAs in chemical engineering plant design.
Evolutionary Algorithms
EAs attempt to find the best solution for a problem by applying basic principles of natural
evolution processes. The term EA include genetic algorithms (GA), Evolutionary
Programming (EP) and evolution strategies (ES). In this paper we mainly refer to Evolution
Strategies (ES) that have been mainly developed by Rechenberg and Schwefel (c.f. Schwefel
1995) to solve numerical optimization problems.
The basic principles of EAs were originated in the late fifties and early sixties. In the last two
decades a broad acceptance and use for engineering applications can be reported. Several
hundred articles are published now annually. There are at least twenty regularly scheduled
conferences on the topic of Evolutionary Algorithms.
There are some main elements common to practically any
existing EA.
Modern Evolution Strategies operate with several instead of just one search point in the
search space at one time. Search points, together with some additional information related to
them, are called individuals. A set of individuals existing at a given time is called population.
Although it is not the correct biological term, the objective function is called fitness
function.
An individual consists of two main components. The object variables that are parameters of
the objective function, and the strategy parameters that are only used to control the
evolutionary operators. In vector representations strategy parameters often encode standard
deviations for the variation of the object variables.
In EAs the variation of individuals is carried out by mutation and recombination operators.
A mutation operator changes the object variables of a search point randomly. The strength of
the mutation is controlled by strategy parameters. The recombination operator generates a
new individual from two randomly chosen individuals, by mixing their information. Today,
there are plenty of variation procedures that have been developed and their choice strongly
depends on characteristics of the specific problem, like the (assumed) topology of the fitness
function and the types of variables and constraints.
In this paper EAs are applied that are similar to the Evolution Strategy (ES) as published in
Schwefel 1995. The main loop of a so-called (+)-Evolution Strategy is described in the
following:
t := 0
initialize 3RSXODWLRQ P( 0 )RI LQGLYLGXDO V
evaluate 3RSXODWLRQ P( 0 )
while WHUPLQDWLR QFULWHULRQQRWIXOILOOHG
recombine RIIVSULQJ LQGLYLGXDO VRXWRI SDUHQWV
mutate RIIVSULQJ
evaluate RIIVSULQJ
select WKH EHVWLQGLYLGXDO VIRU P (t + 1)
t = t + 1
By introducing a maximal life span for the individuals we get the so-called (, ,)Evolution Strategy.
After a random initialization and evaluation of the first individuals, offspring individuals are
generated step-wise with a recombination operator and a mutation operator. Then the fitness
of the offspring is evaluated. The selection operator chooses the best individuals out of the
offspring individuals and the parental individuals not exceeding the maximal life-span (age).
As long as the termination criterion is not fulfilled, these selected individuals form the
parental generation for the next iteration loop.
For a more sophisticated description of this algorithm the reader is referred to Schwefel 1995.
Many other variants of this algorithms scheme are developed, e.g. multi-population Evolution
Strategies and Evolution Strategies that work on spatial distributed populations. For a
overview the reader is referred to Bck et al. 1996.
Notice, that the Evolution Strategy needs no explicit formulation of the target function, which
allows the use of numerical simulator models (black-box evaluation).
The power and flexibility of this optimum seeking method has been demonstrated in many
papers and books. Schwefel compared a standard evolution strategy with the most common
numerical optimization algorithms (c.f. Schwefel 1995). In this investigations they proved to
be a powerful optimization methodology especially for difficult target functions with a high
number of decision variables.
A similar comparison for mixed integer optimization benchmark problems has been carried
out by Gro 1999. The 21 problems have been taken from Ryoo and Sahinidis 1995 and most
of them are related to chemical engineering. They contain up to 10 constraints and decision
variables. Although a simple ES approach has been applied, all problems could be solved in a
few seconds on a PENTIUM 133 processor. Other optimization strategies may need less
target function evaluations but it should be noted that no problem specific knowledge and no
starting point variations had to be introduced into the solution process. In addition to this, the
time consumption of EAs increases less with problem size than that of many other
algorithms.
Program package
In order to solve chemical engineering problems with the numerical simulator ASPEN PLUS
we developed a program package that can be defined as follows:
The optimizer generates solutions that are translated into input files for the simulator (see
Figure 1). The simulator then calculates thermodynamic states of the process. The results are
used to calculate equipment sizes and costs of the flowsheets. The optimizer generates
solution propositions by specifying respective sets of system design variables that are
forwarded to the ASPEN PLUS simulator by an interface program creating the required
Aspen input files.
In the following step simulation errors and constraints are evaluated in a further custom
program module. Whenever restrictions and constraints are violated or simulator errors occur,
a penalty value is added to the process cost, that grows with the severity of the violations or
errors. Furthermore, infeasible solutions are always assessed with higher costs than feasible
one. As the underlying optimization framework is derivation-free this does not impose any
definition problem --- in practice there is always a sensible upper limit on realistic process
costs that can easily be surmounted by the respective penalty value.
Finally the target function value is returned to the evolutionary optimizer, which interprets it
as the fitness value of the proposed set of system values.
3UREOHP'HILQLWLRQ
7UDQVODWLRQLQWRWKH$VSHQ
,QSXWILOH
ASPEN PLUSTM
6LPXODWLRQ&DOFXODWLRQRI
WKHUPRG\QDPLFDOVWDWHV
(YROXWLRQDU\
2SWLPL]HU
&RVW&DOFXODWLRQV
([DPLQDWLRQRI6LPXODWLRQ
(UURUVDQG&RQVWUDLQWV
)ORZVKHHW(YDOXDWLRQ
Figure 1: The interaction between the flowsheet simulator and the optimizer. After the time expensive
simulation a cost function and a penalty term for the simulation errors is calculated. See Gro, Roosen 1999 for
more details concerning the interface to the simulator and the target function.
A very small heat exchanger network, comprising two cold streams and two hot streams, is
illustrated in Figure 2. Aiming at an internal recovery of energy the heat from the hot streams
should be transferred to the cold streams as far as possible. Only the remaining needs must be
supplied by utilities, like saturated steam or cooling water.
320 o F [16666 .7 BTU
]
hF
o
480 F [ 20000
BTU
h oF
63/,75$7,2
+($7(;&+$1*(5
1.3525
= 0 .8
140 o F [14444
320 F
+
&2/'675($0
+27675($0
&
200 o F
280 o F
500 o F
MBTU
h
1.247 MBTU
h
BTU
h oF
2 .751 MBTU
h
(17+$/3<
Figure 2: Example HEN with stream splitting. The red arrows and circles symbolize hot streams and heaters
whereas the blue arrows and circles symbolize cold streams and coolers. Heat exchangers are symbolized by a
blue and a red circle linked with a thick pink line. A more detailed problem description for the 4SP1 problem can
be taken from e.g. Linnhoff 1978.
The cost of the heat exchanger network, that shall be minimized, depends on the operating
costs and on the investment costs of the HEN. Whereas the operating costs mainly depend on
the amount of energy that is exchanged between the streams, the investment costs grow with
the number of processing units and the size of heat transfer areas of the heat-exchangers. This
area-size grows exponentially with the decreasing difference between the temperature levels
of a coupled stream.
As an internal representation of a heat exchanger network for the evolutionary optimizer we
chose a graph representation, with typed nodes and edges to which real valued parameters are
assigned. Each node of the graph represents an unit operation and the edges represent physical
or logical links between the apparatus. The advantage of a graph representation is, that it
allows for the coding of arbitrary structures.
The initial population consists of empty HENs, i.e. networks in which heating and cooling
demands are completely supplied by external sources. There are two kinds of mutation
operators applied for the variation of the HEN graphs. The first one modifies the real-valued
parameters, which are the enthalpies and split-ratios, similar to the standard mutation operator
(c.f. Schwefel 1995) and the other one accomplishes small structural changes of the graph.
The modification of the topology is carried out by adding and deleting heat exchangers,
swapping their position on a stream, inserting or deleting a split. Two modification operators
for the flowsheet are illustrated in Figure 3. A recombination operator has not been applied,
because it is quite impossible to find corresponding parameters in this variable dimensional
problem.
Figure 3: The topology modification of HENs is carried out by graph modification operators that are
applied randomly to the HEN graph and lead to small changes of the graph structure. For every graph
modification a reversible operator needs to be defined to avoid a drift induced by the mutation. The set of
modification operators is chosen in a way that every possible topology can be reached from every search point
by iterated application of the mutation operator (c.f. Emmerich 1999).
There are several constraints that need to be kept in order to get feasible solutions, e.g.
crossing temperature fields have to be avoided. For the optimization of the HENs we applied
a quadratic penalty function, which grows with the severity of the constraint violation.
The structure of the best solution found for the small 4SP1 problem is shown in Gro, Roosen
96. The structure is identical to the solution found in several papers (e.g. Hesselmann 1985,
Linnhoff 1978). In the reported cases, a minimum temperature approach of 20F was assumed
as part of the problem description. If this restriction is realized in our program package the
heat loads of the optimal solution are also identical to the cited literature results. The costs of
this solution are 13587 $/yr. For the application of the ES optimizer this Tmin fixing is not
necessary, since the exponentially growing cost for the heat exchanger areas are automatically
calculated by the target function. If the Tmin constraint is relaxed to 0F an alternative
solution with a temperature difference of only 1.6F and the heat loads given in Figure 2 is
found automatically with total annual costs decreased by 22%. Obviously the additional
expenditure for the increased heat exchanger surfaces is by far surpassed by the savings in the
annualized utility consumptions.
Several further test problems (5SP1, 6SP1, 7SP1, 8SP1 and 10SP1) were solved with very
good target function values, close to the best literature solutions. Again improved solutions
could be found as soon as Tmin restrictions mandatory for the literature solution
methodologies were dropped. A more detailed description of the given algorithm can be
found in Gro et. al. 1996.
As a typical example for a small chemical plant, the well documented Hydrodesalkylation
(HDA) process (c.f. Douglas 1988) has been optimized with the program package, to prove
the applicability of the methods. The aim of the HDA process is the production of benzene
from toluene. Diphenyl is formed in an undesired side reaction.
63/,7
&2035(6625
+<'52*(1
0(7+$1(
72/8(1(
+<'52*(1
0(7+$1(
)/$6+
+($7(5
)851$&(
5($&725
&22/(5
3803
0(7+$1(
%(1=(1(
9$325
7:23+$6(
/,48,'
67$%,/,6(5&2/801
',67,//$7,21
',3+(1</
67$1'$5'&21),*85$7,21
Figure 4: Flowsheet of the HDA process network (taken from Douglas 1988) that has been chosen to
demonstrate the applicability of overall ASPEN PLUS flowsheet optimization with EAs.
Figure 4 shows a possible flowsheet for the HDA process which serves in the following as
reference configuration. There are many other possible configurations for the HDA process.
Structural alternatives of the process can be represented by binary decision variables,
selecting alternative subsystems (see e.g. Gro, Roosen 1998). The application of such
decision variables is demonstrated in Figure 5.
Before an optimization of the flowsheet in question can be carried out, all parameters that are
to be optimized and their feasible range have to be defined. The candidates as parameters can
be divided into three main classes:
metric real parameters (e.g. pressure and temperature levels, conversion rate of the
reactor),
metric integer parameters which denote for example the number of trays in a
distillation column or which define at which tray a feed stream enters the distillation
column.
nominal discrete parameters which represent parameters for which no reasonable
order can be defined (e.g. structural alternatives, type of utility streams used for
heating and cooling).
The evolutionary algorithm we defined is able to deal with all these parameter types, using
special variation procedures for each of them, like the geometric mutation for integer
parameters developed by Rudolph 1994. The parameters we applied in the optimization of the
HDA process are described in a following section.
How to define an optimization problem ?
The initialization of an optimization problem which comprises a complete chemical plant is
too tedious to formulate in a data file by hand. Furthermore, it is troublesome to check the
consistency, because of the quantity of data. ASPEN WORKBENCH is a comfortable
graphical interface,
which provides a better overview about the process structure. In the developed program
package, the ASPEN WORKBENCH is utilized for the problem definition.
The advantage beside the more comfortable initialization capabilities are the check functions
of the ASPEN WORKBENCH in every template, supported by a built-in expert system.
Furthermore existing plant flowsheets defined for the simulator can be extended with some
competing process alternatives by just adding them.
The ASPEN WORKBENCH can be used as usual with some specialties.
Splitters and mixers representing binary decision variables have to get a special marker. Every
continuous parameter to be treated by the evolutionary optimization has to be defined as if it
would be treated by an internal optimization algorithm. Limits of the variable range need to
be specified, and a marker identifies which algorithm is responsible for its optimization. This
allows the user a comfortable definition of the optimization parameters and the integration of
the built-in optimization algorithms to perform hybrid optimizations.
Finally, the target function has to be defined. The user can choose one of the four cost
calculation levels: equipment costs, capital costs, operating and capital costs and profit
calculations. The cost calculations are performed with ASPEN PLUS 9.3 costing modules.
For the optimization of the HDA Process besides the binary decision variables, which are
illustrated in Figure 5, the following set of process parameters has been specified:
&217,18286 3$5$0(7(56
p STABILIZER [1.0,5.0] 6WDELOL]HU &ROXPQ3UHVVXUH
p LIQUSEP ,1 , p LIQSEP , 4 [5.0,12 .0] 'LVWLOODWL RQ&ROXPQV3UHVVXUH
x CONVERSION , x SPLIT [0.0,1.0] 5HDFWRU &RQYHUVLRQ 6SOLWUDWLR *DV 5HF\FOH
,17(*(5 3$5$0(7(56
n FEED ,1 n FEED , 4 [3,50] IN 1XPEHU RI )HHGVWDJH RI /LTXLG6HSDUDWLRQ&ROXPQV
n STAGE ,1 n STAGE , 4 [5,100 ] IN 1XPEHU RI 6WDJHVRI /LTXLG6HSDUDWLRQ&ROXPQV
',6&5(7( 3$5$0(7(56
Notice that all parameters have a bounded domain that need to be specified by the user.
Furthermore FORTRAN Blocks that check the constraints of the process, e.g. checking the
reactor temperature, have to be introduced by the user. The constraint data is written in a
special file that is accessed by the external optimizer. For a more detailed description of the
software the reader is referred to Gross 1999 and Gross, Roosen 1998.
Results:
The optimization has been carried out with a (15,5,100)-ES, applying a random initialization
for the flowsheet variables in their domains. This is a kind of standard ES for difficult target
functions. In order to avoid unnecessary fitness evaluations, the mutation operator and
recombination operator generate parameters that are in their feasible domains, by applying a
transformation function.
0(7+$1(
b1
b2 { 0 ,1}
b3
0,;(5
5(&<&/(
%(1=(1(
',3+(1</
Figure 5: A partial view of the optimized flowsheet of the HDA Plant : Only the lower part of the entire
flowsheet is described here, because the configuration of the other part is like that one in Figure 4. The black
lines represent the active parts of the flowsheet. The small colored boxes are the splitters to which the binary
decision variables b1 , b2 and b3 are assigned. Their values determine the flowsheet configuration. The small
black boxes symbolize mixers. Notice that besides the three binary decision variables 18 further variables have
been adjusted by the optimizer.
As a target function for this example problem the total costs per mole product was chosen,
including investment and operating costs. The results of the optimization are shown in Figure
5. The found sequence for the liquid separation is in agreement with other documented
solutions. The flash instead of the stabilizer column is found to be the better solution for the
removal of the light ends. This takes care of the optimal trade-off between product losses
and the investment costs for the flash. For more details please see Gro 1999 or Gro, Roosen
1998.
The given results were obtained on a SPARC 20 processor each optimization run taking
between 1 and 7 days, which is quite a long time caused by the high accuracy simulation
processing. When comparing to alternative approaches the easy formulation of the
optimization task, saving significant preparation time, must also be taken into consideration.
As stated before no requirements of steadiness or derivativity have to be met.
The necessary total time needed for an optimization run can be strongly reduced by exploiting
the parallelization potential of the EA paradigm. As there is no information exchange between
the evaluations of co-existing individuals, they can be treated completely in parallel resulting
in an approximately linear speedup of this robust optimization process with increasing
numbers of processors.
parameters like the choice of utilities, the choice of a cut-sequence and the number of trays, as
well as the usual continuous parameters could be included as decision parameters. This led to
better results than calculated by the built-in optimization procedures, as can be seen in Gro,
Roosen 1998.
Furtheron we developed graph representations for the overall plant optimization and graph
operators. They allow the formulation of a combined knowledge based search improved by
evolutionary computation. The interested reader is referred to Emmerich 1999.
Future research will aim at the integration of heuristic knowledge to reduce the parameter
space to be examined. This can be achieved by excluding physically or thermodynamically
unreasonable solutions before time consuming simulation calculations are executed.
Furthermore the integration of the C++ program package into the PC environment of ASPEN
PLUS 10.0 including the integration of Windows interoperability features like ActiveX
Automation or OLE and a support of user defined or commercial costing modules shall be
enforced in the future. Another important research direction will be the development of
improved and more powerful EAs for the given applications, both by theoretical and
empirical investigations.
Acknowledgements
We gratefully acknowledge the financial support of the German Volkswagenstiftung enabling
us to perform the reported research and development.
LITERATURE:
Bck, Th. et. al. Handbook of Evolutionary Computation. Oxford Univ. Press, New York, and Inst. Of Physics
Publ., Bristol, 1997
Douglas, J.M., Conceptual Design Of Chemical Processes, McGraw-Hill, NY, 1988
Emmerich, M. Optimierung verfahrenstechnischer Prozestrukturen mit evolutionren Algorithmen, Technical
report, Dept. Computer Science, University of Dortmund, 1999
Gro, B., Roosen, P. (1996), Optimization Strategies and their Application to Heat Exchanger Network
Synthesis Chem. Eng. Tech. 2 175
Gro, B. et. al., (1996) Optimization of heat exchanger networks by means of evolution strategies. In H.-M.
Voigt, W. Ebeling, I. Rechenberg and H.P. Schwefel, Eds., Parallel Problem Solving from Nature PPSN IV,
Springer, Berlin, LNCS 1141
Gro, B. (1999) Gesamtoptimierung verfahrenstechnischer Systeme mit Evolutionren Algorithmen,
Dissertation, Fortschritts Berichte VDI, Reihe 3, Nr. 608, VDI Verlag
Gro, B., Roosen P. (1998) Total process optimization in chemical engineering with evolutionary algorithms,
Comp. Chem. Eng. Suppl. Vol 22, ESCAPE 8, pp. 229-236
Grossmann I. et. al. (1996) New Trends in optimization-based approaches to process synthesis, Comp. Chem.
Eng. Vol. 20, No. 6/7, pp. 665-683
Hesselmann, K. (1985) Wrmeaustauscher-Netzwerke: Eine exergokonomische Betrachtung. Dissertation,
RWTH-Aachen
Linnhoff, B, Flower, J.R., (1978) Pinch Analysis a state-of-the-art overview, AIChE Jl, 24 (4): 642-654
Maia, L.O.A. et. al. (1995) Comp. Chem. Engng. 19 481
Rudolph, G. , An evolutionary algorithm for integer programming, InY. Davidor, H.-P. Schwefel and R.Mnner,
Eds. Paralel Problem Solving from Nature PPSN III, Springer, Berlin, 1994, LNCS 866, 138-148
Ryoo, H.S. and Sahinidis, N.V. (1995), Comp. Chem. Engng. 19 551
Schwefel, H.-P. , Evolution and Optimum Seeking. John Wiley & Sons, New York, 1995