## 1. Introduction

A framework or skeletal structure is one of the most used structures in engineering applications. Using such a structure is said to be advantageous since they are simple and inexpensive to construct. It can be employed in many engineering purposes e.g. transmission towers, wind turbine towers, communication towers, civil engineering structures, and mechanical parts. In the past, the design of such a structure was usually carried out in such a way that the conceptual design stage was somewhat ignored and the initial shape of a structure had been formed by an experienced engineering designer. The classical civil engineering design approach is then applied in the preliminary and detailed design stages.

Recently topology optimisation, an efficient tool for structural conceptual design, has been studied and used in a wide range of engineering applications. Considerable research work on the conceptual design of truss and frame structures by means of topology optimisation has been conducted in the last two decades (e.g. Hajela & Lee,1995; Ohsaki, 1995; Zhou, 1996; Yunkang & Xin, 1998; Hasancebi & Erbatur, 2002; Kawamura et al., 2002; Stolpe & Svanberg, 2003; Ohsaki & Katoh, 2005; Achtziger & Stolpe, 2007; Svanberg & Werme, 2007; Hagishita & Ohsaki, 2009). By the use of such design technology, the optimised layout of a skeletal structure can be accomplished. Having a structural layout from this design process, the shape and sizing design processes are then implemented. Nevertheless, a more efficient design approach can be achieved if a design problem is posed to have topological, shape, and sizing design variables at the same time. Much research work has been conducted towards combining and performing topology, shape, and sizing optimisation at the same optimisation run (e.g. Wang et al., 2004; Zhou et al., 2004; Tang et al., 2005; Shea & Smith, 2006; Achtziger, 2007; Martínez et al., 2007; Chan & Wong, 2008; Noiluplao & Bureerat, 2008; Zhu et al., 2008). The optimum design problem of the simultaneous topology/shape/sizing optimisation is usually structural weight minimisation subject to stress and other safety constraints. The optimisers employed can be a derivative-based method (Martínez et al., 2007), an evolutionary algorithm (Tang et al., 2005; Shea & Smith, 2006; Noiluplao & Bureerat, 2008), and a hybrid method (Chan & Wong, 2008). From the literature, most if not all of the optimisation problems are formulated to have one objective function, whereas in a practical design of an engineering system, there usually are several design objectives for decision making. If we have an optimiser that can deal with the design problem with such multiple design objectives, it could be useful and beneficial for designers.

Evolutionary algorithms (EAs) have been developed and used for several decades. The methods can be viewed as optimisation methods that search for optimum solutions by imitating some kinds of natural behaviours, and relying on a random process. Most EA searches are based upon using a population or a group of design solutions; consequently, they are often called population-based optimisers. Compared to their gradient-based counterparts, the EAs have advantages in that their searching process can be performed without using function derivatives. The methods are robust and simple to use, and can deal with almost all kinds of optimum design problems. As they are dependent on randomisation, EAs however have some drawbacks since they have low convergence rates and lack search consistency.

Evolutionary algorithms can deal with both single- and multiple-objective optimisation problems. The best known single-objective evolutionary optimiser is the genetic algorithm. Some other established single-objective EAs are: evolution strategies, particle swarm optimisation, and population-based incremental learning. EAs that can deal with multiple objective functions are called multiobjective evolutionary algorithms (MOEAs). As they use a set of design solutions for searching, on each generation those solutions can be sorted to find the so-called non-dominated design solutions. The set of non-dominated design solutions are improved iteratively, and the final set is taken as an approximate Pareto optimal set. This is the most outstanding and attractive ability of MOEAs as they can search for a group of non-dominated design solutions within one optimisation run. While the other conventional gradient-based optimisers need to be run as many times as the number of Pareto optimal solutions required.

The use of EAs for structural optimisation has been investigated and studied by many researchers around the world. By using EAs, any kind of design variables, constraints, as well as objective function can be defined. With additional strong points i.e. the capability of reaching near global optima of EAs and their robustness, they are even more popular than their gradient-based counterparts. The methods can be used for structural weight or cost minimisation (Benage & Dhingra, 1994, Bureerat & Cooper, 1988), passive vibration alleviation (Keane, 1995; Ton et al., 1998; Moshrefi-Torbati et al., 2003; Alkhatib et al., 2004; Kanyakam & Bureerat, 2007), and other performance maximisation (Xu et al., 2003; Achtziger & Kocvara, 2007).

