Optimization of structures under load uncertainties based on hybrid genetic algorithm

This paper describes a technique for design under uncertainty based on hybrid genetic algorithm. In this work, the proposed hybrid algorithm integrates a simple local search strategy with a constrained multi-objective evolutionary algorithm. The local search is integrated as the worst-case-scenario technique of anti-optimization. When anti-optimization is integrated with structural optimization, a nested optimization problem is created, which can be very expensive to solve. The paper demonstrates the use of a technique alternating between optimization (general genetic algorithm) and anti-optimization (local search) which alleviates the computational burden. The method is applied to the optimization of a simply supported structure, to the optimization of a simple problem with conflicting objective functions. The results obtained indicate that the approach can produce good results at reasonable computational costs.


Introduction
Today with the development of industry and research, reliability is more and more emphasized. However in some analysis of engineering problem, uncertainty inevitably exists. For example, wind loading on structures often varies over certain range instead of one determined value. The design of such systems should take into account the intervals of deviation of variables from their nominal values. In some uncertainty cases, the variety in the system parameters is neglected for the convenience of analysis. In certain situation it is possible to obtain the response which is valid and reasonable. However, sometimes the uncertainty analysis of system can not be neglected because the uncertainty would significantly affect the system performance. Real structures are designed for numerous combinations of loads and load cases. In structure design it is typical that there may be several hundred important design load cases, which implies uncertainty in the loading conditions. A number of approaches have been developed to tackle such problems. According to the different types of uncertainty, several distinct methods were introduced to handle the analysis: fuzzy sets and fuzzy logic, interval analysis, and probabilistic approach. Fuzzy set theory was firstly introduced in (Zadeh 1978). The application of fuzzy sets concept into structural analysis was studied in (Ayyub 1997) systematically. In his work, a fuzzy set was utilized to describe every objective function. Some mechanical systems can be analyzed through this method (Qiu and Rao 2005). Fuzzy set theory has several advantages including the fact that the mathematics is formally defined and it can provide a potentially simple approach. The disadvantages of fuzzy set theory include validation necessary, justification of each step necessary, complexity when using many variables and the fact that the membership function may be heuristically defined. When less information is given and only range of each uncertainty quantity can be specified, interval method is widely used. Bounded uncertainty approach (Ben-Haim and Elishakoff 1990), a convex model of uncertainty was an extension to the interval analysis (Moore 1966). Some disadvantages of interval analysis are: (1) Interval arithmetic computations are slower than floating point operations, roughly by a factor of three, though there are problems that are solved faster when implemented using interval arithmetic; (2) There are no additive or multiplicative inverses for intervals; (3) There is no strict distributive law of multiplication over addition, only sub-distributivity; and (4)  For anti-optimization, the robustness of the objective function can be achieved by maximizing the minimum value of the objective functions. It tends to move the uncertainty parameters to the desirable range and results in higher reliability with respect to uncertainty. Therefore an anti-optimization procedure implemented in this work by local search is employed to search for the values of the worst objective functions. Consider a problem subject to u uncertain variables u x and nu − normal design variables n x . The formulation can then be written as below: where C N is a particular set of solutions which are related to generation number N. The formulation considered leads to a nested optimization problem which will be solved by means of a GA for optimization and a local search for anti-optimization. Generally speaking, this is a sort of minmax search, where the "max" part is dealt with by a local search algorithm, and the "min'' part is realized by a GA.

