A genetic algorithmic framework for process design and optimization

Computers them. Engng. Vol. 15, No. 4, pp. 217-228. 1991 Printed in &eat Britain.All rightsreserved

A GENETIC

ALGORITHMIC DESIGN AND

0098-1354/91 $3.00 + 0.00

Copyright0 1991 Pergamon pTc88plc

FRAMEWORK OPTIMIZATION

FOR

PROCESS

I. P. ANDROULAICIS and V. VxNrcATAsuBRAMANL4Nt Laboratory for Intelligent Process Systems, School of Chemical Engineering, Purdue University, West Lafayette, IN 47907, U.S.A. (Received

27 November 1989;final revision received 19 November received for publication 5 December 1990)

1990;

Abstract-A general optimization framework for discrete and continuous problems based on genetic algorithmic techniques is presented. The proposed framework exhibits a structured information exchange leading to an efficient search procedure. The utility of the method is demonstrated by applying it to the design of an unsplit heat exchanger network. We also extend the basic idea of genetic algorithms to develop a new method, called the extended genetic search (EGS), for optimization in continuous spaces. The performance of the extended genetic search for several test cases of nonlinear unconstrained and constrained optimization problems is also presented. We also compare the performance of this technique with others, including simulated annealing. The proposed framework turns out to be quite robust and is able to locate global optimum solutions for problems where gradient-based algorithms fail. The results obtained were comparable to the ones obtained with simulated annealing. The implementational simplicity of the EGS framework makes it particularly useful as a tool for problems where approximate solutions are needed quickly.

1. INTRODUCTION

Computer-aided process synthesis is an important area in process systems engineering and has commanded a considerable amount of attention from researchers over the past two decades or so. Many important research milestones have been accomplished over the years, largely by the use of mathematical programming techniques (Nishida et al., 198 1). In the mathematical programming framework, many of the process design problems are formulated as mixed-integer non-linear programming (MINLP) models that have continuous variables and integer decisions. Currently available algorithms in this framework can lead to sub-optimal solutions in the presence of non-convexities. Most such approaches are deterministic in nature, resulting in convergence to the closest local optimum solution without much effort to identify the global optimum. In such cases, the choice of starting points is a very crucial factor in the performance of the algorithms. A further disadvantage is the fact that the first- and secondorder derivatives are needed in this framework. As a result the overall computational effort is increased, and furthermore, the same techniques cannot be used for both differentiable and non-differentiable problems. Some recent developments have shown that better results can be obtained using stochastic methods. An interesting work in this context is the use of simulated annealing for the design of heat exchanger networks (Dolan et al., 1987, 1989).

In this paper, we present a new optimization framework based on the concept of genetic aigorithms. Genetic algorithms are a class of optimization algorithms that seek improved performance by sampling areas of the parameter space that have a high probability for leading to good solutions. The algorithms are called genetic because the manipulation of the possible solutions resembles the mechanics of natural selection. We demonstrate the utility of this framework by applying to to the design of unsplit heat exchanger networks. We also extend the basic idea of genetic algorithms to develop a new method, called the extended genetic search (EGS), for optimization in continuous spaces. An important aspect of this method is the development of a unified framework that can deal with both continuous and discrete variables at the same time. The performance of the extended genetic search for several test cases of non-linear unconstrained and constrained optimization problems is also presented. 2.

GENETIC

ALGORITHMS

Genetic algorithms (GAS) are a general methodfor searching a solution space in a manner analogous to the natural selection procedure in biological evolution (Holland, 1975). An an optimiztechnique, it examines and manipulates ation simultaneously a set of possible solutions. Each candidate solution is represented by a string of symbols. The set of solutions Pi, is referred to as the population of the jth generation. The population

ology

TAuthor to whom all correspondence should be addressed. 217

218

I.

P. ANDROULAKI~

and V. VBWUTA.WBIUMANIAN

evolves for a prespecified total number of generations, under the application of evolutionary rules called genetic operators. The basic structure processed by GAS is the strhg. Strings are compo_sd of a sequence of characters of finite length 1 composed over some alphabet V. Strings may be represented as follows: A=ala2...aL. Strings of the current population are then manipulated to generate a new population in the next time step. This is done by the use of the transition rules, namely, reproductive plan and genetic operators. Reproductive plan determines the number of copies of an existing string to be made during a reproductive cycle. Genetic operators determine the modification and combination of these strings which will form the strings of the next generation. A simple GA may be composed of one reproductive plan and two genetic operators such as: l Fitness proportionate reproduction 0 Crossover 0 Mutation