The work in this chapter, consisting of two parts, is concerned with the demonstration of implementing MOEAs on the optimisation of skeletal structures. In the first part, the performance comparison of a number of MOEAs including non-dominated sorting genetic algorithm (NSGAII), strength Pareto evolutionary algorithm (SPEA2), population-based incremental learning (PBIL), Pareto archive evolution strategy (PAES), and multiobjective particle swarm optimisation (MPSO) are employed to solve seven simultaneous shape and sizing design problems of two structures. The comparative studies are carried out by using the so-called *C*-indicator. Also, a new indicator *C* is proposed to be used along with the *C* parameter. The Pareto optimum results obtained from using the various MOEAs are compared and discussed while some of the best performers are obtained. The second part involves the application of MOEAs to a simultaneous topology/shape/sizing optimisation of a skeletal structure. The design demonstrations are performed and illustrated.

## 2. Multiobjective evolutionary algorithms

A typical multiobjecitve optimisation problem can be defined as:

find x such that

Subject to

*g*

_{ i }(x) 0

*h*

_{ i }(x) = 0

where x is a vector of *n* design variables, *f*
_{
i
} are the *m* objective functions, *g*
_{
i
} are inequality constraints, and *h*
_{
i
} are equality constraints.

If the problem has one objective, there will be one (global) optimum solution. In cases that the problem has more than one objective function, there will be a set of optimal solutions, which is called a Pareto optimum set. Figure 1 illustrates a particular 2-objective design problem where all the feasible solutions are plotted in the objective domain. The Pareto optimal solutions are the points on the bold frontier, which is called a Pareto frontier or Pareto front. Similarly to the 2-objective case, the Pareto front of a 3-objective design problem is a 3-dimensional surface and so on.

The basic concept of exploring Pareto optimum points via MOEA search is that, in each generation while a new population is created, non-dominated solutions are classified and carried on to the next generation. The term non-dominated solutions define the local Pareto solutions among the members of the current population during evolutionary search. The definitions associated with non-domination (for minimisation) are given as:

*Definition 1* Dominance Given *f*
_{
i
}(x) for *i* = 1,…*m* are objective functions, if *f*
_{
i
}(x_{1}) *f*
_{
i
}(x_{2}) for every *i* {1,…,*m*} and there is *j* such that *f*
_{
j
}(x_{1}) < *f*
_{
j
}(x_{2}), then x_{1} dominates x_{2}.

*Definition 2* Non-Dominated Solutions (Approximate Pareto Set) Given a set of solutions or population G size *N*, a solution x_{
e
} G is a non-dominated solution in G if there does not exist x G such that x dominate x_{
e
}.

Figure 2 displays the plot of 7 design solutions on a 2-objective domain. The non-dominated solutions are x_{1}, x_{3}, and x_{7}.

In this work, five MOEAs are used including PAES, NSGAII, SPEA2, PBIL and MPSO. The methods are briefly reviewed as follows:

### 2.1. Pareto archive evolution strategy

PAES was proposed by Knowles & Corne, 2000. The simplest version of PAES is the combination of (1+1)-ES with a Pareto archiving strategy. The algorithm starts with an initial solution called a parent and an external Pareto archive. A candidate solution is then created by mutating the parent. The parent will be replaced by the candidate based on a selection strategy. The decision to discard or add the candidate to the Pareto archive is made. In cases that the archive is full, the adaptive grid algorithm is activated to remove one of the non-dominated solutions from the archive. As the procedure continues, the Pareto archive is iteratively updated until reaching the termination criterion.

### 2.2. Non-dominated sorting genetic algorithm

NSGA was proposed by Srinivas & Deb, 1994, and later the upgraded version NSGAII was presented by Deb et al., 2002. Starting with an initial population, selection, crossover and mutation operators are then used for exploring Pareto solutions. After crossover and mutation are operated with the given probabilities, the children are obtained and a new population is selected by performing non-dominated and crowding distance sorting algorithms on a union set of the parents and their children. The population is updated iteratively until the termination criterion is met.

### 2.3. Strength Pareto evolutionary algorithm