Tournament selection for local search
Since a local search strategy requires a tournament selection between an initial solution and its neighboring solution, a comparison strategy is needed. For multi-objective optimization without constraints, a single objective function converted from multiple objectives can be used. For constrained optimization, constraint handling mechanisms should be given first.
In most applications, penalty functions using static (Srinivas & Deb 1994), dynamic or adaptive concepts (Michalewicz & Schoenauer 1996) have been widely used. The major problem is the need for specifying the right value for penalty parameters in advance. The method from (Ray, Tai et al. 2001) incorporates a Pareto ranking of the constraint violations, so it does not involve any aggregation of objectives or constraints and thus the problem of scaling does not arise. Note that Pareto ranking is not well suited for hybridization with local search. For an anti-optimization problem, when Pareto ranking is used, the current solution x is replaced with its neighboring solution y (i.e., the local search move from x to y is accepted) only when x dominates y (i.e., x is better than y). That is, the local search move is rejected when x and y are non-dominated with respect to each other. However, change of the rank of a given solution may require significant changes of the objective/constraint values, thus, many local moves will not degenerate the rank. A constraint handling method (Deb 2000) was proposed which is also based on the penalty function approach but does not require the prescription of any penalty parameter. The main idea of this method is to use a tournament selection operator and to apply a set of criteria in the selection process. For an anti-optimization problem, it can be easily changed to: i. Any infeasible solution is preferred to any feasible solution.
ii. Between two feasible solutions, the one having worse objective function value is preferred. iii. Between two infeasible solutions, the one having bigger constraint violation is preferred.
According to these criteria, the constrained optimization can be constructed as where () f x is the artificial unconstrained objective function, F is the feasible region of the design domain, max f is the objective function value of the worst feasible solution in the population, and () vio x is the summation of all the violated constraint function values. However, this approach is only suitable for single-objective constrained optimization problem if no further handling mechanisms for multiple objectives are given. And () vio x as summation of constraint values cannot reflect the real relative comparison between them because of different orders of magnitude among the constraints, and in this sense, it is also based on the penalty functions where all the penalty parameters are set to 1. Extending the basic idea of Deb's method, a technique combining Pareto ranking and weighted sum is suggested in this work for the local search selection process. There are only 3 combinations for the two solutions: both feasible, both infeasible, and one feasible and the other infeasible. The main idea of the technique is to use a tournament selection operator and to apply a set of criteria in the selection process. For an anti-optimization procedure, any infeasible solution is preferred to any feasible solution. When both solutions are feasible, Pareto ranking based on objectives is calculated. The one with bigger rank value is preferred. If the situation still ties, a more sophisticated acceptance rule is used for handling the situation. The fitness function of the solution x is calculated by the following weighted sum of the m objectives: where () f x is a combined objective, and 12 ,,, m ww w are nonnegative weights for the objectives set according to different orders of magnitude among them. Constant weight values are used in this work to fix the search direction based on user's preference. The solution with a bigger () f x will survive. When both solutions are infeasible, Pareto ranking based on constraints is calculated. The one with bigger rank value is preferred. If the rank is same, the one with worse fitness value survives. A tournament selection criterion can be described as below to decide whether a current solution x s h o u l d b e r e p l a c e d b y a neighboring solution y: i. If x is feasible and y is infeasible, replace the current solution x with y (i.e., let x = y).
ii. If both x and y are feasible, then if RankObj RankObj iii. If both x and y are infeasible, then if RankCon RankCon

Selection of initial solutions
Local search applied to all solutions in the current population in the algorithm is inefficient, as shown in (Ishibuchi, Yoshida et al. 2002). In the proposed algorithm, the computation time spent by local search can be reduced by applying local search to only selected solutions in selected generations. If n is the number of decision variables, the best n number of solutions from the current population (based on Pareto ranking) are selected. These n number of mutated solutions and elites from N-th generation after local search are then put into the next population. The generation update mechanism in the proposed algorithm is shown in Fig. 1. The implementation of the anti-optimization part is modularized.

Local search procedure
As explained in the above, a local search procedure is applied to elite individuals and new solutions generated by the mutation in selected generations. Generally, a local search procedure can be written as follows: Step 1. Specify an initial solution and its corresponding design variable under uncertainty.
Step 2. Apply Hooke and Jeeves Method to determine the search path using the tournament selection criteria stated above as the function values.
Step 3. If the prescribed condition is satisfied, terminate the local search.

Main algorithm
The overall algorithm uses a framework which combines the method stated in ) and the local search proposed above. The algorithm is given below: Step 1. Generate random initial population P of size M .
Step 2. Evaluate objective as well as constraint functions for each individual in P.
Step 4. Select elite individuals. Elite individuals carried from the previous generation preserve the values of their objective and constraint functions.
Step 5. Select n best individuals from P, mutate and apply local search procedure in specified generation, then put them into new population P'.
Step 7. If a prescribed stopping condition is satisfied, end the algorithm. Otherwise, return to Step 2.