According to “fitness proportionate reproduction”, the probability that a particular string will participate in the process of creating the new generation is set proportional to its fitness, where fitness is defined as some non-negative measure of merit. Clearly high fitness individuals will have a higher expected number of offsprings. By definition, fitness proportionate re&duction guides the search towards the most promising direction. If we consider populations with only one member, then the procedure reduces to gradient descent. Crossover is the most important operator of a genetic-based technique. First, members of the current generation are mated at random. Then, an integer position k along the string is selected at random, with uniform distribution, in the interval, [l, 1 - I]. Then two new strings are created at random, by exchanging all the characters between positions k and II inclusively. As an example consider: A = aIa2 a3a4,

B = b,b2b3b4, For k = 3, the resulting crossover will be: A’=a,a2a,b4, B’ = b,b2b3a4. Crossover results in a randomized yet structured information exchange. Each solution created combines the characteristics of both parents. In this case, crossover is like linear combination. In every search procedure, there is a trade-off between creating new knowledge and exploiting the already existing knowledge. Crossover permits us to exploit already existing knowledge. New knowledge can he introduced in the system by mutation.

Basically, it is the occasional and random alteration of a string position, and as a result it introduces a sort of a random walk in the search space. The main characteristics of this genetic-based search procedure are the following:

(9 it is a multiple point search technique that

attempts to identify “good regions” rather than “good points”; (ii) the operators make use of a coding of the parameter space; only objective function information is used. This results in a simpler implementation; (iii) although the approach has a stochastic flavor, it makes use of all the information that has been obtained during the search. It also permits the structured exchange of that information. Holland (1975) defmed the concept of schemata as the subset of the set of strings of given length consisting of those strings which have specific values at specific positions. He further proceeded to show that through a fitness proportionate reproduction rule, an exponentially increasing number of trials is allocated to the above average schemata. This led to what is called an optimal allocation of trials, implicitly showing that the selection of candidates for the reproductive step is performed optimally. In the following sections, we will present an application of this framework for a combinatorial optimization problem, namely, the heat exchanger network synthesis. In later sections of the paper, we will present a new genetic-based framework that we developed that can handle not only discrete decision problems but continuous ones as well. 3. COMBINATORIAL OPTIMIZATION USING GENETIC SEARCH: APPLICATION TO HEAT EXCHANGER NETWORK SYNTHKSIS

The general problem of HEN defined as follows:

synthesis, can be

A set of “cold” streams, initially at supply temperature T,,Sare to be heated at target temperatues T,,, . Another set of “hot streams”, at supply temperatures Th,Sare to be cooled to target temperatures T,,,. The problem is to determine the structure of a heat exchange network, which will bring all the streams to their target temperatures, while minimizing the cost of the plant. The complete definition of the problem, as well as the assumptions, constraints, objective function and data used can be found in Su and Motard (1984) and Pho and Lapidus (1973). Optimal HEN synthesis is traditionally associated with two characteristics, that of maximum energy recovery (MER) and minimum exchanger number (MEN). Linnhoff and Flower (1978) were the first to

Genetic algorithmicframework

recognize that it was possible to predict both MER and MEN. Later, Linnhoff and Hindmarsh (1983) developed the pinch design technique based on an understanding of the thermodynamics behind the process of heat exchange. Papoulias and Grossmann (1983) solved the HEN synthesis problem by using a transshipment model (MILP). Floudas ef al. (1986) extended this work by optimizing a superstructure developed for the network. Lately, Dolan et al. (1987, 1989) applied simulated annealing to the synthesis of optimal HEN successfully. In this section, we describe the application of the genetic-based search technique to the synthesis of optimal unsplit HENS. 3.1. Application of the genetic algorithmic technique Every synthesis of an optimum operational system is basically a search in two different spaces: (a) the space of the distinct structural alternatives, as defined by the topology and the nature of the interactive units; and (b) the space of the alternative designs for the operating units composing the system. At this point we propose a way for attacking the first of the two objectives, namely, the structural optimization of the system. Once the structure of the system has been specified, the actual performance needs to be determined. This usually requires the solution of a nonlinear problem. This parametric optimization, is performed using the GRG2 non-linear optimizer. Therefore, the genetic search is used as a master driver, searching for the optimum values of the discrete variables. For the application of the genetic search, we need a consistent representation of the search space in the form of a string. Since the independent variable used in our optimization is the structure of the network, we assign to each bit-value a unit, representing a connec-

219

tion between a cold and a hot stream, as shown in Fig. 1. A collection of such matches form a string as in: [{1,21,

(2321,

{1,111-

The grid representation of Flower and Linnhoff (1980) is adopted here. The alphabet that we use contains nh x n, characters, corresponding to the total number of distinct matches (nh and a= are the number of hot and cold streams, respectively). Therefore, the representation scheme is expanded considerably, in the sense that each bit contains the attribute of a particular module. Given this representation the genetic operators are applied in a manner as shown in Fig. 2. The examples presented in Fig. 2 are explained in the following manner. Reproduction is demonstrated through the first example where we have two matches between streams {1,3} and {2,3}; the actual nature of those streams is unimportant for the sake of our argument. Both those matches are propagated to the children solution. In the mutation example, the parent solution contained matches between streams { 1,3} and {2,3). Mutation, in our context is defined as the random alteration of a particular match. It is defined as the change of either one of the two streams specifying a particular match. In that sense, we first select a match to be mutated (namely, match { 1,3}) and then we assign to it a new value for one of its constituents streams (namely, stream 3 was changed to stream 2). Needless to mention that at this point we might have to impose constraints as far as the allowed matches are concerned, i.e. we might have agreed to consider only matches between hot-cold streams. In the last example, where crossover is presented, we assume that the two parent solutions contain, respectively, the following set of matches: [Cl,3}, {2,4}, {L4)1 and [{1,4), {1,3), (2,411. Again the nature of the streams is unimportant as long as it

I. P. ANLBROIJLAKIS and V. %NKATASUBRAMANlAN

220

Genetic

Operators I

Fig. 2. Application of geneticoperators. satisfies the possible restrictions imposed. The crossover point selected was position one and the two structures that are on the left of the crossover point (i.e. the shaded ones) are interchanged to produce the children solutions defined as: [( 1,4}, {2,4}, { 1,4}] and [{ 1,3j, (1,3}, (2,411. The main idea of our implementation is that the actual structure of the system is used as the independent variable and not some transformation of it through the use of decision (O-l) variables.

14t&ooo

h _\

In a different context, structures such as {1,3}, {1,4} etc. could have different meanings. A critical point in this approach is how to determine the string length A. This is addressed in Section 3.2. 3.2, Computational results All results presented in the following sections concern the lOsp1 HEN synthesis problem, for which the optimum solution is known. From the definition of

140.000

rat. = 1.0 rat* - 0.0

Mutation

Cr0ss0v.r .A”-‘\

,*,*,.

Ll,. optmtum

O.F.

Mutation

* j

value .

rate

= 0.0

Crossover rat* I

\

locx000

\.
O.F.

- 1.0

value

i;

5

--

80,000 t

q

OS. VQLUC

0 No.

of

Qenerotions

Fig. 3a. Cost evolution for pure mutation.

40

20 No.

of

60

80

generations

Fig. 3b. Cost evolution for pure crossover.

100

Genetic algorithmicframework

221

Llnwr

dacrogse mutation

No. of generations

No. of generotions

Fig. 4. Cost evolution for intermediaterate values. crossover and mutation, it is clear that they play a very important role in the search. Figures 3a and 3b illustrate this. When only mutation is applied to the system the search is equivalent to a random one. Gn the other hand when only crossover is applied, the system very quickly identifies tha sub-optimal solution that exists within the initial population, thus resulting in a premature convergence. Here no new knowledge is introduced. The trade-off between a crossover and mutation can be seen in Fig. 4, where the behavior has improved but still the convergence to a sub-optimal solution occurred. To avoid the problem of premature convergence, we introduced a generation-dependent mutation rate. The search starts with a high mutation rate and a low crossover rate, so that a random search is performed leading to the location of regions of good performance. We then increase the crossover rate and decrease the mutation rate in order to investigate these regions in depth. The mutation schedules tested are shown in Fig. 5. Simulation with both the mutation schedules gave the same optimum, which is also the global optimum structure found in the literature (Papoulias and Grossmann, 1983). Thus the generation-dependent mutation-crossover rates technique combines the benefits of mutation and crossover at appropriate points of the search leading to an optimal

0

eo

100

150

Of the

rate

200

No. of generations

Fig. 5. Linear and sigmoid mutation rate schedules.

Fig. 6a. Cost evolution for linear schedule.

“Slgmold”

k

dgcmasa

mutation

.

120+000

of thm

rota

2; = .a

1OQ000

1

P ~SoKlo -

1.. e

I l

l

l l

svl

.

,k\4___,__*_q 0

100

200

No. of generations Fig. 6b. Cost evolution for sigmoid schedule.

solution. The results are shown in Figs 6a and 6b. The results presented so far were based on the assumption that the number of units needed, expressed as the string length, was known. But in a synthesis problem this is often not the case. As a result the string length could also be introduced as a variable in the genetic search. In the genetic algorithmic literature, this has not been considered to any great extent. Although such messy genetic algorithms (Goldberg, 1989) have been developed to some extent, it would be very difficult for such an implementation to give satisfactory results since the dimensionality of the space is constantly changing. The reason for that is that we are not simply looking, even for simple (0 1) genetic algorithms, for the best distribution of the OS and the 1s within a string, but we are also looking for the optimum size of the string. Clearly this implies an optimization procedure taking place in a space whose dimension is constantly changing. Unless the fitness value is a monotonic function of the string size, which would guarantee that strings with rZ= n are better than strings with 1 = n + 1 for example, no conclusions can be drawn when two strings of different sixes are compared. In our computational experience on this matter, we found that the algorithm fails very quickly due to the trap of premature convergence, where a super-individual

I. P. ANDROULAJUS and V. VENKATASUBRAMANUN

222

Table 1. Effect of string length and population size

CaseNo.

A

P

Schedule

optimum (S yr)_

1

10

s 4 5 6 7 8

;: 20 10

10 IO 10 10 20 20 20 20

Linear Sigmoid Linear Sigmoid Linear Sigmoid Linear Sigmoid

49,392 44.400 45.918 43,979 44,278 44,929 44,556 44,627

; 20

’

dominates the rest of the population. In our implementation, crossover was designed such that different positions on the parent strings within their length are selected. This way we could create legltimate strings of different sizes. The algorithm failed to converge due to the fact that the dimensionality of the problem was not tixed and because no monotonicity existed between the performance of the solution and the size of the string (e.g. we can very easily obtain a network with 3 units that is much better than a network with 7 or 8 units). We then tried the idea of defining some sort of a super-string as an upper bound on the number of units that can be used. In Table 1, we present the results based on this idea for the initial populations that had members having 10 and 20 units. One other free parameter of the genetic search is the population size. Table 1 again presents typical results. In all cases the genetic search was able to locate near-optimal solutions. For the GA framework described above for discrete optimization problems, the following points can be made:

(9 The nature of GAS makes formal theoretical

assertions about their convergence properties very difficult. Therefore, it appears that the computational properties can only be established through empirical studies. A set of sufficient, but not necessary, conditions for the convergence of GAS was derived by Bethke (1980). Our empirical investigations show that the idea of incorporating a generationdependent mutation rate increases the efficiency of the search. (ii) The proposed method seems promising for the synthesis of structures that can be represented as a sequence of modules. If there are no restrictions in the way these modules can he aligned in order to perform the task, the overall approach is domain independent. If this is not the case, then we need to redefine the operators so that the resulting structures will always satisfy the constraints (an instructive example is the traveling salesman problem). In such a case, it is clear that domain-specific knowledge is needed as with any stochastic or heuristic optimization technique (Dolan ef al., 1989). (iii) An interesting aspect of the algorithm is the fact that it uses populations of solutions and therefore the final answer forms a family of oossible alternatives. If one allows for vew,

(iv)

long runs, these members will become very close to each other. Having populations of alternative solutions to a problem can become important in cases where we are not simply interested in the most economical design, but we are also interested in satisfying other objectives, such as safety, controllability etc. It has been shown previously (Huang and Fan, 1988) that incorporation of controllability to the design can modify the final solutions. As a consequence it might be helpful to have a collection of solutions to choose. from. It is also interesting to study the quality of solutions that belong to the final population and establish possible relations between them. The two reported optimum solutions for the 1Ospl HEN synthesis problem are those by Flower and Linnhoff (1980) and by Papoulias and Grossmann (1983). In both these solutions almost 50% of the total heat that is exchanged between the streams is exchanged between hot stream 3 and cold stream 5 and between hot stream 5 and cold stream 4. In a sample run of our genetic algorithms having a population of 10 members and run for a total of 10 generations, the optimum that we found was a structure having a total cost of $45,604 yr-’ achieving an energy recovery of 6115 kW (compared to the optimal value of $43,934 yr-’ and 6130 kW). Interestingly, the match between hot stream 5 and cold stream 4 appeared in all 10 members of the final population, the match between hot stream 3 and cold stream 5 appeared in 8 of the 10 members. Further, the match between hot stream 4 and cold stream 1, which exchanges almost 13% of the total heat was found in 8 out of the 10 members. The above results indicate that the algorithm is able to identify very quickly important modules of the solution and propagate them through the search (note that a 10 member population for 10 generations, creates and tests at most 100 networks). The nature of the overall algorithm, as far as the state transitions used and the selection procedure are concerned, leads to the identification of sub-optimal solutions in some cases. This is a typical behavior of any stochastic search algorithm, and therefore results can he appreciated only through the form of collective statistics. Nevertheless, as mentioned earlier, the important modules that constitute the optimal solution were retrieved very successfully.

4. CONTINUOUSOPTIMIZATIONUSING THE GENETIC

SEARCH

EXTENDED

FRAMEWORK

4.1. Non-linear unconstrained optimization The problem considered here is that of the search for a minimum of a real function f of n variables

Genetic algorithmicframework x, each of which can take any value within Xl,..., Dc R”. It is thus a question of finding a point X* of D such that: vx E o:f(x*)


(1)

is a global minimum off in D. While necessary and sufficient conditions that guarantee local optimality exist [vf(x *) = 0 and V*f(x *) b 01, there is no definite way to identify the global optimality of a general non-linear function except for very specific cases (e.g. convex problems). Deterministic methods attempt to create sequences of successive approximations that converge to points that satisfy the local optimality criteria. A simple implementation of such an idea is the Cuuchy’s method:

dk+l)=X(*)+eVf,

L = +1.

(2)

These methods are also known as the trajectory methods as after the algorithm has been completed, a trajectory connecting the starting point to the final point can be identified. Once a point that satisfies the optimality criteria is identified, the algorithm stops. As a result, such approaches are very useful when the starting point belongs the region of attraction of the global optimum. This is due to the fact that the techniques used in mathematical programming create the sequences by strictly following the decreasing gradient of the function. As a result only moves that improve the objective function are accepted, and that leads to the identification of the closest optimum to the initial point. Stochastic methods, except for the pure random search, try to identify either convergent sequences of points whose limit satisfies some property indicating an optimum solution (e.g. simulated annealing), or a collection of points (e.g. clustering techniques) in an effort to efficiently cover the entire search space, so that all the local optima, and hopefully the global as well, will be identified. In the extended genetic search (EGS) algorithm we propose, we will treat the direction of the search as the independent variable, expressed in terms of vectors with unit Euclidean norm. Our initial population is composed of vectors created randomly with uniform distribution in the space Dc R”. The details are presented in the following sections. 4.1.1. Description of the EGS algorithm. The proposed optimization framework combines the characteristics of a genetic search through recombinations and mutations of the existing solutions and the characteristics of both trajectory and clustering methods. The overall approach is stochastic in nature, which makes theoretical analysis quite difficult, especially issues concerning the convergence. Therefore the behavior is to be determined computationally through a series of experiments. As mentioned earlier, in this approach we use the direction of the search as the independent variable. At the beginning we randomly create a starting point

223

x0 E DE R”. An initial set of “m” directions is then created randomly such that all vectors have unit Euclidean norm and are uniformly distributed. Such vectors can be created as: j=l,...,

D]=_d’,

n,

(3)

& d: where di _ iV(0, 1) and n is the dimensionality of the space. A set of “m” points are then generated as: (4)

Xj=Xo+pDj,

with p being the step size. The objective functionfis evaluated at each one of these points and a weight pJ is assigned to each point depending on how “good” that particular point was. The new point x0 is now determined as:

Furthermore, this weight is also assigned to the direction vector D1 that corresponds to each of the “m” points. The weight of each vector will be used as a probability of selecting that vector during the process of generating the new generation of search vectors. In the next step, the vectors are selected with the aforementioned probability of selection and paired. The convex combinations of the old pairs create the new search directions as: &=vDk+(l D,=(I

-v)D,, -v)D,+vD,

(64 (6b)

where 0 < v 6 1. The convex combinations will play the role of a generalized crossover. The generalized mutation is the random alteration of none, one, or more coordinates of the new direction vectors. After the direction vectors have been created they are normalized. The updated center point and the new directions are used through equation (4) to create a new set of points. The procedure is repeated for a prespecified number of times. Thus, this approach combines characteristics of both the trajectory and the grid (clustering) methods. At every major iteration a collection of points is used to determine a new center of mass. At the end of the algorithm, the trajectory of these centers will determine the quality of the performed search. 4.1.2. Parameter selection. The number of vectors needed is expressed in terms of the dimensionality of the problem, based on the following heuristic. A 2-D space can be covered by use fo 2’ = 8 vectors. A general n-D space can be covered by (2% - n) vectors. This number gives an approximate estimate for an upper bound. We observed that equally good results could be obtained with a smaller number of vectors, especially for the cases where the number of independent variables was greater than 4.

I.P. ANDROULAKISand V. VENKATASU~IRAMANIAN

224

4.1.3. Computational results. We applied the proposed search procedure to a set of six test functions, introduced by Dixon and Scerzo (1978). First, we present a brief summary of the various functions used.

3

2

l

stepsize

Goldstein

and Price:

J-(x, y) = [l + (x +y + 1);(19 - 14x +3x2 - 14~ + 6xy + 3~3 x

1

x [(30 + 2x - ~JJ)~(18 - 32x + 12x2 + 48~ - 36xy + 27~31. .a

Fig. 7. Generation-dependent step size. A key parameter in the proposed search is the selection of the step size used. Very small, as well as very large step sizes, result in equally bad searches. We used a time dependent step size, similar to the one used by Vanderbilt and Louie (1984) as shown in Fig. 7. The maximum value of the step size is proportional to the volume of the search space in order to ensure an effective search of the entire search space. This is an assumption valid only for the cases where the search domain is at least bounded. For more practical applications this is a valid assumption. Both these selections are based on the experience gained through several experiments and are found to behave adequately. The EGS algorithm can be summarized as: Step 0 Create x0, Create D,, i=l,..., Step 1 Evaluate x, = x0 + @,. Step 2 Evaluate _&).

n.

Step 3 Evaluate p, = n.

tx,fw Step 4 Evaluate x0 = $r p, xi. Step 5 Apply operators to create new Di, Go to Step 1 As mentioned in the previous section, the procedure will terminate after a prespecified number of iterations. Table 2. Comparison SRSA’ GP

-

of algorithm etTidencics OSA’ 1.00

EC&

i.00

0.78 1.00 Hartman 3 ::z 1.00 0.72 Hartman 6 0.54 0.58 0.61 Shekel 5 0.64 0.47 0.80 Sbckel 7 0.81 Shekel 10 0.47 0.76 ‘Self-regulating simulated annealing (Vanderbilt and Louie. 1983). ‘Gmeralized simulated annealing (Brooks and Vcrdini. 1988). ‘Extended genetic search (our method). -Not reported.

In the region of interest (- 2 Q x, y & 2) there exist four minima. The global minimum is located at (0, -1) and/=3. a Hartman’s

f(x) =

family:

1,-ia,,tx, 1

5 ci exp

i-

1

-i:

-pij12 .

Two members were tested, having 3 and 6 variables, and n = 4. For both these cases there exist 4 local optima, one of which is the global. The domain of interest is 0 Q xi < 1. o Shekel’s family:

f(x)= +*

(Xj -

a,,)$, -

Ufj) +

ci.

It is a four-variable problem, and m = 5, 7, 10. The local optima are 5, 7 and 10, respectively. The region of interest is 0 < xj < 10. More details about the functions are given in Dixon and Scerzo (1978). The probability of obtaining the global optimum using our approach is shown in Table 2. Our results are compared with two different versions of simulated annealing reported in the literature (Vanderbilt and Lottie, 1983; Brooks and Verdini, 1988). Based on our results the following observations can be made: (i) The EGS algorithm is quite efficient in identifying the global optimum solutions. The results of Table 2 represent the number of times the algorithm located the global optimum out of 100 independent runs. The results are comparable to the ones obtained using simulated annealing. The number of function evaluations needed in our approach is up to 20% more than the number of function evaluations nsecl in self-regulated simulated annealing. For generalized simulated annealing, the authors do not provide the number of function evaluations that were needed in order to adjust the values of the parameters. In our approach, the only adjustable parameter is the maximum step size, which is set proportional to the volume of the search space as explained earlier. The choice of the functional form of the parameters used is not the optimal one, as it cannot be so

Genetic algorithmicframework

optimum solution. As a representative example of the behavior of our genetic search, we studied the 2-D Shekel function. For this function, in the region of interest there exist 6 local minima, as shown in Figs 8 and 9 (one of the six is unimportant and not shown). For this function, the genetic search was run 100 times. For 70 times it converged to the global minimum (f= -11.0122), for 18 times the second best (f = -5.38392) and for 12 times to the third best (f-5.22521). None of the runs terminated in any of the other two minima. The above results indicate that the EGS algorithm, even when it failed, was able to locate sub-optimal points of good quality. For a deterministic algorithm, the results are proportional to the magnitude of the region of attraction of each one of the optimal points and as a result there is a uniform distribution of the optima obtained when multiple local deterministic searches are preformed (a quasiNewton method with finite difference gradient, implemented in IMSL, was used for the comparison).

Fig. 8. Two-dimensionalShekel function. determined due to the nature of the search. Their selection was based on a qualitative analysis of the characteristics of those parameters, aa was explained earlier. (ii) As mentioned above, the proposed algorithm (and the simulated annealing procedures for comparison) was run 100 times and the results were evahiated. It is important to determine what kind of solutions were found during the runs that were not able to locate the global

For a stochastic algorithm where the step size is an adjustable parameter the shape of the function is very critical. The most difficult functions to be optimized are those of the Shekel’s family and this result has

10.0

-5.38392 -3.14032 (3.00822.6.99727) 1.s

so -11.0122 ~4.00262.4.002141

2.8

3.3

225

5.0

6.7

Fig. 9. Minima of the 2-D Shekel function.

I. P. ANnRouLAxIsand V.

226

also been reported by Vanderbilt and Louie (1983) and by Brooks and Verdini (1988). The reason for this is the narrowness of the region of attraction of the global minimum, compared with the entire domain. In other words:

VENICATASUBRAMANfAN

4

AbaRlnutoa

mW,=) < ~ m(V) ’ where m(V) represents the measure of region V. Figure 10 gives an example of functions that are easy and hard to optimize. Functions with deep and well-isolated minima +re really difficult to optimize. The results presented above concerning the application of our EGS framework in continuous optimixation demonstrate that this search is competitive with the best alternative methods available for global optimization. In the next section we move one step further, studying the application of our algorithm for a general non-linear constrained optimization problem. 4.2. Non -linear constrained optimization The problem of non-linear constrained optimixation is very important in chemical engineering process design and control since many of the problems fall into this category. Such problems can be stated as: minf(x) subject to: xoDsR”, g-(x) 2 0,

j = 1,.

h,(x)=O,

k-1

. . , J,

,...,

K,

The basis of our approach is to use a transformation method (Reklaitis et al., 1983) in order to convert the constrained problem to an unconstrained one and then apply our technique. The transformation used is of the form:

P(x, R) -f(x)

+ m (R, g, h ),

(7)

where f is the objective function, g the inequality constraints, h the equality constraints, R the penalty parameter and o a predefined function of the above. For our application the following penalty terms o were used: Equality constraints: w - Rh (x)2 Inequality constraints: w = R*, where

:-

-I

g(x)* 0

3

t?(X) -= 0, g(x) 3 09

The above transformation equally penalizes positive and negative derivations of the equality

Fig.

10.

Hard and

easy functions.

constraints, and penalizes only deviations of the inequality constraints. This kind of transformation has a common weakness, namely, the distortion of the penalty contours. In cases where this happens the unconstrained search method might fail to converge. But, this is a property of the unconstrained searches used in classical optimization. We believe that such problems should not appear in our method. Our implementation proceeds as follows: step 0

Step 1 step2 Step 3

get t=o, Set R=&, Set X =x#J. Optimize using genetic search P(x, R) =f(x) + MI&g, h). stop. Else go to Step 3. Ifr=t_, and go to Step 1. Set R,,, =R,+SR

Once R is fixed a search is performed. At the end of each step a new centroid is identified. The values of R are increased gradually. At the first stage we do not penalize the constraints very much, so that the algorithm can perform a more efficient overall search, even in the non-feasible region. 42.1. Computational dts. As in the case of unconstrained optimization, we compare our results with the existing results in the literature. We present the results for one of the standard test problems (Timmer, 1984). Results of two constrained NLP problems from Floudas et al. (1989) are also presented.

Genetic algorithmic framework The first problem is from Timmer (1984). The form of the objective function is very similar to the general family of the Shekel functions: mm=-?

1 I- 14(x

-P,Y(X

--Pi)+“,’

subject to XI + x* - 5 = 0, -3
-4
In the region of interest, this problem has two local minima. Our method identified the global optimum in all of the runs. The problem combines a linear equality constraint with two non-linear inequality constraints. The second example from Floudas et al. (1989) is a small problem of two variables with one non-linear equality constraint: min-

12x, -7x*+x;,

subject to --2x:+2-xX,=0, O
-3u,

=o,

XI + 2U, C 4, x2 + 2u2 < 4, x,<3andu2<

1,

xlv2and ule22 0. The algorithm again located the global optimum, as reported by the authors, in all runs from any randomly created initial point. Thus as the empirical studies indicate the proposed approach, although quite simple in its ideas and implementation, was able to locate the global optimum solutions of the problems that were tested. These results allow us to argue that this technique can be used at least in the initial phase of the optimization procedure in order to obtain good initial feasible solutions. Having feasible initial points is not always a trivial task, and several existing optimization techniques (e.g. generalized reduced gradient) will not converge if the initial points are in the infeasible region. Finally, the genetic search avoids several numerical difficulties that might be encountered, such as non-differentiability, matrix singularity, etc.

227 5. CONCLUSIONS

In this paper, we have presented a new solution framework for optimization problems based on the concept of GAS. It is shown that this framework is well-suited for problems which are wmbinatorial in nature and are often formulated in MINLP. As an example of this category, the heat exchanger network synthesis problem was explored, where GAS were used as a master driver for the selection of the decision (O-l) variables. We have also extended the basic idea of GAS to develop a new method, called the extended genetic search (EGS), for optimization in continuous spaces. The performance of the EGS for several test cases of non-linear unconstrained and constrained optimization problems is also presented. The proposed EGS framework turns out to bc quite robust and is able to locate global optimum solutions for problems where gradient based algorithms fail. The results obtained were comparable to the ones obtained with simulated annealing. An important aspect of this approach is that it generates a set of possible altemative solutions. Furthermore, the EGS algorithm, even when it failed to discover the global optimum, was able to locate sub-optimal solutions of good quality. The implementational simplicity of the EGS framework makes it particularly useful as a tool for problems where good approximate solutions are needed

quickly. REFERENCES

Allti-Pentini F., V. Parisi and F. Zirilli, Global optimization and stochastic differential equations. J. Optimization Theorv Avvlic. 47. I-16 (1985). Bethke A. D.: G%etic algorith&s as ‘function optimizers. Ph.D. Thesis, University of Michigan (1980). Bosworth J., N. Y. Foo and B. P. Zeigler, Comparison of genetic algorithms with conjugate gradient methods. Technical Report 003120-l-T, University of Michigan (1972). Brooks G. D. and W. A. Verdini, Computational experience with generalized simulated annealing over continuous variables. Am. J. Mathi Mgmt Sci. 8, 425-l49 (1988). Davis L. and M. Steenstrup, General Afgorithms and Simulated Annealing: An Overview. Morgan-Kaufman, Los Altos. Calif. 11987). DeJong. K. A., ‘An analysis of the behavior of a class of genetic adaptive systems. Ph.D. Thesis, University of Michigan (1975). Dixon L: C. W. and G. P. &ergo, Towards Global Optimization 2. North-Holland, Amsterdam (I 978). Dolan W. B.. P. T. Cummings and M. D. LeVan, Heat exchange network de&n-bv simulated annealina. Foundations Gf Computer-Ai&d besign and Opera&& Reklaitis and H. D. (G. V. EdS). Sprints, CAChE/Elsevier. New York (1987). _ Dolan W. B., P. T. Cummi ngs and M. D. LeVan, Process optimization via simulated annealing. AIChE JI 35, 725-736 (1989). Floudas C. A., A. Aggarwal and A. R. Chic, Global optimum search for nonconvex NLP and MINLP problems. Cumpurers &em. Engng 13, 1117-l 132 (1989). Flower J. R. and B. Linnhoff, A thermodynamic-combinatorial approach to the design of optimum heat exchanger networks. AIChE JI 26, 1-9 (1980).

228

I. P. ANDROULAK~Sand V.

Floudas C. A., A. R. Ciric and I. E. Grossmann. Automatic synthesis of optimum heat exchanger configurations. AZChE JI 35, 725-736 (1986). Foo N. Y. and J. L. Bosworth, Algebraic, geometric and stochastic aspects of genetic operators. Technical Report 003120-2-T, University of Michigan (1972.) Goldberg G. E., Computer-aided gas pipeline operation using genetic algorithms and rule learning. Ph.D. Thesis, University of Michigan (1983). Goldberg G. E., Genetic Algorithms in Search, Optimization and Machine Learning. Addison-Wesley, Reading, Mass. (1989). Goldberg G. E., Messy genetic algorithms: motivation, analysis, and first results. Technical Report 89003, University of Alabama (1989). Himmelblau D.M.. Applied Nonlinear Programming. McGraw-Hill, New York (1972). Holland J. H., Adaptation in Natural and Arti$ciaf Systems. University of Michigan Press, Ann Arbor (1975). Huang Y. L. and L. T. Fan, Strategy for enhancing structural controllability in heat exchanger network synthesis. Presented at the AIChE National Meeting, Washington, D.C. (1988). Kirkpatrick S., D. Gellat and M. Vecci, Optimization by simulated annealing. Science 220, 671-680 (1983). Linnhoff B. and J. R. Flower, Synthesis of heat exchanger networks-I. Systematic generation of energy optimal networks. AiChE JI 24, 633-642 (1978). Linnhoff 8. and E. Hindmarsh, The pinch design method for heat exchange networks. Chew. Engng Sci. 38, 745-763 (1983). Luenberger D. G.. Linear and Nonlinear Programming. Addison-Welsey, Reading, Mass. (1984).

VENKATASUBRAMANIAN Minoux M., Mathematical Programming: Theory and Applications. Wiley, New York (1986). Nishida M., G. Stephanopoulos and A. W. Westerberg, A review of process synthesis. AfChE 31 27, 321-351 (1981). Papoulias S. A. and I. E. Grossmann, A structural optimiiation approach in process synthesis IT. Computers c/rem. Engng 7, 707-721 (1983). Pho T. K. and L. Lapidus, Topics in computer-aided design: Part II. Synthesis of optimal heat exchanger networks by tree searching algorithms. AIChE JC 19, 1182-1189 (1973). Ponton J. W. and R. A. Donaldson, A fast method for the synthesis of optimal heat exchanger networks. Chem. Engng Sci. J. 29, 375-385 (1974). Reklaitis G. V., A. Ravindran and K. M. Ragsdell, Engineering Optimization MethodF and Applications. Wiley, New York (1983). Rodrigo F. R. B. and J. S. Seader, Synthesis of separation sequences by ordered branch search. AIChE JI 21, 885-894 (1975). Stephanopoulos G. and A. W. Westerberg, Studies in process synthesis. Chem. Engng Sri. 31, 195-205 (1976). Su J.-L. and R. L. Motard, Evolutionary synthesis of heat exchanger networks. Computers them. Engng 8, 67-86 (1984). Timmer G. T., Global Optimization: A Stochastic Approach. Centrum Voor Wiskunde Informatica, Amsterdam (1984). Vanderbilt D. and S. G. Lottie, A Monte-Carlo simulated annealing approach to optimization over continuous variables. J. Comput. Phys. 56, 259-271 (1984).