SPEA was proposed by Zitzler & Thiele, 1999, whereas its second version SPEA2 was presented in Zitzler et al., 2002. The search procedure starts with an initial population and an external Pareto set. Fitness values are assigned to the population based upon the level of domination. A set of solutions are then selected to a mating pool by performing the binary tournament selection operator. The new population are created by means of crossover and mutation on those selected solutions. The new external Pareto solutions are the non-dominated solutions of the combination of the old external Pareto set and the new population. In cases that the number of non-dominated solutions in the Pareto archive exceeds its predefined size, the truncation operator also known as the nearest neighbourhood technique is executed to delete some design solutions from the archive. The Pareto archive is updated repeatedly until the termination criterion is met.

### 2.4. Multiobjective particle swarm optimisation

The particle swarm optimisation method uses real codes and searches for an optimum by mimicking the movement of a flock of birds, which aim to find food (Reyes-Sierra & Coello Coello, 2006). The search procedure herein is based on the particle swarm concepts combining with the use of an external Pareto archiving scheme. Starting with an initial set of design solutions as well as their corresponding particle velocities and objective function values, an initial Pareto archive is filled with the non-dominated solutions obtained from sorting the initial population. A new population is then created by using the particle swarm updating strategy where the global best solution is randomly selected from the external Pareto archive. Afterwards, the external archive is updated by the non-dominated solutions of the union set of the new population and the previous non-dominated solutions. In cases that the number of non-dominated solutions is too large, the adaptive grid algorithm (Knowles & Corne, 2000) is employed to properly remove some solutions from the archive. The Pareto archive is repeatedly improved until fulfilling the termination criteria.

### 2.5. Population-based incremental learning

Population-based incremental learning, probably the simplest form of an estimation of distribution algorithm (EDA), was first proposed by Beluja, 1994, for single-objective optimisation. The method was then extended for multiobjective optimisation (Bureerat & Sriworamas, 2007; Kanyakam & Bureerat, 2007). The algorithm of multiobjective PBIL starts with an external Pareto archive and initial probability matrix having all elements set to be 0.5. A set of binary design solutions are then created corresponding to the probability matrix while their function values are evaluated. The Pareto archive is updated by replacing its members with the non-dominated solutions of the new population and the previous Pareto archive. In cases that the number of non-dominated solutions exceeds the predefined archive size, the normal line method is activated to remove some members from the archive. The probability matrix and the Pareto archive are iteratively updated until the termination criteria are fulfilled.

## 3. Structural modelling

A structural dynamic model can be described as a structure being in a dynamic equilibrium state. It is the state at which the system has minimum potential energy (potential energy herein includes structural elastic potential energy, the work done by external forces, and the work due to inertial forces). The equations of motion basically comprise of kinetic energy, structural restoration (spring and damping) and external dynamic forces. By using the finite element approach, the structural dynamic model is represented by

where is the vector of structural displacements, M is a structural mass matrix, K is a structural stiffness matrix, and F is the vector of dynamic forces acting on the structure.

The computation is traditionally achieved by means of finite element analysis. With the given boundary conditions (say _{
b
} = 0) being given, Equation (2) can be partitioned as

where the subscript *b* indicates the known displacements and unknown reactions at the boundary conditions, and the subscript *a* denotes unknown displacements and predefined external forces.

Equation (3) can be rearranged leading to 2 systems of differential equations as:

andIn the cases of free vibration analysis, Equation (4) can be written as

By substituting

Solving such a system of equations leads to *N* natural frequencies
*N* is the size of the square mass and stiffness matrices. The orthogonality conditions can be written as

By using the proportional (Rayleigh) damping concept, a damping matrix can be introduced to the model yielding

where C_{
aa
} =M_{
aa
} + K_{
aa
}, and and are damping coefficients to be defined.

From equation (9), by substituting

where is the frequency of the input force and its output displacement, and

The relation between displacement response and input force can be expressed as

Figure 3 illustrates how to measure H(r,s) which represents the ratio of a displacement response at the r^{th} degree of freedom to a harmonic excitation at the s^{th} degree of freedom.

Furthermore, by defining force transmissibility (FT), denoted by T(), as the ratio of output harmonic reaction forces to the input external harmonic forces, it can be written as

By letting the reaction force be