Enhanced geometric representation scheme for structural topology optimization
developed in (Tai & Akhtar 2005;Tai & Prasad 2007). As in any structural topology optimization procedure, the geometry of the structure has to be represented and defined by some form of design variables. The enhanced morphological representation efficiently cast structure topology as a chromosome code that makes it effective for solution via a GA. In the proposed scheme, the connectivities and the number of curves used are made variable and to be optimized in the evolutionary procedure. The process of the scheme definition is illustrated as follows.
A square design space shown in Fig. 2(a) is discretized into a 50 by 50 mesh of identical square elements. While it is initially unknown how the design space will be occupied by the structure, there must exist some segments of the structure such as the support and the loading that have functional interactions with its surroundings. The support point is some segment of the structure that is restrained (fixed, with zero displacement) while the loading point is where some specified load (input force) is applied to deform the structure. Collectively, the support and loading points represent the input points of the structure. There is also usually an output point which is some segment of the structure where the desired output behavior is attained. As shown in Fig. 2(a), the problem is defined with four I/O locations, each made up of one element in black. Six connecting curves in the illustration of Fig. 2(b), three of which are active and three of which are inactive, are used such that there is one connecting curve between any two points (i.e. every I/O point is directly connected to the other three). Before continuing, it is important to make a clear distinction between the active and inactive curves. The active curves are the curves which are in the 'on' state. The structure is generated based only on the active curves. Although the inactive curves, which are in the 'off' state, temporarily contribute nothing to the current structure, they are still very important in subsequent generations because they may be turned 'on' later through the crossover or mutation operations. In Fig. 2(b), the active curves are marked with thick lines and the inactive with thin dotted lines. The connectivity of the I/O points is based on all connecting active curves joining one point to another. Each curve is a Bezier curve defined by the position vector which can be derived from the element number of control point. The set of elements through which each active curve passes form the 'skeleton' (Fig. 2 (c)). Some of the elements surrounding the skeleton are then included to fill up the structure to its final form ( Fig. 2 (d)) based on the skeleton's thickness values. Each curve is defined by three control points, and hence each curve has four thickness values. And the union of all skeleton, surrounding elements and I/O elements constitute the structure while all other elements remain as the surrounding empty space. In order to use a GA for the optimization, the topological/shape representation variables have to be configured as a chromosome code. Hence the structural geometry in Fig. 2 (d) can be encoded as a chromosome in the form of a graph as shown in Fig. 3. Each curve is represented by a series of nodes connected by arcs in the sequence of start element number, thickness values alternating with control element number and end point. For identification purpose, the active curves are shown by solid lines and the inactive curves are represented by dotted lines. Altering the curve states can vary the connectivity of the I/O regions, and therefore the representation scheme can automatically decide the connectivity. The resulting scheme therefore increases the variability of the connectivity of the curves and hence the variability of the structure topology. Two of the important operations in a GA are the crossover and mutation. In this implementation, the crossover operator works by randomly sectioning any single connected subgraph from a parent chromosome and swapping with a corresponding subgraph from another parent (as shown in Fig. 3.). As a result, two offsprings are produced which have a mix of the topological/shape characteristics of the two parents, and the advantages of the representation (such as no checkerboard patterns and single-node hinge connections) are maintained in the offspring. The 'on' and 'off' state of different curves which are crossed by the loop are also swapped. If the 'on' variables dominate a curve, i.e. when the number of 'on' variables is more than the number of 'off' variables, the curve of in the child chromosome will be active. Otherwise, the child curve will be inactive. As for mutation, the mutation operator works by randomly selecting any vertex of the chromosomal graph and altering its value to another randomly generated value within its allowable range. Mutation about the on-off state is simple, which is altering the state of curves. When the selected curve is active, it will be inactive after mutation, and vice verse.  In summary, this morphological representation scheme uses arrangements of skeleton and surrounding material to define structural geometry in a way that will not render any undesirable design features such as disconnected segments, checkerboard patterns or singlenode hinge connections because element edge connectivity of the skeleton is guaranteed, even after any crossover or mutation operation. Any chromosome-encoded design generated by the evolutionary procedure can be mapped into a finite element model of the structure accordingly.

