- Email: [email protected]
Optimization of finite element bidimensional models: an approach based on genetic algorithms
Finite Elements in Analysis and Design 29 (1998) 231—257
Optimization of finite element bidimensional models: an approach based on genetic algorithms W. Annicchiarico, M. Cerrolaza* Bioengineering Center, Faculty of Engineering, Central University of Venezuela, PO Box 50.361, Caracas 1050-A, Venezuela
Abstract This paper deals with the optimization of 2D finite element shapes using the very promising methods based on genetic algorithms. The codification of the design variables is carried out by generating series of strings in binary code. Classical genetic operators such as crossover, mutation and reproduction are used for the optimization process. More refined operators needed to improve the performance of the process are used as well. Some illustrative examples are presented and discussed ( 1998 Elsevier Science B.V. All rights reserved. Keywords: Shape optimization; Genetic algorithms; Finite elements
1. Introduction In the past, many design processes were made by the designer’s intuition and experience instead of an intensive application of optimization theory. This fact has recently changed since optimal designs became relevant in order to reduce costs in manufacturing. All this has been possible thanks to the development of new numerical methods such as the finite element method (FEM) boundary element method (BEM) and new optimization techniques, along with more powerful computers. The goal in optimum design is to obtain a solution for a given engineering problem which must satisfy all the limitations and constraints and under a previous selected criteria. The optimization techniques have been developed in several ways, such as the determination of the cross-sectional areas of a structure, or by considering variation in the structure geometry, i.e., both joints position and members geometry are considered as design variables. Other problems are the optimization of the restraint conditions and the shape optimization, a type of problem which has recently received much effort.
* Corresponding author. E-mail: [email protected]. 0168-874X/98/$19.00 ( 1998 Elsevier Science B.V. All rights reserved PII S 0 1 6 8 - 8 7 4 X ( 9 8 ) 0 0 0 2 2 - 5
232
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
The shape optimization problem can be solved by using several methods. One of them involves the fundamental equation written in a discrete form. Then, this equation is differentiated in order to get the function gradients. Other approaches need the differentiation of the fundamental equation and, by using the definition of an adjunct problem, the function gradients are calculated as the variations of the domain occur, as shown in Braibant and Fleury [1]. It is usually necessary to use a numerical method, which allows the user to get the values of the fundamental variables and their derivatives with respect to the design variables. Several works have been developed using the FEM to solve shape optimization problems. Zienkiewicz and Campbell [2] used such design variables (joint positions of the finite element mesh). They derived the stiffness matrix and loading vector with respect to the design variables and then they used a sequential linear programming method as the optimization technique. Ramakrishnan and Francavilla [3] used a formulation based on finite elements, but they used the penalty method as the optimization technique. Francavilla et al. [4] employed the same formulation in order to minimize the stress concentration around the concave zone of the union between two different thickness. Bhavikatti and Ramakrishnan [5,6] used a polynomial function, their coefficients being the design variables, to characterize the contour’s shape. Their goals were to minimize the stress concentration factor as well as the volume model. They obtained a uniform stress distribution along the contour shape. Again, they used a linear programming method to solve the optimization problem. Chun and Haug [7] carried out their analyses based on the gradient projection method. The goal again was to minimize the model weight, with restrictions involving the Von-Mises stresses along its moving boundary. Yoo et al. [8] and Yang et al. [9] optimized the shape of several elements in order to get the minimum weight, taking into account some stress constraints. The main drawback of these approaches to the shape optimization problems lies on the loss of accuracy when computing the function derivatives, which is an essential factor for using conventional optimization techniques. The search for a robust method, which presents a good balance between efficiency and efficacy, is of utmost concern to get good results in many different environments. In this work, genetic algorithms (GAs) have been selected due to the fact that they are theoretically and empirically proven to provide a robust search in complex space [10—12]. Evolution programs started in the 1960s, when a group of biologists [13—15] used digital computers to simulate genetic evolutionary systems. But it was Holland [10], who in 1975 published the work entitled “adaptation in natural and artificial systems”, who layed down the two main principles for current GAs: the ability of simple representation (bit strings) to encode complicated structures and the power of simple transformations to improve such structures. In his schemata theory, this author explained how, by using an appropriate control, a rapid improvement in the bit strings is obtained, as it is observed in animal populations. Their application to engineering problems is recent and it is due to Goldberg [16—18].
2. Some basic concepts in genetic algorithms Nowadays, engineers have to optimize more and more their designs, due to the high costs involved in their production. The first phase in optimal design consists in the definition of the basic
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
233
characteristics of the final product that will be optimized. These features could be the size, thickness, shape or topological configuration of their members. In genetic algorithms, a particular design (a) is represented by a set of its characteristics a , called phenotype and defined as real i numbers c A"a *a *a *2*a " < a , a 3R. (1) 1 2 3 c i i i/1 The phenotype of a structure is formed by the interaction of the total genetic package with its environment and his genotype is obtained by encoding each a into a particular code (binary in this i case). The transformation of the phenotype structure into a string of bits leads to the so-called chromosomes, and it represents, like in natural systems, the total genetic prescription for the construction and operation of some organism c c A" < a " < (e: a PM0, 1N). (2) i i i/1 i/1 Genetic algorithms operate on populations of strings (structures) and progressively (t"0, 1, 22) modifies their genotypes to get the better performance of their phenotype environment E. The adaptation process is based on the mechanics of natural selection and natural genetics. They combine the survival of the fittest together with string structures, with a structured yet randomized information, which is exchanged to form a search algorithm. In each generation, a new set of artificial creatures (strings) is generated by using bits and pieces of the fittest of the previous generation. They efficiently exploit historical information to speculate on new search points with expected improved performance. In order to use GAs we have to define an objective function, or fitness function, that measures the behavior of each individual into its environment. This function provides a direct indication of the performance of each individual to solve the optimization problem subjected to the imposed constraints of the environment. With the population ranked according to fitness, a group of chromosomes are selected from the population. There exist several methods to select parents. In this work, the following methods have been used [19]: f Stochastic sampling with replacement. f Remainder stochastic sampling without replacement. The selected chromosomes are then reproduced through the crossover and mutation operators. The crossover operator consist in taking two selected chromosomes as parents. Then, they are either crossed by using a certain probability value in order to create two new chromosomes (children) or they are directly included into the new population. It can be found many procedures to carry out this task in the technical literature [20]. The present work used both one break point (simple crossover) and two break points (double crossover), since these methods lead to simple and reliable solutions. In simple crossover, the position n along the string is randomly selected between 1 and the 1 string length less one [1, 1, !1]. Two new strings are then created by swapping all the characters between positions n #1 and 1 (see Fig. 1). 1
234
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 1. Graphic scheme of simple crossover procedure.
The two break-points method is based in the interchange of the sub-string placed between bits n and n of the parents. Both numbers are randomly chosen, in the same way as simple crossover. 1 2 The mutation operator gives each bit in a chromosome the opportunity of changing (if it is 1 to 0 or vice versa). The selection according to the fitness, combined with the crossover, provides genetic algorithms the bulk of their processing power. Mutation plays a secondary role in the operation of the genetic algorithms and it is needed because, eventhough selection and crossover effectively search and recombine extant notions, occasionally they may become overzealous and they can lose some potentially useful genetic material (1’s or 0’s at particular locations) [11,20]. Thus, the mutation operator protects against such irrecoverable premature loss of important notions. Due to the secondary importance of this operator, a low mutation probability value is usually considered. What is the power involved behind these simple transformation over a random population of n strings, that allow genetics algorithms to find the optimum point (or nearly the optimum point) in complex and non-linear situations? The answer of this question was found by John Holland and it is exposed in the ‘Schema Theorem’ or ‘The Fundamental Theorem of Genetic Algorithms’ [21,22]. A schema is a similarity template describing a subset of string displaying similarities at certain string positions. It is formed by the ternary alphabet M0.1,*N, where * is simply a notational symbol, that allows the description of all possible similarities among strings of a particular length and alphabet. In general, there are 21 different strings or chromosomes of length 1, but schemata display an order of 31. A particular string of length 1 inside a population of n individuals into one of the 21 schemata can be obtained from this string. Thus, in the entire population the number of schemata present in each generation is somewhere between 21 and n21, depending upon the population diversity. But, how many of they are processed in a useful way?. Clearly, not all of them have the same possibility to carry on its genetic information through each generation, since genetic operators could destroy some of them. Holland [10] estimated that in a population of n chromosomes, the GAs process O(n3) schemata into each generation. He gives this a special name: ‘Implicit Parallel Process’, which can be observed in Fig. 2. Eventhough at each generation we perform a proportional computation to the size of the population n, the algorithm carries out a more efficient process on the n3 schemata with no further information than that contained into the population itself. Holland has shown that ‘good schemata’, which are related to high-fitness values, are propagated generation by generation,
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
235
Fig. 2. Real and implicit parallel processes.
whereas ‘bad schemata tend to disappear. Moreover, the best schemata will reproduce exponentially. At the end of the generation process, the chromosomes in the population will be only constituted from good schemata, which he called ‘Building Blocks’.
3. Some refined operators Several new operators exists [11,20], that can be used together with the basic reproduction, crossover and mutation operators, explained before. In the present work, it has been used the following refined operators: rebirth, elitism and sharing functions to induce niche exploitation and speciation throughout the population, when dealing with multicriteria optimization [11]. 3.1. Rebirth [23] Rebirth is not really a genetic operator, because it does not apply directly over the chromosome. This operator works on the procedure to obtain the best chromosomes. The basic idea of this operator is to obtain a good fit, in order to get the desired optimum. The Genetic algorithm method is a stochastic and dynamic search process, which achieves a near optimal solution of the objective scalar function in every evaluation. This near-optimal solution means: DO !O D i bj )e , 0)e @1, k"1, 2,2, n , k k % O i
(3)
236
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
where O is the ideal objective scalar function, O is the best objective scalar function in run j, n is I bj % the number of runs. Also, in the phenotype structure DI !b D k jk )e , 0)e @1, k"1, 2,2, n , (4) k k # I k where I is the phenotype k of the ideal chromosome I. b is the phenotype k of the best k jk chromosome in the process j, n is the number of chromosomes. # When the genetic algorithm converges to the near-optimal solution, the process is stopped and it is made the ‘rebirth’ of a new population, thus creating random chromosomes in a subspace of the initial phenotype space. This subspace is defined by taking a range of variation (r) for each phenotype, which is less than the one used during the initial step, and by taking the best last value obtained before rebirth as the starting of the interval for each new chromosome. So, if P is the best last chromosome, then the new phenotype structures are bound in the following way: P !r)a )P #r, i"1, 2,2, n . (5) i i i # In this way, it is possible to fit the searching of the optimum, because a new process is started with a new random population created into a reduced phenotypical space of the best chromosomes belonging to the initial step. 3.2. Elitism The operations of crossover and mutation can affect the best member of the population, without producing offspring in the next generation. The elitist strategy [20] fixes this potential source of losing by copying the best member of each generation into the succeeding generation. The elitist strategy may increase the domination speed exerted by a super individual on the population. However, it seems to improve the genetic algorithm performance [20]. 3.3. Niche and speciation Sometimes, when several criteria are used simultaneously there exist no possibility to combine them by using a single number. When dealing with this situation, the problem is said to be a multiobjective or multicriteria optimization problem. This sort of problems do not have a unique optimal solution in the strict sense of the word, since the algorithm provides the best way to combine the criteria in order to obtain a plausible solution, but it is not the only one. The final selection depends on the designer’s criterion. When using genetic algorithms, the way to deal with multicriteria problems is by employing niche-like and species-like behavior in the members of the population [11,24,25]. In order to induce stable sub-population of strings (species) serving different sub-domains of a function (niches), it is then necessary to define a sharing function to determine the neighborhood and degree of sharing for each string in the population [25]. First of all, it is necessary to count the number of individuals belonging to the same niche. This is carried out by accumulating the total number of shares, which is determined by summing up the
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
237
sharing function values provided by all other individuals in the population, in a factor called niche count (m ) i m " + Sh[d(i, j)], (6) i j|P01 where d(i, j) is the distance between individuals i and j and Sh[d] is the sharing function. Sh[d] is a decreasing function of d(i, j), so that Sh[0]"1 and Sh[d*p ]"0. Usually, a triangular 4)!3% sharing function is used [25]: d Sh[d]"1! for d)p , 4)!3% p (7) 4)!3% Sh[d]"0 for d'p , 4)!3% here p is the niche radius, defined by the user at some estimates of the minimum separation 4)!3% desired or expected between the desired species. Then, the selection of the best no-dominating individuals to go on reproduction will be made over the one that has the less niche count m . Thus, the pareto frontier will be formed by a wide i variety of individuals.
4. Bidimensional and tridimensional models optimization The optimization process carried out by genetic algorithms can be applied to a different optimization engineering problems. In order to test its robustness, efficiency and flexibility, we will discuss herein the application of GAs to 3D-truss structures, due to their easy interpretation compared with other optimization techniques. Following, the optimization of bidimensional finite element models will be discussed. 4.1. The elastic problem In elasticity, design specifications and restrictions must be imposed in a model by defining limit values for stresses and displacements. In order to obtain those values in each design a basic finite element formulation is used to analyse 2D FEM models, as well as a classic displacement-based matrix method is used to analyze 3D truss structures. The basic formulation of the FEM method in elasticity can be found elsewhere [26,27]. 4.2. Chromosome representation In the examples presented herein, 3D truss and 2D finite-element models of mechanical pieces, the nodal restraints, loading conditions and the topological configuration are known. Hence, the design variables are f For truss structures: 1. the cross section type (t ) of each group of bars. j 2. the increment value (D ) of each joint coordinate (x ) h h
238
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
f For finite elements models: 1. the increment value (D ) of each joint coordinate (x ) in finite elements along the moving h h boundary. Thus, the structure characteristics or phenotypes can be represented as a chain of decimal parameters, as shown below: Structure"St t 2t 2t aD D 2D 2D mT, 12 j n 1 2 h n
(8)
where t is the the integer representing the cross sectional type, (D ) the real number representing j h the increment of the displacements, n the number of group of bars, n the number of mobile joint a m coordinates. The genotype structure is then obtained by encoding each parameter in binary format, in the following way: Individual"S1 0 0 0 1 0 1 1 0 0 1 120 1 1T.
(9)
4.3. The objective functions and constraints Once the population has been created, it is necessary to define the merit functions, which permit the measurement of the success or failure of the chromosomes of the problem to be solved. In the case of 3D-truss optimization, the problem is to find the minimum weight without violating the constraints M*/ N nj Min(¼)" + A + o ¸ j i i j/1 i/1
(10)
p !p61)0, i"1, 2,2, n , l"1,2, n , il "!3 #!4%4 p61"pt ° pc, i i
(11)
S !S61)0, S61"S5 ° S#, i
(12)
u !u61)0, k"1, 2,2, n , kl k r
(13)
s.t.
where A is the area of the section type j, o , ¸ : the density and length of the member i, p the stress j i i il of the member i under the load case l, pt, pc the upper bounds of tensile and compression stress of i i member i, S the slenderness ratio of bar i, S5, S# the upper slenderness bounds for tensile and i compression bars, º the elastic displacement of the k degree of freedom under load case l, u61 the kl k upper o lower limit of the elastic displacement of k d.o.f., and n the number of restricted r displacements. In the case of 2D finite element shape optimization, the problem consists in finding the minimum area of the model without violating the stresses and displacements limit values. The objective
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
239
function is then stated by N %-
Min(¼)"+ » o . i ii 1
(14)
In multiobjective optimization, the above objective function is used together with the minimization of the stress concentration factor (k) p Min(K)" .!9 , p .%$
(15)
where » , o are the volume and specific weight of element i. N the number of finite elements of the i i %model, p , p the maximum and mean Von-Mises stress of the model. subjected to .!9 .%$ (a) Stress restrictions: the Von-Mises stresses, calculated at Gauss points of the finite element, must not exceed the limit value p !$. p i!p )0. 70/ !$.
(16)
The Von-Mises stress is calculated as p "[(p !p )2#3q2 ]1@2. 70/ x y xy
(17)
(b) Nodal coordinates restrictions: nodal coordinates of some nodes of the mesh should not move beyond certain limit values in X and ½ directions in order to maintain the mesh topology. (c) Restrictions on the shape of the elements of the moving boundary, in order to avoid singular or negative Jacobians. In order to transform the constrained problem in an unconstrained problem and due to the fact that GAs do not depend on continuity and existence of the derivative, penalty methods have been used in this paper [11,18,23,31].Thus, the evaluation function is the minimization of the model weight, penalized with the constraints g "¼#jD*pD2#kD*ulD2#gD*xlD2, i
(18)
where g is the penalized weight, i the structured index, D*pD the increment over allowed stress, D*ulD the increment over allowed displacement, D*xlD the increment over allowed coordinate, (l"1,)2,)3). The parameters j, k and g are adjusted by trial and, in this paper, they have been evaluated in such a way that a 10% of increase in stress or displacement increases the structure weight by about 5% [18]. 4.4. Fitness function and normalization The duality of cost minimization and profit maximization is well-known. In normal operations, the transformation of a minimization problem into a maximization problem is simply carried out by multiplying the cost function by minus one.
240
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
In this paper, the minimization of the evaluation function has been transformed into a maximization search by using Constant fitness(i)" . g(i)
(19)
This new evaluation function is called fitness function. The relative rate between chromosomes is maintained as defined by the original evaluation function (Eq. (18)). When working with GAs, an early trend on few members to dominate the selection process is observed. therefore, the optimum solution cannot be reached. As a consequence, when the population is very convergent the average fitness is close to the best. This fact leads to a strong competition between many individuals displaying very similar fitness, thus generating a somewhat erratic behaviour of the algorithm. In both cases, it is necessary to normalize the fitness function in order to keep the competition more or less the same in generation after generation. This situation can be controlled by the normalization number f N " .!9 , (20) n f .%$ where f is the maximum fitness, f the averages fitness. Goldberg [18] suggested that these .!9 .%$ values should be taken between 1.2 and 2 for small population (n "50—100). Here, N "2 has 101 n been taken for all generations, and the Forrest transformation has been adopted in order to obtain it [28].
5. Illustrative examples 5.1. Three-dimensional truss steel tower The topology and dimensions of a three-dimensional 25-bar-truss [29,30] are shown in Fig. 3. The goal is to optimize the weight of the structure by considering the following cases: 1. Continuous optimization of the design variables such as cross section types and geometry configuration ( joint coordinates). 1.1. optimization of cross-sectional types with fixed geometry configuration, 1.2. optimization of cross-sectional types and variable geometry configuration. 2. Discrete optimization of the cross-sectional areas of the members. Table 5 displays the sections used to optimize the structure [30]. Tables 1—4 contain the design properties, loading cases and the design variables chosen for optimization with the genetic algorithm (see Table 5). Fig. 4 illustrates the best three runs of the weight evolution of the truss versus the number of generations. Those analyses have considered stress, displacement and buckling restrictions, by assuming a fixed geometry. Note that the initial weight was about 2333.77 N and the final one was 1696.45 N.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
241
Fig. 3. Three dimensional 25-bar transmission tower. Node and member numeration. (dimensions in mm). Table 1 Design properties (three-dimensional 25-bar truss tower) Young modulus Density Maximum tensile or compression stress Maximum displacement (nodes 3 and 4) Maximum buckling stress
68 950 Mpa 0.02713 N/cm3 275.8 Mpa $0.89 cm Following the AISC methodology
Table 2 Loading cases (three-dimensional 25-bar truss tower) Hypothesis
Node
F (kN) x
F (kN) y
F (kN) z
1
1 2 3 4
2.224 2.224 4.448 4.448
0.0 0.0 44.48 44.48
0.0 0.0 !22.24 !22.24
2
3 4
0.0 0.0
88.96 88.96
!22.24 !22.24
Fig. 5 also shows the weight evolution of the truss versus the number of generations, but in this case the variation of the geometry of the truss was included into the optimization criteria. It can be observed that the final weight, 974 N, reflects a reduction of more than 700 N with respect to the fixed geometry case.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
242
Table 3 Physical design variables (three dimensional 25-bar truss tower) Variable
Element
Minimum value (cm2)
1 2 3 4 5 6 7 8
1 2345 6789 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
0.645 0.645 0.645 0.645 0.645 0.645 0.645 0.645
Table 4 Geometric design variables..(three dimensional 25-bar truss tower) Variable
Coordinate (nodes)
» ($cm) .!9
» ($cm) .*/
9 10 11 12 13
Z(1,2,5,6) X(1,2,5,6) Y(1,2,5,6) X(7,8,9,10) Y(7,8,9,10)
254 95.25 95.25 254 254
25.4 25.4 25.4 25.4 25.4
Table 5 Member types used for discrete optimization (three-dimensional 25-bar truss tower) Type
Area (cm2)
Type
Area (cm2)
Type
Area (cm2)
1 2 3 4 5 6 7 8 9 10
0.645 1.290 1.935 2.580 3.225 3.870 4.516 5.161 5.806 6.451
11 12 13 14 15 16 17 18 19 20
7.096 7.741 8.387 9.032 9.677 10.322 10.967 11.612 12.258 12.903
21 22 23 24 25 26 27 28 29 30
13.548 14.193 14.838 15.483 16.129 16.774 18.064 19.356 20.645 21.935
The results of the three best runs when dealing with discrete optimization are shown in Fig. 6. In this case, the optimization was carried out including stress and displacement restrictions. In the Fig. 7, the buckling stress of the bars are compared against the maximum value allowed for them. As it can be noted, all stresses in the bars are below this limit.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
243
Fig. 4. Genetic history with fixed geometry (Three dimensional 25-bar truss tower).
Fig. 5. Genetic history with variable geometry (Three dimensional 25-bar truss tower).
Table 6 summarizes the results obtained in the different cases under study. Finally, in Fig. 8 we can observe the initial and final geometry configuration of the tower. Note that in the final ground-plan view the algorithm has reduced the distance between supports along the x-direction and how in the final lateral view, the first level of the tower has moved upwards 5.2. 2D finite element models The following examples illustrate the ability of the algorithm to optimize finite element model shapes, in order to get the minimum area and an uniform stress distribution.
244
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 6. Genetic history with discrete optimization (Three dimensional 25-bar truss tower).
Fig. 7. Comparison between buckling stress in members and the maximum allowed stress (three dimensional 25-bar truss tower).
5.2.1. Fillet bar The first example corresponds to a fillet bar subjected to uniform loading. The dimensions and geometry of the model are shown in Fig. 9. The loading and boundary conditions are displayed in Fig. 10. The goal is to find the best shape of the boundary AB which produces a Von-Mises stress less than or equal to p "620.55 KPa. In Fig. 11 we can see the final shape of the model, where it can # be noted the optimized AB contour. The comparison between the initial and final stresses is shown in Figs. 12 and 13. Note how the final stress is very close to the allowed value.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
245
Table 6 Final results of the different runs (three-dimensional 25-bar truss tower) Variable active restriction
1 2 3 4 5 6 7 8
Fixed geometry (stress/displ/buckling) (cm2) 3.193 8.909 15.987 1.742 0.935 6.561 4.219 4.361
Final Weight (N) 1698.18
Variable geometry (stress/displ/buckling) (cm2) 2.380 14.006 8.335 1.0064 1.735 2.632 4.206 3.703 973.88
Fig. 8. Final three dimensional 25-bar truss tower configuration.
Discrete optimization (stress/displ/buckling) (cm2) 0.645 10.967 5.806 0.645 1.290 0.645 0.645 0.645 787.22
246
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 9. Fillet bar geometry and dimensions.
Fig. 10. Loading and boundary conditions.
Fig. 11. Final shape of the optimized piece.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
247
Fig. 12. Initial stress distribution (stresses in KPa).
Fig. 13. Final stress distribution (stresses in KPa).
The stress concentration factor changed from K "1.830 in initial configuration to K "1.5635 * & in the final one. This fact can be clearly observed in the final stress distribution. Fig. 14 shows the evolution of the maximum stress (p ) divided by the allowed stress (p ) and .!9 # the evolution of the weight of the model divided by the initial weight, versus the number of
248
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 14. Genetic evolution of stress and weight.
Fig. 15. Left window: smooth shape. Right window: final smooth stress.
generations. Those quantities factors are below one, which is an indication of the quite confidence of the reached solution. In Fig. 15, a smooth shape of the optimized piece is shown. It was obtained modifying the final shape shown in Fig. 11. The aim is to show how the designer can act interactively with the program
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
249
and modify the contour given by the optimization process. Moreover, it can be seen how the final stress in the smoothed piece are similar to those obtained in the optimized model (Fig. 11). 5.2.2. Cross section of a pipe Fig. 16 shows a quarter part of a cross section of a pipe, which is loading with an internal constant pressure of P "9810 KPa. The goal is to find the best cross section with the minimum 1 weight, by providing that the following constraints are not violated: 1. Inner radius. 2. Thickness of the pipe’s wall. 3. Maximum allowed Von-Mises stress. In the Fig. 17, the final mesh of the model is displayed. It corresponds to the best design obtained by the algorithm. It can be observed how the optimization process has smoothed the piece’s contour, thus getting a circular shape, as expected. Figs. 18 and 19 display a comparison between the initial stresses and the final ones. It can be noted an important reduction in the stress values. The maximum stress moved down from 2374.02 to 1265.50 MPa, which represents almost a reduction of 50% in the maximum stress. Fig. 20 shows the weight evolution of the model, versus the number of generations. The evolution of the maximum stress (p ) divided by the allowed stress (p ) is depicted in .!9 # Fig. 21. Note how this factor is lower than one. 5.2.3. Plate subjected to tractions As a demonstration of effectiveness of shape optimization using GAs, the optimization of the surface area of a plate is considered. The Fig. 22 contains the geometry, dimensions and the moving boundary AB of the plate. The loading and boundary conditions are depicted in Fig. 23
Fig. 16. Initial geometry, loading and boundary conditions (cross sectional pipe problem).
250
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 17. Top window: Zoom on the final geometry reached by the algorithm; bottom left window: Initial geometry of the model; bottom right window: final geomety of the model.
Fig. 18. Initial Von Mises stress distribution (stresses in MPa).
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
251
Fig. 19. Final Von Mises stress distribution (stresses in MPa).
Fig. 20. Genetic history of cross sectional pipe problem.
The objective function of the problem is the optimization of the surface area of the plate. The constraints are the Von-Mises stresses, which must not exceed 34 475 kPa. The initial area of the shape was 619.35 cm2, as shown in Fig. 22. The optimized shape is shown in Fig. 24. It can be noted how the plate area was reduced significantly. The final area obtained was 456.65 cm2, which means a 26.27% reduction. The comparison between the initial and final stresses is shown in Figs. 25 and 26. Note that the final shape is the result of the interaction between allowed stress and the elements stresses in the narrow zone of the final shape.
252
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 21. Evolution of the tension factor (p /p ) versus number of pipe problem). .!9 #
Fig. 22. Plate with moving boundary (A-B), geometry and dimensions.
Finally, the Fig. 27 illustrates the evolution curves of the structure stresses and weight. It can be noted the slope decreasing in the adimensional weight factor as the generation process progresses. This factor is always below than one, which clearly shows that the algorithm is stable and will converge to the optimal solution in terms of minimum weight. Another subject dealing with the relation between both curves should be remarked herein. When the adimensional stress curve displays its minimum value (see Fig. 27, generation number 17), it is a clear indication that the area
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
253
Fig. 23. Loads and boundary conditions.
Fig. 24. Final shape of the optimized piece.
of the piece can be further reduced. However, beyond generation number 40 or so, the adimensional stress curve displays a nearly constant value close to one, which involves that the area is very close to its optimal value. This explain why the adimensional weight curve behaves in very stable manner after generation 40.
254
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 25. Initial stress distribution (plate subjected to traction).
Fig. 26. Final stress distribution (plate subjected to traction).
6. Conclusions The deterministic classic resolution in optimization methods starts from a predefined point and, from it the search is performed until the nearest minimum value is achieved. However, there is the possibility that such a minimum value is not the required point. When using genetic algorithms, the search starts from a set of points which form the initial population and from them the exploration
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
255
Fig. 27. Genetic evolution of stresses and weight (plate subjected to traction).
of the design space of plausible solutions is started. It is carried out by the random application of crossover and mutation operators. One of the greatest advantages of the method is that it is a discrete method. The power of the GA, based in the Darwin mechanics and having certain random behavior, is due to the exponential increment of the short and low degree “schemata” associated to good merit functions. Its discrete formulation and its ability in searching, suggest that this method is a powerful tool in engineering applications. New algorithms to simulate the process artificially have been developed and applied to various types of engineering optimization problems. Studying the algorithms from the optimization point of view, they surpass the other existing optimization techniques in some points, and they require no calculation of sensitivities and comparatively to attain the global optimum point. The 2D finite element models optimized herein have shown that the proposed approach is suitable to carry out this task. The results indicate a good performance even in the presence of stress concentrations caused by corners and geometry changes. A promising field in structural optimization is currently under research. The application of GA methods to finite element models suggest the exciting possibility of optimizing complex 3D problems, such as the shape optimization of gravity dams.
Acknowledgements The authors wish to acknowledge the financial support provided by Consejo de Desarrollo Cientifico y Humanistico (CDCH) and Consejo Nacional de Investigaciones Cientificas y Tecnolo´gicas (CONICIT). Also, to H. Bricen8 o for her careful revision of the english manuscript.
256
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
References [1] V. Braibant, C. Fleury, Sensitivity analysis in shape optimal design, in: De. C.A. Mota (Ed.), Computer Aided Optimal Design, vol. 2, Troı´ a Portugal, 1986, pp. 7—28. [2] O.C. Zienkiewicz, J.S. Campbell, Sequential linear programming, in: De. R.H. Gallagher, O.C. Zienkiewicz (Eds.), Optimum Structural Design, Wiley, New York, 1973, pp. 109—126. [3] C.V. Ramakrishnan, A. Francavilla, Structural shape optimization using penalty functions, J. Struct. Mech. 3 (4) (1975) 403—432. [4] A. Francavilla, C.V. Ramakrishnan, O.C. Zienkiewicz, Optimization of shape to minimize stress concentration, J. Strain Anal. 10 (2) (1975) 63—70. [5] S. Bhavikatti, C.V. Ramakrishnan, Optimum design of fillets in flat and round tension bars, ASME Paper 77-Det-45, 1977. [6] S. Bhavikatti, C.V. Ramakrishnan, Shape optimization of structural systems using finite elements and sequential linear programming, in: A. Wexley (Ed.), Large Engineering Systems, Pergamon Press, Oxford, 1977, pp. 224—225. [7] Y.W. Chun, E.J. Haug, Two dimensional shape optimal design, Int. J. Numer Meth. Eng. 13 (5) (1978) 311—336. [8] Y.M. Yoo, E.J. Haug, K.K. Choi, Shape optimal design of multiple loaded structural components, in: H. Eschenaver, N. Olhoff (Eds.), Optimization Methods in Structural Design. pp. 250—255, 1983. [9] R. Yang, K.K. Choi, E.J. Haug, Finite element computation of structural design sensitivity analysis. Report no. 3, University of Iowa, 1984. [10] J. H. Holland, Adaptation in Natural and Artificial Systems, University of Michigan Press, Michigan, 1975 and First MIT Press ed., Cambridge, Ch. 7, 1992, p. 121. [11] D.E. Goldberg, Genetic Algorithms in Search Optimization and Machine Learning, Chaps. 3 and 4, AddisonWesley, Reading, MA, p. 89, 1989. [12] L. Davis, Genetic Algorithm and Simulated Anneling, Morgan Kaufmann, Los Altos, CA, 1987. [13] N.A. Barricelli, Numerical testing of evolution theories, ACTA Biotheoretica 16 (1962) 69—126. [14] A.S. Fraser, Simulation of genetic systems, J. Theoret. Biol. 2 (1962) 329—346. [15] F.G. Martin, C.C. Cockerham, High speed selection study, in: O. Kemptthorne (Ed.), Biometrical Genetics, Pergamon Press, London, 1960, pp. 35—45. [16] D.E. Goldberg, Computer-aided gas pipeline operation using genetic algorithm and rule learning. Doctoral Dissertation, University of Michigan. Dissertation Abstracts International 44 (10), 317413 (University mycrofilms num. 8402282), 1983. [17] D.E. Goldberg, C.H. Kuo, Genetic algorithms pipeline optimization, J. Comput. Civil Eng. 12 (1987) 128—141. [18] D.E. Goldberg, M.P. Santani, Engineering optimization via genetic algorithm. Proc. 9th Conf. On Elect. Comp., ASCE, New York, 1986, 471—482. [19] A. Brindle, Genetic algorithms for function optimization, Doctoral Dissertation, University of Alberta, Edmonton, 1981. [20] L. Davis, Handbook of Genetic Algorithms, Van Nostrand Reinhold, New York, 1991. [21] J.H. Holland, Processing and processors for schemata, in: E.L. Jacks (Ed.), Associative information Processing, Elsevier, New York, 1971, pp. 127—146. [22] J.H. Holland, Schemata and intrinsically parallel adaptation. Proc. NSF Workshop of Learning System Theory and its Applications, Gainesville, University of Florida, 1973, pp. 43—46. [23] M. Galante, Genetic algorithms as an approach to optimize real-world trusses, Int. Numer. Meth. Eng. 39 (1996) 361—382. [24] J. Horn, N. Nafpliotis, Multiobjective optimization using the niched pareto genetic algorithms. Illinois genetic laboratory, Report No. 93005, July 1993. [25] D.E. Goldberg, J. Richardson, Genetic algorithms with sharing for multimodal function optimization, Proc. 2nd Int. Conf. on Genetic Algorithms, 1987. [26] O.C. Zienkiewicz, The Finite Element Method in Engineering Science, McGraw-Hill, New York, 1979. [27] S.P. Timoshenko, J.N. Goodier, Theory of Elasticity, McGraw-Hill, New York, 1951.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
257
[28] S. Forrest, A study of parallelism in the classifier system and its application to classification in KL-ONE semantic networks, Doctoral dissertation, University of Michigan, Ann Arbor, 1985, unpublished. [29] P. Marti, P. Company, M. Sanchis, Acoplamiento elementos finitos-te´cnicas de Optimizacio´n en el Sistema Disseny, II Cong. Met. Num. Ing., Espan8 a, 1987. [30] S. Rajeev, C.S. Krishnamoorthy, Discrete optimization of structures using genetic algorithms. J. Struct. Eng. 118 (5) (1992). [31] M. Cerrolaza, W. Annicchiarico, Structural optimization algorithms based on Genetic Simulation. Ed.CDCHUCV, 1996 (in Spanish) pp. 163.
Optimization of finite element bidimensional models: an approach based on genetic algorithms W. Annicchiarico, M. Cerrolaza* Bioengineering Center, Faculty of Engineering, Central University of Venezuela, PO Box 50.361, Caracas 1050-A, Venezuela
Abstract This paper deals with the optimization of 2D finite element shapes using the very promising methods based on genetic algorithms. The codification of the design variables is carried out by generating series of strings in binary code. Classical genetic operators such as crossover, mutation and reproduction are used for the optimization process. More refined operators needed to improve the performance of the process are used as well. Some illustrative examples are presented and discussed ( 1998 Elsevier Science B.V. All rights reserved. Keywords: Shape optimization; Genetic algorithms; Finite elements
1. Introduction In the past, many design processes were made by the designer’s intuition and experience instead of an intensive application of optimization theory. This fact has recently changed since optimal designs became relevant in order to reduce costs in manufacturing. All this has been possible thanks to the development of new numerical methods such as the finite element method (FEM) boundary element method (BEM) and new optimization techniques, along with more powerful computers. The goal in optimum design is to obtain a solution for a given engineering problem which must satisfy all the limitations and constraints and under a previous selected criteria. The optimization techniques have been developed in several ways, such as the determination of the cross-sectional areas of a structure, or by considering variation in the structure geometry, i.e., both joints position and members geometry are considered as design variables. Other problems are the optimization of the restraint conditions and the shape optimization, a type of problem which has recently received much effort.
* Corresponding author. E-mail: [email protected]. 0168-874X/98/$19.00 ( 1998 Elsevier Science B.V. All rights reserved PII S 0 1 6 8 - 8 7 4 X ( 9 8 ) 0 0 0 2 2 - 5
232
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
The shape optimization problem can be solved by using several methods. One of them involves the fundamental equation written in a discrete form. Then, this equation is differentiated in order to get the function gradients. Other approaches need the differentiation of the fundamental equation and, by using the definition of an adjunct problem, the function gradients are calculated as the variations of the domain occur, as shown in Braibant and Fleury [1]. It is usually necessary to use a numerical method, which allows the user to get the values of the fundamental variables and their derivatives with respect to the design variables. Several works have been developed using the FEM to solve shape optimization problems. Zienkiewicz and Campbell [2] used such design variables (joint positions of the finite element mesh). They derived the stiffness matrix and loading vector with respect to the design variables and then they used a sequential linear programming method as the optimization technique. Ramakrishnan and Francavilla [3] used a formulation based on finite elements, but they used the penalty method as the optimization technique. Francavilla et al. [4] employed the same formulation in order to minimize the stress concentration around the concave zone of the union between two different thickness. Bhavikatti and Ramakrishnan [5,6] used a polynomial function, their coefficients being the design variables, to characterize the contour’s shape. Their goals were to minimize the stress concentration factor as well as the volume model. They obtained a uniform stress distribution along the contour shape. Again, they used a linear programming method to solve the optimization problem. Chun and Haug [7] carried out their analyses based on the gradient projection method. The goal again was to minimize the model weight, with restrictions involving the Von-Mises stresses along its moving boundary. Yoo et al. [8] and Yang et al. [9] optimized the shape of several elements in order to get the minimum weight, taking into account some stress constraints. The main drawback of these approaches to the shape optimization problems lies on the loss of accuracy when computing the function derivatives, which is an essential factor for using conventional optimization techniques. The search for a robust method, which presents a good balance between efficiency and efficacy, is of utmost concern to get good results in many different environments. In this work, genetic algorithms (GAs) have been selected due to the fact that they are theoretically and empirically proven to provide a robust search in complex space [10—12]. Evolution programs started in the 1960s, when a group of biologists [13—15] used digital computers to simulate genetic evolutionary systems. But it was Holland [10], who in 1975 published the work entitled “adaptation in natural and artificial systems”, who layed down the two main principles for current GAs: the ability of simple representation (bit strings) to encode complicated structures and the power of simple transformations to improve such structures. In his schemata theory, this author explained how, by using an appropriate control, a rapid improvement in the bit strings is obtained, as it is observed in animal populations. Their application to engineering problems is recent and it is due to Goldberg [16—18].
2. Some basic concepts in genetic algorithms Nowadays, engineers have to optimize more and more their designs, due to the high costs involved in their production. The first phase in optimal design consists in the definition of the basic
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
233
characteristics of the final product that will be optimized. These features could be the size, thickness, shape or topological configuration of their members. In genetic algorithms, a particular design (a) is represented by a set of its characteristics a , called phenotype and defined as real i numbers c A"a *a *a *2*a " < a , a 3R. (1) 1 2 3 c i i i/1 The phenotype of a structure is formed by the interaction of the total genetic package with its environment and his genotype is obtained by encoding each a into a particular code (binary in this i case). The transformation of the phenotype structure into a string of bits leads to the so-called chromosomes, and it represents, like in natural systems, the total genetic prescription for the construction and operation of some organism c c A" < a " < (e: a PM0, 1N). (2) i i i/1 i/1 Genetic algorithms operate on populations of strings (structures) and progressively (t"0, 1, 22) modifies their genotypes to get the better performance of their phenotype environment E. The adaptation process is based on the mechanics of natural selection and natural genetics. They combine the survival of the fittest together with string structures, with a structured yet randomized information, which is exchanged to form a search algorithm. In each generation, a new set of artificial creatures (strings) is generated by using bits and pieces of the fittest of the previous generation. They efficiently exploit historical information to speculate on new search points with expected improved performance. In order to use GAs we have to define an objective function, or fitness function, that measures the behavior of each individual into its environment. This function provides a direct indication of the performance of each individual to solve the optimization problem subjected to the imposed constraints of the environment. With the population ranked according to fitness, a group of chromosomes are selected from the population. There exist several methods to select parents. In this work, the following methods have been used [19]: f Stochastic sampling with replacement. f Remainder stochastic sampling without replacement. The selected chromosomes are then reproduced through the crossover and mutation operators. The crossover operator consist in taking two selected chromosomes as parents. Then, they are either crossed by using a certain probability value in order to create two new chromosomes (children) or they are directly included into the new population. It can be found many procedures to carry out this task in the technical literature [20]. The present work used both one break point (simple crossover) and two break points (double crossover), since these methods lead to simple and reliable solutions. In simple crossover, the position n along the string is randomly selected between 1 and the 1 string length less one [1, 1, !1]. Two new strings are then created by swapping all the characters between positions n #1 and 1 (see Fig. 1). 1
234
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 1. Graphic scheme of simple crossover procedure.
The two break-points method is based in the interchange of the sub-string placed between bits n and n of the parents. Both numbers are randomly chosen, in the same way as simple crossover. 1 2 The mutation operator gives each bit in a chromosome the opportunity of changing (if it is 1 to 0 or vice versa). The selection according to the fitness, combined with the crossover, provides genetic algorithms the bulk of their processing power. Mutation plays a secondary role in the operation of the genetic algorithms and it is needed because, eventhough selection and crossover effectively search and recombine extant notions, occasionally they may become overzealous and they can lose some potentially useful genetic material (1’s or 0’s at particular locations) [11,20]. Thus, the mutation operator protects against such irrecoverable premature loss of important notions. Due to the secondary importance of this operator, a low mutation probability value is usually considered. What is the power involved behind these simple transformation over a random population of n strings, that allow genetics algorithms to find the optimum point (or nearly the optimum point) in complex and non-linear situations? The answer of this question was found by John Holland and it is exposed in the ‘Schema Theorem’ or ‘The Fundamental Theorem of Genetic Algorithms’ [21,22]. A schema is a similarity template describing a subset of string displaying similarities at certain string positions. It is formed by the ternary alphabet M0.1,*N, where * is simply a notational symbol, that allows the description of all possible similarities among strings of a particular length and alphabet. In general, there are 21 different strings or chromosomes of length 1, but schemata display an order of 31. A particular string of length 1 inside a population of n individuals into one of the 21 schemata can be obtained from this string. Thus, in the entire population the number of schemata present in each generation is somewhere between 21 and n21, depending upon the population diversity. But, how many of they are processed in a useful way?. Clearly, not all of them have the same possibility to carry on its genetic information through each generation, since genetic operators could destroy some of them. Holland [10] estimated that in a population of n chromosomes, the GAs process O(n3) schemata into each generation. He gives this a special name: ‘Implicit Parallel Process’, which can be observed in Fig. 2. Eventhough at each generation we perform a proportional computation to the size of the population n, the algorithm carries out a more efficient process on the n3 schemata with no further information than that contained into the population itself. Holland has shown that ‘good schemata’, which are related to high-fitness values, are propagated generation by generation,
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
235
Fig. 2. Real and implicit parallel processes.
whereas ‘bad schemata tend to disappear. Moreover, the best schemata will reproduce exponentially. At the end of the generation process, the chromosomes in the population will be only constituted from good schemata, which he called ‘Building Blocks’.
3. Some refined operators Several new operators exists [11,20], that can be used together with the basic reproduction, crossover and mutation operators, explained before. In the present work, it has been used the following refined operators: rebirth, elitism and sharing functions to induce niche exploitation and speciation throughout the population, when dealing with multicriteria optimization [11]. 3.1. Rebirth [23] Rebirth is not really a genetic operator, because it does not apply directly over the chromosome. This operator works on the procedure to obtain the best chromosomes. The basic idea of this operator is to obtain a good fit, in order to get the desired optimum. The Genetic algorithm method is a stochastic and dynamic search process, which achieves a near optimal solution of the objective scalar function in every evaluation. This near-optimal solution means: DO !O D i bj )e , 0)e @1, k"1, 2,2, n , k k % O i
(3)
236
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
where O is the ideal objective scalar function, O is the best objective scalar function in run j, n is I bj % the number of runs. Also, in the phenotype structure DI !b D k jk )e , 0)e @1, k"1, 2,2, n , (4) k k # I k where I is the phenotype k of the ideal chromosome I. b is the phenotype k of the best k jk chromosome in the process j, n is the number of chromosomes. # When the genetic algorithm converges to the near-optimal solution, the process is stopped and it is made the ‘rebirth’ of a new population, thus creating random chromosomes in a subspace of the initial phenotype space. This subspace is defined by taking a range of variation (r) for each phenotype, which is less than the one used during the initial step, and by taking the best last value obtained before rebirth as the starting of the interval for each new chromosome. So, if P is the best last chromosome, then the new phenotype structures are bound in the following way: P !r)a )P #r, i"1, 2,2, n . (5) i i i # In this way, it is possible to fit the searching of the optimum, because a new process is started with a new random population created into a reduced phenotypical space of the best chromosomes belonging to the initial step. 3.2. Elitism The operations of crossover and mutation can affect the best member of the population, without producing offspring in the next generation. The elitist strategy [20] fixes this potential source of losing by copying the best member of each generation into the succeeding generation. The elitist strategy may increase the domination speed exerted by a super individual on the population. However, it seems to improve the genetic algorithm performance [20]. 3.3. Niche and speciation Sometimes, when several criteria are used simultaneously there exist no possibility to combine them by using a single number. When dealing with this situation, the problem is said to be a multiobjective or multicriteria optimization problem. This sort of problems do not have a unique optimal solution in the strict sense of the word, since the algorithm provides the best way to combine the criteria in order to obtain a plausible solution, but it is not the only one. The final selection depends on the designer’s criterion. When using genetic algorithms, the way to deal with multicriteria problems is by employing niche-like and species-like behavior in the members of the population [11,24,25]. In order to induce stable sub-population of strings (species) serving different sub-domains of a function (niches), it is then necessary to define a sharing function to determine the neighborhood and degree of sharing for each string in the population [25]. First of all, it is necessary to count the number of individuals belonging to the same niche. This is carried out by accumulating the total number of shares, which is determined by summing up the
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
237
sharing function values provided by all other individuals in the population, in a factor called niche count (m ) i m " + Sh[d(i, j)], (6) i j|P01 where d(i, j) is the distance between individuals i and j and Sh[d] is the sharing function. Sh[d] is a decreasing function of d(i, j), so that Sh[0]"1 and Sh[d*p ]"0. Usually, a triangular 4)!3% sharing function is used [25]: d Sh[d]"1! for d)p , 4)!3% p (7) 4)!3% Sh[d]"0 for d'p , 4)!3% here p is the niche radius, defined by the user at some estimates of the minimum separation 4)!3% desired or expected between the desired species. Then, the selection of the best no-dominating individuals to go on reproduction will be made over the one that has the less niche count m . Thus, the pareto frontier will be formed by a wide i variety of individuals.
4. Bidimensional and tridimensional models optimization The optimization process carried out by genetic algorithms can be applied to a different optimization engineering problems. In order to test its robustness, efficiency and flexibility, we will discuss herein the application of GAs to 3D-truss structures, due to their easy interpretation compared with other optimization techniques. Following, the optimization of bidimensional finite element models will be discussed. 4.1. The elastic problem In elasticity, design specifications and restrictions must be imposed in a model by defining limit values for stresses and displacements. In order to obtain those values in each design a basic finite element formulation is used to analyse 2D FEM models, as well as a classic displacement-based matrix method is used to analyze 3D truss structures. The basic formulation of the FEM method in elasticity can be found elsewhere [26,27]. 4.2. Chromosome representation In the examples presented herein, 3D truss and 2D finite-element models of mechanical pieces, the nodal restraints, loading conditions and the topological configuration are known. Hence, the design variables are f For truss structures: 1. the cross section type (t ) of each group of bars. j 2. the increment value (D ) of each joint coordinate (x ) h h
238
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
f For finite elements models: 1. the increment value (D ) of each joint coordinate (x ) in finite elements along the moving h h boundary. Thus, the structure characteristics or phenotypes can be represented as a chain of decimal parameters, as shown below: Structure"St t 2t 2t aD D 2D 2D mT, 12 j n 1 2 h n
(8)
where t is the the integer representing the cross sectional type, (D ) the real number representing j h the increment of the displacements, n the number of group of bars, n the number of mobile joint a m coordinates. The genotype structure is then obtained by encoding each parameter in binary format, in the following way: Individual"S1 0 0 0 1 0 1 1 0 0 1 120 1 1T.
(9)
4.3. The objective functions and constraints Once the population has been created, it is necessary to define the merit functions, which permit the measurement of the success or failure of the chromosomes of the problem to be solved. In the case of 3D-truss optimization, the problem is to find the minimum weight without violating the constraints M*/ N nj Min(¼)" + A + o ¸ j i i j/1 i/1
(10)
p !p61)0, i"1, 2,2, n , l"1,2, n , il "!3 #!4%4 p61"pt ° pc, i i
(11)
S !S61)0, S61"S5 ° S#, i
(12)
u !u61)0, k"1, 2,2, n , kl k r
(13)
s.t.
where A is the area of the section type j, o , ¸ : the density and length of the member i, p the stress j i i il of the member i under the load case l, pt, pc the upper bounds of tensile and compression stress of i i member i, S the slenderness ratio of bar i, S5, S# the upper slenderness bounds for tensile and i compression bars, º the elastic displacement of the k degree of freedom under load case l, u61 the kl k upper o lower limit of the elastic displacement of k d.o.f., and n the number of restricted r displacements. In the case of 2D finite element shape optimization, the problem consists in finding the minimum area of the model without violating the stresses and displacements limit values. The objective
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
239
function is then stated by N %-
Min(¼)"+ » o . i ii 1
(14)
In multiobjective optimization, the above objective function is used together with the minimization of the stress concentration factor (k) p Min(K)" .!9 , p .%$
(15)
where » , o are the volume and specific weight of element i. N the number of finite elements of the i i %model, p , p the maximum and mean Von-Mises stress of the model. subjected to .!9 .%$ (a) Stress restrictions: the Von-Mises stresses, calculated at Gauss points of the finite element, must not exceed the limit value p !$. p i!p )0. 70/ !$.
(16)
The Von-Mises stress is calculated as p "[(p !p )2#3q2 ]1@2. 70/ x y xy
(17)
(b) Nodal coordinates restrictions: nodal coordinates of some nodes of the mesh should not move beyond certain limit values in X and ½ directions in order to maintain the mesh topology. (c) Restrictions on the shape of the elements of the moving boundary, in order to avoid singular or negative Jacobians. In order to transform the constrained problem in an unconstrained problem and due to the fact that GAs do not depend on continuity and existence of the derivative, penalty methods have been used in this paper [11,18,23,31].Thus, the evaluation function is the minimization of the model weight, penalized with the constraints g "¼#jD*pD2#kD*ulD2#gD*xlD2, i
(18)
where g is the penalized weight, i the structured index, D*pD the increment over allowed stress, D*ulD the increment over allowed displacement, D*xlD the increment over allowed coordinate, (l"1,)2,)3). The parameters j, k and g are adjusted by trial and, in this paper, they have been evaluated in such a way that a 10% of increase in stress or displacement increases the structure weight by about 5% [18]. 4.4. Fitness function and normalization The duality of cost minimization and profit maximization is well-known. In normal operations, the transformation of a minimization problem into a maximization problem is simply carried out by multiplying the cost function by minus one.
240
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
In this paper, the minimization of the evaluation function has been transformed into a maximization search by using Constant fitness(i)" . g(i)
(19)
This new evaluation function is called fitness function. The relative rate between chromosomes is maintained as defined by the original evaluation function (Eq. (18)). When working with GAs, an early trend on few members to dominate the selection process is observed. therefore, the optimum solution cannot be reached. As a consequence, when the population is very convergent the average fitness is close to the best. This fact leads to a strong competition between many individuals displaying very similar fitness, thus generating a somewhat erratic behaviour of the algorithm. In both cases, it is necessary to normalize the fitness function in order to keep the competition more or less the same in generation after generation. This situation can be controlled by the normalization number f N " .!9 , (20) n f .%$ where f is the maximum fitness, f the averages fitness. Goldberg [18] suggested that these .!9 .%$ values should be taken between 1.2 and 2 for small population (n "50—100). Here, N "2 has 101 n been taken for all generations, and the Forrest transformation has been adopted in order to obtain it [28].
5. Illustrative examples 5.1. Three-dimensional truss steel tower The topology and dimensions of a three-dimensional 25-bar-truss [29,30] are shown in Fig. 3. The goal is to optimize the weight of the structure by considering the following cases: 1. Continuous optimization of the design variables such as cross section types and geometry configuration ( joint coordinates). 1.1. optimization of cross-sectional types with fixed geometry configuration, 1.2. optimization of cross-sectional types and variable geometry configuration. 2. Discrete optimization of the cross-sectional areas of the members. Table 5 displays the sections used to optimize the structure [30]. Tables 1—4 contain the design properties, loading cases and the design variables chosen for optimization with the genetic algorithm (see Table 5). Fig. 4 illustrates the best three runs of the weight evolution of the truss versus the number of generations. Those analyses have considered stress, displacement and buckling restrictions, by assuming a fixed geometry. Note that the initial weight was about 2333.77 N and the final one was 1696.45 N.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
241
Fig. 3. Three dimensional 25-bar transmission tower. Node and member numeration. (dimensions in mm). Table 1 Design properties (three-dimensional 25-bar truss tower) Young modulus Density Maximum tensile or compression stress Maximum displacement (nodes 3 and 4) Maximum buckling stress
68 950 Mpa 0.02713 N/cm3 275.8 Mpa $0.89 cm Following the AISC methodology
Table 2 Loading cases (three-dimensional 25-bar truss tower) Hypothesis
Node
F (kN) x
F (kN) y
F (kN) z
1
1 2 3 4
2.224 2.224 4.448 4.448
0.0 0.0 44.48 44.48
0.0 0.0 !22.24 !22.24
2
3 4
0.0 0.0
88.96 88.96
!22.24 !22.24
Fig. 5 also shows the weight evolution of the truss versus the number of generations, but in this case the variation of the geometry of the truss was included into the optimization criteria. It can be observed that the final weight, 974 N, reflects a reduction of more than 700 N with respect to the fixed geometry case.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
242
Table 3 Physical design variables (three dimensional 25-bar truss tower) Variable
Element
Minimum value (cm2)
1 2 3 4 5 6 7 8
1 2345 6789 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25
0.645 0.645 0.645 0.645 0.645 0.645 0.645 0.645
Table 4 Geometric design variables..(three dimensional 25-bar truss tower) Variable
Coordinate (nodes)
» ($cm) .!9
» ($cm) .*/
9 10 11 12 13
Z(1,2,5,6) X(1,2,5,6) Y(1,2,5,6) X(7,8,9,10) Y(7,8,9,10)
254 95.25 95.25 254 254
25.4 25.4 25.4 25.4 25.4
Table 5 Member types used for discrete optimization (three-dimensional 25-bar truss tower) Type
Area (cm2)
Type
Area (cm2)
Type
Area (cm2)
1 2 3 4 5 6 7 8 9 10
0.645 1.290 1.935 2.580 3.225 3.870 4.516 5.161 5.806 6.451
11 12 13 14 15 16 17 18 19 20
7.096 7.741 8.387 9.032 9.677 10.322 10.967 11.612 12.258 12.903
21 22 23 24 25 26 27 28 29 30
13.548 14.193 14.838 15.483 16.129 16.774 18.064 19.356 20.645 21.935
The results of the three best runs when dealing with discrete optimization are shown in Fig. 6. In this case, the optimization was carried out including stress and displacement restrictions. In the Fig. 7, the buckling stress of the bars are compared against the maximum value allowed for them. As it can be noted, all stresses in the bars are below this limit.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
243
Fig. 4. Genetic history with fixed geometry (Three dimensional 25-bar truss tower).
Fig. 5. Genetic history with variable geometry (Three dimensional 25-bar truss tower).
Table 6 summarizes the results obtained in the different cases under study. Finally, in Fig. 8 we can observe the initial and final geometry configuration of the tower. Note that in the final ground-plan view the algorithm has reduced the distance between supports along the x-direction and how in the final lateral view, the first level of the tower has moved upwards 5.2. 2D finite element models The following examples illustrate the ability of the algorithm to optimize finite element model shapes, in order to get the minimum area and an uniform stress distribution.
244
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 6. Genetic history with discrete optimization (Three dimensional 25-bar truss tower).
Fig. 7. Comparison between buckling stress in members and the maximum allowed stress (three dimensional 25-bar truss tower).
5.2.1. Fillet bar The first example corresponds to a fillet bar subjected to uniform loading. The dimensions and geometry of the model are shown in Fig. 9. The loading and boundary conditions are displayed in Fig. 10. The goal is to find the best shape of the boundary AB which produces a Von-Mises stress less than or equal to p "620.55 KPa. In Fig. 11 we can see the final shape of the model, where it can # be noted the optimized AB contour. The comparison between the initial and final stresses is shown in Figs. 12 and 13. Note how the final stress is very close to the allowed value.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
245
Table 6 Final results of the different runs (three-dimensional 25-bar truss tower) Variable active restriction
1 2 3 4 5 6 7 8
Fixed geometry (stress/displ/buckling) (cm2) 3.193 8.909 15.987 1.742 0.935 6.561 4.219 4.361
Final Weight (N) 1698.18
Variable geometry (stress/displ/buckling) (cm2) 2.380 14.006 8.335 1.0064 1.735 2.632 4.206 3.703 973.88
Fig. 8. Final three dimensional 25-bar truss tower configuration.
Discrete optimization (stress/displ/buckling) (cm2) 0.645 10.967 5.806 0.645 1.290 0.645 0.645 0.645 787.22
246
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 9. Fillet bar geometry and dimensions.
Fig. 10. Loading and boundary conditions.
Fig. 11. Final shape of the optimized piece.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
247
Fig. 12. Initial stress distribution (stresses in KPa).
Fig. 13. Final stress distribution (stresses in KPa).
The stress concentration factor changed from K "1.830 in initial configuration to K "1.5635 * & in the final one. This fact can be clearly observed in the final stress distribution. Fig. 14 shows the evolution of the maximum stress (p ) divided by the allowed stress (p ) and .!9 # the evolution of the weight of the model divided by the initial weight, versus the number of
248
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 14. Genetic evolution of stress and weight.
Fig. 15. Left window: smooth shape. Right window: final smooth stress.
generations. Those quantities factors are below one, which is an indication of the quite confidence of the reached solution. In Fig. 15, a smooth shape of the optimized piece is shown. It was obtained modifying the final shape shown in Fig. 11. The aim is to show how the designer can act interactively with the program
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
249
and modify the contour given by the optimization process. Moreover, it can be seen how the final stress in the smoothed piece are similar to those obtained in the optimized model (Fig. 11). 5.2.2. Cross section of a pipe Fig. 16 shows a quarter part of a cross section of a pipe, which is loading with an internal constant pressure of P "9810 KPa. The goal is to find the best cross section with the minimum 1 weight, by providing that the following constraints are not violated: 1. Inner radius. 2. Thickness of the pipe’s wall. 3. Maximum allowed Von-Mises stress. In the Fig. 17, the final mesh of the model is displayed. It corresponds to the best design obtained by the algorithm. It can be observed how the optimization process has smoothed the piece’s contour, thus getting a circular shape, as expected. Figs. 18 and 19 display a comparison between the initial stresses and the final ones. It can be noted an important reduction in the stress values. The maximum stress moved down from 2374.02 to 1265.50 MPa, which represents almost a reduction of 50% in the maximum stress. Fig. 20 shows the weight evolution of the model, versus the number of generations. The evolution of the maximum stress (p ) divided by the allowed stress (p ) is depicted in .!9 # Fig. 21. Note how this factor is lower than one. 5.2.3. Plate subjected to tractions As a demonstration of effectiveness of shape optimization using GAs, the optimization of the surface area of a plate is considered. The Fig. 22 contains the geometry, dimensions and the moving boundary AB of the plate. The loading and boundary conditions are depicted in Fig. 23
Fig. 16. Initial geometry, loading and boundary conditions (cross sectional pipe problem).
250
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 17. Top window: Zoom on the final geometry reached by the algorithm; bottom left window: Initial geometry of the model; bottom right window: final geomety of the model.
Fig. 18. Initial Von Mises stress distribution (stresses in MPa).
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
251
Fig. 19. Final Von Mises stress distribution (stresses in MPa).
Fig. 20. Genetic history of cross sectional pipe problem.
The objective function of the problem is the optimization of the surface area of the plate. The constraints are the Von-Mises stresses, which must not exceed 34 475 kPa. The initial area of the shape was 619.35 cm2, as shown in Fig. 22. The optimized shape is shown in Fig. 24. It can be noted how the plate area was reduced significantly. The final area obtained was 456.65 cm2, which means a 26.27% reduction. The comparison between the initial and final stresses is shown in Figs. 25 and 26. Note that the final shape is the result of the interaction between allowed stress and the elements stresses in the narrow zone of the final shape.
252
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 21. Evolution of the tension factor (p /p ) versus number of pipe problem). .!9 #
Fig. 22. Plate with moving boundary (A-B), geometry and dimensions.
Finally, the Fig. 27 illustrates the evolution curves of the structure stresses and weight. It can be noted the slope decreasing in the adimensional weight factor as the generation process progresses. This factor is always below than one, which clearly shows that the algorithm is stable and will converge to the optimal solution in terms of minimum weight. Another subject dealing with the relation between both curves should be remarked herein. When the adimensional stress curve displays its minimum value (see Fig. 27, generation number 17), it is a clear indication that the area
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
253
Fig. 23. Loads and boundary conditions.
Fig. 24. Final shape of the optimized piece.
of the piece can be further reduced. However, beyond generation number 40 or so, the adimensional stress curve displays a nearly constant value close to one, which involves that the area is very close to its optimal value. This explain why the adimensional weight curve behaves in very stable manner after generation 40.
254
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
Fig. 25. Initial stress distribution (plate subjected to traction).
Fig. 26. Final stress distribution (plate subjected to traction).
6. Conclusions The deterministic classic resolution in optimization methods starts from a predefined point and, from it the search is performed until the nearest minimum value is achieved. However, there is the possibility that such a minimum value is not the required point. When using genetic algorithms, the search starts from a set of points which form the initial population and from them the exploration
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
255
Fig. 27. Genetic evolution of stresses and weight (plate subjected to traction).
of the design space of plausible solutions is started. It is carried out by the random application of crossover and mutation operators. One of the greatest advantages of the method is that it is a discrete method. The power of the GA, based in the Darwin mechanics and having certain random behavior, is due to the exponential increment of the short and low degree “schemata” associated to good merit functions. Its discrete formulation and its ability in searching, suggest that this method is a powerful tool in engineering applications. New algorithms to simulate the process artificially have been developed and applied to various types of engineering optimization problems. Studying the algorithms from the optimization point of view, they surpass the other existing optimization techniques in some points, and they require no calculation of sensitivities and comparatively to attain the global optimum point. The 2D finite element models optimized herein have shown that the proposed approach is suitable to carry out this task. The results indicate a good performance even in the presence of stress concentrations caused by corners and geometry changes. A promising field in structural optimization is currently under research. The application of GA methods to finite element models suggest the exciting possibility of optimizing complex 3D problems, such as the shape optimization of gravity dams.
Acknowledgements The authors wish to acknowledge the financial support provided by Consejo de Desarrollo Cientifico y Humanistico (CDCH) and Consejo Nacional de Investigaciones Cientificas y Tecnolo´gicas (CONICIT). Also, to H. Bricen8 o for her careful revision of the english manuscript.
256
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
References [1] V. Braibant, C. Fleury, Sensitivity analysis in shape optimal design, in: De. C.A. Mota (Ed.), Computer Aided Optimal Design, vol. 2, Troı´ a Portugal, 1986, pp. 7—28. [2] O.C. Zienkiewicz, J.S. Campbell, Sequential linear programming, in: De. R.H. Gallagher, O.C. Zienkiewicz (Eds.), Optimum Structural Design, Wiley, New York, 1973, pp. 109—126. [3] C.V. Ramakrishnan, A. Francavilla, Structural shape optimization using penalty functions, J. Struct. Mech. 3 (4) (1975) 403—432. [4] A. Francavilla, C.V. Ramakrishnan, O.C. Zienkiewicz, Optimization of shape to minimize stress concentration, J. Strain Anal. 10 (2) (1975) 63—70. [5] S. Bhavikatti, C.V. Ramakrishnan, Optimum design of fillets in flat and round tension bars, ASME Paper 77-Det-45, 1977. [6] S. Bhavikatti, C.V. Ramakrishnan, Shape optimization of structural systems using finite elements and sequential linear programming, in: A. Wexley (Ed.), Large Engineering Systems, Pergamon Press, Oxford, 1977, pp. 224—225. [7] Y.W. Chun, E.J. Haug, Two dimensional shape optimal design, Int. J. Numer Meth. Eng. 13 (5) (1978) 311—336. [8] Y.M. Yoo, E.J. Haug, K.K. Choi, Shape optimal design of multiple loaded structural components, in: H. Eschenaver, N. Olhoff (Eds.), Optimization Methods in Structural Design. pp. 250—255, 1983. [9] R. Yang, K.K. Choi, E.J. Haug, Finite element computation of structural design sensitivity analysis. Report no. 3, University of Iowa, 1984. [10] J. H. Holland, Adaptation in Natural and Artificial Systems, University of Michigan Press, Michigan, 1975 and First MIT Press ed., Cambridge, Ch. 7, 1992, p. 121. [11] D.E. Goldberg, Genetic Algorithms in Search Optimization and Machine Learning, Chaps. 3 and 4, AddisonWesley, Reading, MA, p. 89, 1989. [12] L. Davis, Genetic Algorithm and Simulated Anneling, Morgan Kaufmann, Los Altos, CA, 1987. [13] N.A. Barricelli, Numerical testing of evolution theories, ACTA Biotheoretica 16 (1962) 69—126. [14] A.S. Fraser, Simulation of genetic systems, J. Theoret. Biol. 2 (1962) 329—346. [15] F.G. Martin, C.C. Cockerham, High speed selection study, in: O. Kemptthorne (Ed.), Biometrical Genetics, Pergamon Press, London, 1960, pp. 35—45. [16] D.E. Goldberg, Computer-aided gas pipeline operation using genetic algorithm and rule learning. Doctoral Dissertation, University of Michigan. Dissertation Abstracts International 44 (10), 317413 (University mycrofilms num. 8402282), 1983. [17] D.E. Goldberg, C.H. Kuo, Genetic algorithms pipeline optimization, J. Comput. Civil Eng. 12 (1987) 128—141. [18] D.E. Goldberg, M.P. Santani, Engineering optimization via genetic algorithm. Proc. 9th Conf. On Elect. Comp., ASCE, New York, 1986, 471—482. [19] A. Brindle, Genetic algorithms for function optimization, Doctoral Dissertation, University of Alberta, Edmonton, 1981. [20] L. Davis, Handbook of Genetic Algorithms, Van Nostrand Reinhold, New York, 1991. [21] J.H. Holland, Processing and processors for schemata, in: E.L. Jacks (Ed.), Associative information Processing, Elsevier, New York, 1971, pp. 127—146. [22] J.H. Holland, Schemata and intrinsically parallel adaptation. Proc. NSF Workshop of Learning System Theory and its Applications, Gainesville, University of Florida, 1973, pp. 43—46. [23] M. Galante, Genetic algorithms as an approach to optimize real-world trusses, Int. Numer. Meth. Eng. 39 (1996) 361—382. [24] J. Horn, N. Nafpliotis, Multiobjective optimization using the niched pareto genetic algorithms. Illinois genetic laboratory, Report No. 93005, July 1993. [25] D.E. Goldberg, J. Richardson, Genetic algorithms with sharing for multimodal function optimization, Proc. 2nd Int. Conf. on Genetic Algorithms, 1987. [26] O.C. Zienkiewicz, The Finite Element Method in Engineering Science, McGraw-Hill, New York, 1979. [27] S.P. Timoshenko, J.N. Goodier, Theory of Elasticity, McGraw-Hill, New York, 1951.
¼. Annicchiarico, M. Cerrolaza / Finite Elements in Analysis and Design 29 (1998) 231—257
257
[28] S. Forrest, A study of parallelism in the classifier system and its application to classification in KL-ONE semantic networks, Doctoral dissertation, University of Michigan, Ann Arbor, 1985, unpublished. [29] P. Marti, P. Company, M. Sanchis, Acoplamiento elementos finitos-te´cnicas de Optimizacio´n en el Sistema Disseny, II Cong. Met. Num. Ing., Espan8 a, 1987. [30] S. Rajeev, C.S. Krishnamoorthy, Discrete optimization of structures using genetic algorithms. J. Struct. Eng. 118 (5) (1992). [31] M. Cerrolaza, W. Annicchiarico, Structural optimization algorithms based on Genetic Simulation. Ed.CDCHUCV, 1996 (in Spanish) pp. 163.