In structural vibration design, FRF and FT determine structural merit. The lower FRF or FT at a particular frequency means the better structural vibration suppression design. Therefore, a design objective can be assigned in such a way that frequency responses at a frequency range of interest are minimised. Moreover, maximising structural natural frequency is an alternative criterion for design of structures under dynamic loadings.

Apart from dynamic analysis, structural static analysis can be carried out by using Equation (3) ignoring the structural kinetic energy or the mass matrix. The system of equations then becomes

The solutions of (15) are

andThe displacement vector is used for displacement constraints as well as for computing stresses on structural elements. The reaction F_{
b
} is also an important factor since the reaction force can affect a structural foundation. In an optimisation process, it is also common to use many load cases since, in real-world problems, there are many aspects of applied loads acting on one structure.

## 4. Numerical experiment

### 4.1. Design problems

The optimisation problems assigned in this study are to find a Pareto optimal set that optimises multiple objective functions subject to stress constraints, which can be expressed as

subject to

_{
i
}
_{
a
}, i = 1, …, N_{
e
}

where _{
i
} is a stress on the i^{th} element, _{
a
} is an allowable stress, and N_{
e
} is the total number of elements.

The multiobjective design problems presented here are similar to the work in Kanyakam & Bureerat, 2007, where the objective functions include structural mass, a natural frequency, an FRF crest parameter, and an FT crest parameter. Figure 4 displays the term FRF crest parameter, which determines the shaded area in a frequency range of interest. As we need

to reduce the magnitude of the FRF for vibration suppression, the frequency range is thus the interval that bounds a particular natural frequency. Seven design criteria are defined as:

f

_{1}: structural mass- f
_{2}: - f
_{3}: f

_{4}: FRF crest parameter at_{1}f

_{5}: mean value of FRF crest parameters at_{1},_{2}and_{3}f

_{6}: FT crest parameter at_{1}f

_{7}: mean value of FT crest parameters at_{1},_{2}and_{3}

where _{
i
} is the i^{th} natural frequency.

The sets of objective functions can be arranged as:

P1: min f

_{1}, min f_{2}P2: min f

_{1}, min f_{3}P3: min f

_{1}, min f_{4}P4: min f

_{1}, min f_{5}P5: min f

_{1}, min f_{6}P6: min f

_{1}, min f_{7}P7: min f

_{1}, min f_{5}, min f_{7}.

As a result, we have seven multiobjective design problems with the above sets of objective functions and stress constraints. It should be noted that some of these objective function sets have been implemented and examined in Kanyakam & Bureerat, 2007, and it is shown that using P2, P4 and P6 leads to a more effective design based on vibration alleviation. All seven design problems are however used to benchmark the performance of the MOEAS.

We apply the above-mentioned problems to design three frame structures made of material with E = 20910^{9} N/m^{2}, and = 7000 kg/m^{3}. The structures are named as:

ST1: a 2D bridge (Figure 5). The structure has 12 nodes and 21 elements where nodes 1 and 12 are fixed. The vertical loads are applied to nodes 2, 4, 6, 8, and 10. The design variables consist of element diameters, and the vertical positions of nodes 3, 5, 7, 9, and 11. Both shape and sizing design variables are treated as being symmetric; as a result, there are 14 continuous design variables. The FRF used as design criteria is the point receptance in the vertical direction of node 7. The FT is the vertical direction transmission from node 7 to the fixed support node 1.

ST2: a walking tractor handlebar (Figure 6). The finite element model of the handlebar structure (Kanyakam & Bureerat, 2007) consists of 16 nodes and 27 elements with four nodes being fixed. The structure is subjected to two static load cases, one is the load for turning the tractor and the other is the load for balancing and controlling the tractor. The total number of continuous shape and sizing design variables is 24. The FRF being used is the point receptance in the vertical direction of node 6. The FT is the vertical direction transmission from node 6 to the fixed support node 1.

ST3: a 2D bridge (Figure 7). The ST3 structure is modified from ST1 for the demonstration of simultaneous topology/shape/sizing optimisation. The structure in Figure 7, which is set as a ground structure for topology optimisation, has 12 nodes and 25 elements. The design variables are element diameters, and the vertical and horizontal positions of nodes 3, 5, 7, 9, and 11, which have a total of 18 variables after symmetry treatment. The diameters are discrete and allowed to have nearly zero size so that we can have various structural topologies. The position of FRF and FT measures are the same as that of ST1.