Target matching problem
Before a GA is relied upon for solving a structure design problem with unknown solutions, it is important that the performance of the GA be tested and tuned by using it to solve a problem with known solutions. Various kinds of test problems (Michalewicz, Deb et al. 2000;Schmidt & Michalewicz 2000;Martin, Hassan et al. 2004) have been established for testing multi-objective GAs. They were created with different characteristics, including the dimensionality of the problem, number of local optima, number of active constraints at the optimum, topology of the feasible search space, etc. However, all of these test problems have well-defined objectives/constraints expressed as mathematical functions of decision variables and therefore may not be ideal for evaluating the performance of a GA intended to solve problems where the objectives/constraints cannot be expressed explicitly in terms of the decision variables. In essence, a GA is typically customized to tackle a certain type of problem and therefore `general-purpose' test problems may not correctly evaluate the performance of the customized GA. The test problem should, therefore, ideally suit (or be customized to) the GA being used.
In numerous real-life problems, objectives/constraints cannot be expressed mathematically in terms of decision variables. One of such real-life problems is structural topology optimization, where a procedure (structure geometry representation scheme) first transforms decision variables into the true geometry of the designed structure and then finite element analysis of the designed structure is carried out for evaluating the objectives/constraints. The GA solving such problems may have special chromosome encoding to suit the structure geometry representation used and there may also be specially devised reproduction operators to suit the chromosome encoding used. As such, the structure geometry representation scheme, the chromosome encoding and the reproduction operators introduce additional characteristics to the search space and, therefore, they are very critical to the performance of the GA. The test problem for such GAs, therefore, must use the same structure geometry representation scheme, chromosome encoding and reproduction operators. The conventional test problems found in literature cannot make use of the GA's integral procedures such as structure geometry representation scheme and therefore they are not suitable for testing such GAs.
Ideally, the test problem should emulate the main problem to be solved. The test problem should be computationally inexpensive so that it can be run many times for the GA parameters to be changed or experimented with and the effect thereof can be studied for the purpose of fine-tuning the GA. However, the main problem in the present work, being a structural topology optimization problem under uncertainty, requires structural analysis which consumes a great deal of time. Taking the running time into consideration, the test problem needs to be designed without any need for structural analysis. A test problem emulating structural topology optimization does not necessarily need structural analysis as the main aim of topology optimization is to arrive at an optimal structural geometry. Without using structural analysis, if a GA is successfully tested to be capable of converging the solutions to any arbitrary but predefined and valid `target' structural geometry, then it may be inferred that the GA would be able to converge design solutions to the optimal structural topology when solving an actual topology optimization problem. Based on this inference, a test problem can be designed such that simple geometry-based (rather than structural analysis based) objectives/constraints help design solutions converge towards the predefined target geometry. This type of test problem may be termed as ``Target Matching Problem", which is capable of using exactly the same GA (including structural geometry representation scheme, chromosome encoding and reproduction operators) as that intended for solving the actual topology optimization problem. The present problem is similar to the Target Matching Problem solved in Tai et al. 2008;Wang & Yang 2009).
The target matching problems are defined here as multiobjective optimization problems under uncertainty which are more difficult (e.g. more nonlinear problem) and computationally intensive.

Formulation
The test problem makes use of the design space shown in Fig.5, which has one support point, two loading points and one output point. The problem does not represent a structural analysis problem, but the original terms 'support', 'loading' and 'output' are still used for ease of reference. The loading point 1 is positioned anywhere along the left boundary and loading point 2 is positioned anywhere along the right boundary. The position of output point is fixed as shown in Fig.5. The support point is positioned in a specified area marked as ''under uncertainty'' and its position is random in the area. In this problem, the target geometry is shown in Fig. 6. The aim is therefore to evolve structures that match as closely as possible this target geometry. The problem presented here is more difficult than the original problem described in , since the support point is under uncertainty that makes the geometry more complex and not easy to converge to the target. The problem is formulated with the following two objectives and two constraints: distance objective, material objective, forbidden area constraint and prescribed area constraint. Such a problem is defined with the help of Fig.7. The distance objective is given by where l d is the centroid-to-centroid Euclidean distance between the actual loading point 1 and the actual support point. The material objective is given by where i x is the material density of the i -th element in the design space, with a value of either 0 or 1 to represent that the element is either void or material (solid), respectively. n is the total number of elements in the discretized design space. In other words, this objective function is defined as the summation of the material density of all elements in the current geometry. The forbidden area constraint can be written as where i y is the material density of the i -th element in the forbidden area, and f n is the total number of elements in the forbidden area. In other words, the summation of the material density of elements in the forbidden area is required to be less than or equal to zero. The prescribed area constraint can be written as where i z is the material density of the i -th element in the prescribed area, and p n is the total number of elements in the prescribed area. In other words, the summation of the material density of elements in the prescribed area is required to be more than or equal to the total number of elements in that area.