Finite element analysis is carried out by using the 2-node 3-dimensional 12 d.o.f. beam element. The ST1 and ST2 structures are used for the investigation of a comparative performance of MOEAs, while the ST3 structure is assigned for the demonstration of performing topology/shape/sizing optimisation at the same time.

### 4.2. Performance assessment

Seven multiobjective evolutionary strategies are employed in this study, which are designated as:

BPBIL: PBIL using binary codes with 0.05 mutation rate, and 0.2 mutation shift

BPAES: PAES using binary codes

BNSGA: NSGAII using binary codes with crossover and mutation rates as 1.0 and 0.5 respectively

BSPEA: SPEA2 using binary codes with crossover and mutation rates as 1.0 and 0.5 respectively

RSPEA: SPEA2 using real codes (Sirsomporn & Bureerat, 2008) with crossover and mutation rates as 1.0 and 0.5 respectively

RNSGA: NSGAII using real codes with crossover and mutation rates as 1.0 and 0.5 respectively

RMPSO: MPSO (Reyes-Sierra & Coello Coello, 2006) using real codes with a starting inertia weight, an ending inertia weight, a cognitive learning factor, and a social learning factor as 0.5, 0.01, 0.5, and 0.5 respectively.

For the comparative performance study of optimising ST1 and ST2, starting with the same initial population, the number of generations of MOEAs is set to be 100 while both the population and archive size are set as 50. Each optimiser is applied to solve each design problem for designing each structure with 10 optimisation runs. The last updated Pareto archive is set to be a Pareto optimal set. For the design of ST3, one of the best MOEAs obtained from the first two design case studies is employed to solve the P2, P4, and P6 design problems with the number of generations being 150 while the population and the external archive are sized 100. The non-dominated sorting scheme proposed in Deb et al. (2001) is used to handle the stress constraints.

The performance assessment is accomplished by using the C-parameter (Zitzler et al., 2000). Having two particular non-dominated fronts {A} and {B}, C_{
A,B
} determines the number of solutions in {B} that are dominated or equal to some solution in {A} divided by the total number of solutions in {B}. On the other hand, C_{
B,A
} determines the number of solutions in {A} that are dominated or equal to some solution in {B} divided by the total number of solutions in {A}. From the definition, if C_{
A,B
}> C_{
B,A
}, we can say that front {A} is better than front {B} or vice versa. Figure 8 displays two approximate Pareto fronts {A} and {B} where C_{
A,B
} = 0.5000 and C_{
B,A
} = 0.4000. This implies that {A} is better than {B}. However, we can observe from the shaded areas between the two fronts dominating each other that the dominating area of {B} is larger than {A}. These shaded areas are also meaningful for front comparison. Nevertheless, it is difficult or even impossible to calculate such areas of hundred pairs of Pareto fronts directly. We can laterally calculate them and define a new performance indicator as:

where

andIf V_{
A
} = V_{
B
} = 0, C_{
A,B
} = C_{
B,A
} = 0.5.

_{ A,B }> C

_{ B,A }, it means {A} is better than {B}. From the definition, the C-indicator also has the effects of front extension and the number of dominant points.

From Figure 8, we have C_{
A,B
} = 0.4556 and C_{
B,A
} = 0.5444, which means {B} is the better front than {A}. Based on these two indicators, we have two different views of comparing approximate Pareto fronts. The special characteristic of the C parameter is illustrated in Figure 9, which shows three particular non-dominated fronts {A_{1}}, {A_{2}}, and {A_{3}}. The C and C comparisons of the three fronts are given in Table 1. In the Table, the value in the i^{th} row and j^{th} column represents C_{
Ai,Aj
}. From the C comparison, it can be concluded that the front {A_{3}} is the best method while the second best is {A_{1}}. However, in cases of C comparison, we can see that {A_{2}} is better than {A_{1}}, {A_{1}} is better than {A_{3}}, and {A_{3}} is better than {A_{2}}. This scenario can be called a scissor-paper-stone situation, which can happen in comparing approximate Pareto fronts obtained from using MOEAs. The two indicators, in fact, are not always conflicting to each other but they should be used together to provide different aspects of front comparing.

## 5. Design results