Main results
The target matching problem to be solved is defined by the design space shown in Fig. 5. The optimization was run for 1001 generations with a population size of 100 per generation. The local search procedure is triggered once every ten generations. By the end of evolutionary process 132,113 objective function evaluations have been performed. One of the solutions at the end of 1001 generations is shown in Fig. 8. It is the same as the target solution shown in Fig. 6. As can be seen from the result, the support point is on the extreme point where the element number is 37. The following Fig. 9 illustrates how the solution shown in Fig. 8 is obtained by applying local search. Apply the Hooke and Jeeves method to determine the search path using the tournament selection criteria stated in Section 2. Each data point is labeled with its index where some indexes are coincided. At the start point,  Fig. 11. History of the best distance objective ( material f ). Fig.12 shows a plot in objective space, where the solid shape markers are used to denote the feasible non-dominated solutions at any particular sample generation, viz. the 51st, 101st, 301st, 501st and 1001 st generation. Although Fig.12 shows all the non-dominated solutions at any particular generation, only one or two distinct points (in the objective space) can be seen for that generation. However, a few distinct solutions in the design variable space may have the same objective function values and therefore, such solutions would coincide in the objective space. The number shown in parenthesis next to every point marker indicates the total number of such coincident solutions.

Discussion
For the test problem, the results are summarized in Section 4.2. The hybrid GA proves its efficiency by converging the two objective function values ( tan dis ce f and material f ) to the optimal values. The recurrent fluctuations in Fig. 10 and Fig. 11 show the effect of local search to the hybrid algorithm. Fig. 9 shows how the local search works.

Optimization of structures under load uncertainties
In this work, the structural optimization problem is defined in a design space shown in Fig.  13 using 4-noded quadrilateral elements. The present problem is similar to design problem solved in (Wang, Yang et al. 2008  where w is the weight of the structure and d is the displacement of the loading point 2. p is the vector of loads under uncertainty, that is 1 P in this problem. Local search which is triggered every 10 generations in this work is only applied to selected N C solutions. A constraint on the vertical displacement, d g , is used to prevent big deformations which are supposed likely to occur. A constraint on the maximum stress in the structure (to prevent fatigue or failure) is important. A dimensionless expression for the stress constraint may be written as  Fig. 14. Fig. 14 (a) shows the solution with best weight objective under the worst load case where 1 P is 55.5N. Fig. 14 (c) shows the solution with the best displacement objective under the worst load case where 1 P is 55.4N. One solution with median weight and displacement objective is given in Fig. 14 (b). Fig. 15 shows a plot of the best weight objective ( w ) and the corresponding solution's displacement objective ( d ) versus generation number. w and d values on the plot corresponding to any particular generation number belong to that generation's nondominated feasible solution which has the best weight objective.    Fig. 16 shows plots of the best displacement objective (d) and the corresponding solution's weight objective ( w ) versus generation number. As can be seen from Fig. 15 and Fig. 16, there are some fluctuations because of anti-optimization. Fig. 17 shows a plot in the objective space, where the solid shape markers denote the feasible non-dominated solutions at any particular sample generation.

Conclusion
The versatility and effectiveness of the topology optimization methodology developed in this work rest on three key components: an efficient morphological geometry representation that defines practical and valid structural geometries, a compatible graph-theoretic chromosome encoding and reproduction system that embodies topological and shape characteristics, and a multiobjective hybrid GA with local search strategy as the worst-casescenario technique of anti-optimization. The use of local search strategy helps to direct and focus the genetic search in uncertainty design variable space. A multiobjective target matching problem with known solutions has been formulated and solved to demonstrate the validity of presented algorithm. Simulation results of the structure optimization under load uncertainty are encouraging, which indicates that the hybrid algorithm integrates local search as anti-optimization is applicable. The proposed tournament constrained selection method works well and the computation cost is reasonable.