The performance comparison of 7 MOEAs is displayed in Figures 10 – 17. In Figure 10, the box-plots of C parameters comparing the various MOEAs when solving the P1-P7 design problems for the case of ST1 are illustrated. Each box-plot represents 1010 C values of comparing 10 fronts obtained from using a method X to 10 fronts obtained from using a method Y when solving a particular design problem. The box-plots at the i^{th} row and j^{th} column display the C values comparing non-dominated fronts of the seven design problems obtained from the i^{th} optimiser with respect to those obtained from using the j^{th} optimiser. For example the box-plots at the first row and second column in the figure present the C values of the fronts obtained from using BNSGA to those obtained from using BSPEA. We also need the box-plots in the second row and first column to compare these two optimisers. From the Figure, it can be concluded that BPBIL is the best for all the design problems. The C comparison of these four methods is given in Figure 11. It can be shown that BPBIL is the best method except for the P7 design problem, which BNSGA gives the best results.

The best multiobjective evolutionary algorithm using binary codes BPBIL is taken to be compared with the methods using real codes as shown in Figures 12-13. From the C comparison in Figure 12, the overall best method is BPBIL with RSPEA being the second best. In Figure 13, it can be seen that RSPEA is as good as BPBIL based on C comparison. This occurrence is probably similar to that illustrated in Figure 8. However, when taking an account of both indicators, the best method for designing ST1 is BPBIL.

The comparative performance of MOEAs for solving the P1-P7 optimisation problems for the case of ST2 is illustrated in Figures 14-17. The best method that use binary codes based on the C indicator is BPBIL whereas the second best method is BNSGA as shown in Figure 14. Based on the C indicator, the best method is BPBIL while the close second best is BNSGA. BNSGA even outperforms BPBIL in cases of the P3 and P4 design problems. This situation is similar to that displayed in Figure 8.

The best evolutionary optimiser among the methods using binary codes (BPBIL) is taken to be compared with the optimisers using real codes as shown in Figures 15 and 16. Based on the C indicator, the best optimiser is RSPEA whereas the second best is BPBIL. For the C comparison, the best optimiser is still RSPEA with BPBIL being the second best. From both ST1 and ST2 case studies, it can be concluded that BPBIL tends to be more efficient when dealing with design problems with a smaller number of design variables while RSPEA is the better optimiser for the design problems with a greater number of design variables.

The BPBIL algorithm is chosen to solve the P2, P4 and P6 design problems where the structure to be optimised is ST3. In this design case, some structural elements are allowed to be removed if their diameters are too small. The Pareto optimum results of P2 for optimising ST3 obtained from using BPBIL are plotted in Figure 18. Some selected non-dominated solutions are labelled whereas their corresponding structures are given in Figure 19. It can be seen that the structures have various layouts and shapes as well as various element sizes. The approximate Pareto front of P4 is plotted in Figure 20 whereas some selected structures from the front are illustrated in Figure 21. Similarly to the P2 design case, there consist of various structural topologies, shape and element sizes. The non-dominated front of P6 obtained from using BPBIL is given in Figure 22 whilst the selected solutions are shown in Figure 23. Similarly to the first two cases, the structures have various topologies, shapes, and element sizes. The obvious advantage of using MOEA for solving the simultaneous topology/shape/sizing optimisation is that they are robust and simple to use, and we can have a set of Pareto optimal solutions for decision making within one simulation run.

## 5. Conclusions and discussion

The search procedures of the various MOEAs as well as the static and dynamic analysis of skeletal structures are briefly detailed. Seven multiobjective optimum design problems are posed for benchmarking the performance of the MOEAs. The design problems are applied to three skeletal structures. The first two structures have simultaneous shape/sizing design variables while the third structure is used for design demonstration of simultaneous topology/shape/sizing optimisation. From the optimum results, it can be concluded that SPEA2 using real codes is the best performer for the ST2 structure which has greater number of design variables. The PBIL using binary string is said to be the best method for ST1, which has a smaller number of design variables. From the design demonstration of using topological, shape, and sizing design variables at the same time, it is shown that BPBIL is an efficient and effective optimisation tool for such a type of structural optimum design. The C and C indicators should be used together to provide some insight that is missing from using either of them solely. The application of MOEAs for design optimisation of skeletal structures is said to be advantageous since they are robust and simple to use. The method can cope with all kind of design criteria as demonstrated in this work. Most importantly, we can have a set of non-dominated solutions for decision making within one optimisation.