Multi-objective topology optimization using the Boundary Element Method

Multi-objective topology optimization using the Boundary Element Method

Structures 19 (2019) 84–95 Contents lists available at ScienceDirect Structures journal homepage: www.elsevier.com/locate/structures Multi-objectiv...

3MB Sizes 0 Downloads 117 Views

Structures 19 (2019) 84–95

Contents lists available at ScienceDirect

Structures journal homepage: www.elsevier.com/locate/structures

Multi-objective topology optimization using the Boundary Element Method a,⁎

b

c

Hélio Luiz Simonetti , Valério S. Almeida , Francisco de Assis das Neves , Marcelo Greco

T

d

a

Department of Mathematics, Federal Institute of Minas Gerais (IFMG) Department of Geotechnical and Structural Engineering from the School of Engineering of the University of São Paulo (EPUSP) Department of Civil Engineering, Federal University of Ouro Preto (UFOP) d Department of Structural Engineering, School of Engineering of Federal University of Minas Gerais (UFMG) b c

A R T I C LE I N FO

A B S T R A C T

Keywords: Multi-objective optimization Weighted sum Exponentially weighted criterion Boundary Element Method Evolutionary Structural Optimization

This article aims to explore the application of an evolutionary optimization technique for multi-objective optimization problems using as criteria the minimization of von Mises maximum stress and minimization of the maximum growth of internal structural strain energy. To evaluate the overall effect on the optimal design configuration, due to the removal of inefficient material from the structure by using these two optimization criteria, a goal weighting scheme is adopted, whereby the weight factors emphasize and balance the stress and strain energy criteria. Also considered in this study was the method of the exponentially weighted criterion for multi-objective optimization and the Pareto optimal concept. Thus, a contribution is made to the study of these two methods in the structural optimization procedure using a linear analysis by the Boundary Element Method. Four examples are presented to demonstrate the ability of the proposed method to solve structural design problems using multi-objective optimization.

1. Introduction The traditional methods of structural topology optimization, such as Solid Isotropic Microstructure with Penalty (SIMP), Bendsøe [1] and Rozvany et al. [2], the Homogenization-Based Optimization (HBO) method, Bendsøe and Kikuchi [3] and Ma et al. [4] and the Evolutionary Structural Optimization (ESO) method, Xie and Steven [5,6] are applicable to multi-objective optimization problems related to static or dynamic problems. Multi-objective optimization is the method of forming a solution that satisfies a series of conflicting objectives in the best possible way; that is, multi-objective optimization is the generation of a project that obtains the best performance of the structure when several criteria are taken into consideration. The solution to the multiobjective problem is known as a Pareto optimum. Pareto optima are not unique solutions, but rather a set of optimal solutions that simultaneously satisfy several optimization criteria. This article incorporates the multi-objective optimization problem in the ESO algorithm coupled to the Boundary Element Method (BEM) using the weighted objective technique and method of the exponentially weighted criterion. The ESO method developed by Xie and Steven [5] slowly removes the inefficient material from a structure and the residual shape of the structure evolves into a stationary optimum. Cervera and Trevelyan [7,8] used ESO-BEM coupling for the

structural optimization of 2D and 3D problems. The developed algorithm creates cavities in the structure's domain during the optimization process, based on the von Mises stress criterion. The concept of a topological derivative based on boundary elements was used by Cisilino [9] and Marczak [10] for the optimization of topology in potential problems. These two studies were all theoretically based on the topological derivative presented by Novotny et al. [11], who proposed a new approach to topological sensitivity analysis. Neches and Cisilino [12], proposed the optimization via BEM using a direct formulation. The model is discretized by linear elements and a distribution of internal points in the domain. The total potential energy is selected as the objective function. The evaluation of the topological derivative at internal points is performed as a post-processing procedure. Then, the internal points with lower values of the topological derivative are removed. The triangulation algorithm, capable of detecting cavities, is used to reconstruct the mesh and the new geometry. The procedure is repeated until a certain stop criterion is satisfied. The strategy proved to be flexible and robust. In all the methods discussed, an additional mechanism is always necessary to deal with topological changes; that is, the merging of cavities with each other or with the boundary during the optimization iterations. In addition, the evolving geometry and optimal topology exhibit uneven edges. In the last decade, several researchers have used BEM for Topology



Corresponding author. E-mail addresses: [email protected] (H.L. Simonetti), [email protected] (V.S. Almeida), [email protected] (F. de Assis das Neves), [email protected] (M. Greco). https://doi.org/10.1016/j.istruc.2018.12.002 Received 31 October 2018; Accepted 2 December 2018 Available online 11 December 2018 2352-0124/ © 2018 Published by Elsevier Ltd on behalf of Institution of Structural Engineers.

Structures 19 (2019) 84–95

H.L. Simonetti et al.

where u (q) and p (q) are the prescribed values and Γ = Γu ∪ Γp. The equations of the boundary elements are derived from Eq. (1) by introducing the boundary element into the discretization in u, p, we have the expression:

Optimization (TO) by coupling it to the Level Set Method (LSM) in acoustic problems Isakari et al. [13], heat conduction Jing et al. [14] and electromagnetic Isakari et al. [15]. This method also presents other relevant studies such as Abe et al. [16], Allaire and Jouve [17], Yamada et al. [18,19], Yamasaki et al. [20], Ullah and Trevelyan [21], Isakari et al. [13] and Ullah et al. [22]. However, these approaches treat BEM by coupling it to optimization processes using only one optimization criterion. An exception is the study of Isakari et al. [15], which presents a multi-objective optimization model that can find a material distribution by locally minimizing the intensity of transverse magnetism and transverse waves with multiple frequencies. The basic idea of the proposed method is to use objective weighting or the KreisselmeierSteinhause function (see [23]) of the functional original objectives as a new functional objective. Anflor and Marczak [36] coupled the Topological Derivative (DT) to BEM for simultaneous heat transfer and mass optimization problems using a multi-objective approach via a goalweighting technique. Toledo et al. [24] used an Evolutionary MultiObjective Optimization (EMOO) approach coupled with a Dual Formulation of BEM for shape optimization of a Y noise barrier. The analysis carried out in this study focuses on the simultaneous optimization of two conflicting objectives: a maximizing noise attenuation and minimizing the costs of all barrier material used, represented by the total length of the boundaries. In this article, a multi-objective topology optimization approach is proposed using the ESO-BEM coupling, based of combining a boundary representation and a criterion for creating cavities. As such, the optimization of shape is not performed, therefore, no boundary evolution occurs, that is, the boundary mesh is fixed throughout the optimization process and the optimal topology is obtained from the distribution of objective functions in the solution domain. In addition, cavity fusion is not performed since the representation of the boundary in the hexagonal shape to create the cavity facilitates the direct connectivity of the elements. Thus, the topology is improved for each cavity created in the domain without the need for special techniques like the methods already presented. Thus, to show the efficiency and effectiveness of the methodology used for the EMOO procedure, two methods were implemented: a) the weighted sum method and b) the exponentially weighted criterion method. In the two methodologies explored in this article, the von Mises stress and the strain energy were considered as objective functions. The present article is organized as follows: the basic details of the BEM formulation are introduced in Section 2, the formulation and the mono-objective optimization problem in von Mises stress and strain energy are developed in Section 3. In Section 4, the details of the cavity insertion algorithm in the domain and calculation of the weighted average of the von Mises stress are presented. In Section 5, the problem of multi-objective optimization is presented. The results obtained with the proposed optimization program with the present formulation are presented in Section 6.

[H ]{U } = [G]{P }

(3)

where {U} and {P} are, respectively, the vectors of the nodal values for u and p. The components of the matrices [H] and [G] are given by:

Hij = cij +

∫ p∗ (s, q) φj (q) dΓq Γj

Gij =

(4)

∫ p∗ (s, q) ϕj (q) dΓq (5)

Γj

with j = 1, 2, …, NE. For NE being the number of boundary elements, with Hij and Gij being the 2 × 2 sub-matrices of [H] and [G], s is the set point and ϕj is the j interpolation function of jth of Γj and cij is the submatrix derived from free term. The unknown terms of {U} and {P} are obtained by solving the linear system:

[A]{X } = {B}

(6)

where {X} is a vector given by the unknown sub-vectors of {U} in Γp and {P} in Γu, [A] is a matrix composed of sub-matrices of [H] and [G], and {B} is a vector composed of the prescribed values of {U} and {P} of Eq. (1). 3. Formulation and the problem of mono-objective optimization in von Mises stresses and strain energy 3.1. von Mises stress formulation In a BEM analysis, stresses (or any other required constraint) within the project domain are calculated at internal points. In this article, the von Mises stress for the plane state stress is given by the expression:

σ vm =

2 2 σ 11 +σ 222−σ11 σ22 + 3σ12

(7)

Removal of material from the structure occurs when von Mises stress at an internal point obeys the following inequality

σivm < RR∗ σivm ,max

(8)

with i = 1, 2, …, NPI, where NPI represents the number of internal points, RR is the removal ratio, σivm is the von Mises stress at the internal point i, and σi, maxvm is the maximum von Mises stress of the structure. The initial RR removal rate is prescribed and ranges from 0.01 to 0.05. This ratio increases incrementally with an evolutionary ratio (ER) as the optimization procedure progresses and the update is performed by equation

RRk + 1 = RRk + ER

(9)

2. The Boundary Element Method 3.2. Strain energy formulation Consider a two-dimensional linear elasticity problem. A homogeneous and isotropic material is assumed. Thus, the integral boundary equation for this problem is given by:

cij (s ) ui (s ) = −

∫ pij∗ (s, q) uj (q) dΓ + ∫ uij∗ (s, q) pj (q) dΓ Γ

Γ

Considering the constitutive relationships between stress and strain, such that the strain is a homogeneous function of degree n of stress, then, J*, Eq. (10), must be a homogeneous function of the stress components of degree n + 1.

(1)

J∗ =

where Γ is the boundary of the structural domain Ω, u and p are displacement and traction, u* and p* are the fundamental solutions and c is the free term. Integration in Γ is achieved in relation to q. Thus, the boundary conditions for u(q) and p(q) are given by:

u (q) = u (q) in Γu p (q) = p (q) in Γp

n σij ε ij n+1

(10)

Thus, from Euler's theorem for homogeneous functions, When the stress is proportional to the strain (n = 1) obtains J given by:

J= (2)

1 σij ε ij 2 Knowing that the equation εij can be written in the form:

85

(11)

Structures 19 (2019) 84–95

H.L. Simonetti et al.

εij =

(1 + ν ) ν σij − δij σkk E E

(12)

where the E and ν are Modulus of Elasticity and Poisson's ratio. Using indicial notation and mathematically manipulating Eqs. (11) and (12), we obtain the strain energy J as a function of the maximum stresses given by Eq. (13).

J=

1 ν (1 + ν ) 2 2 2 (σ 11 ) + σ11 σ22 + σ12 +σ22 2E E E

(13)

3.3. The problem of mono-objective optimization in von Mises stress and strain energy In a macrostructural approach, the topology of the structure is modified by the insertion of cavities in the problem domain. ESO can be cited as an example of this TO group, which is based on the solution of the objective function when a domain region is removed. Tanskanen [25] mathematically showed that the ESO material removal heuristic is associated with regions with low strain energy levels. It is a non-explicit minimization procedure in relation to its weight when compared to a deterministic method, such as the algorithm of optimization technique called Sequential Linear Programming (SLP). Furthermore, according to Querin [26], ESO is a subset of the Integer Linear programming method. The main difference is that the inequality constraint functions are not limited by a finite number. Instead, a linear or nonlinear evolutionary function is used. Also, the inequality constraints change or evolve from the value of an initial value zero to the end of the value defined. Thus, the generalized mathematical representation for the optimization problem through the ESO-BEM method with the von Mises stress criterion can be defined as:

Fig. 1. Creation of the hexagonal cavity.

cavities are created around the internal points as shown in Fig. 1. In Fig. 1, it is seen that DY is the distance between the internal point IP5 and element 14, DX is the length of the element 14. Considering that the internal point IP5 satisfies the inequality 10, it is necessary to create the hexagonal cavity around the point IP5, Fig. 3. In this way, elements 11 and 14 have length L = DX and the other elements are calculated using the expression (16):

minimize F (W , Xi ) subject to [H]{U} − [G]{P} = 0 {U} = {U} in Γu {P} = {P} in Γp

L=

d

Ω

(14)

where F is the objective function, V is the volume of the structure and Vmax is the final volume prescribed. The shape of the domain Ω is represented by the design variable Xi and W represents the variation of the vectors {U} and {P}. In this sense, the generalized mathematical representation for the optimization problem through the ESO-BEM method with the strain energy criterion can be defined as:

4.2. Weighted stress average

minimize F (W , Xi ) subject to [H]{U} − [G]{P} = 0 {U} = {U} in Γu {P} = {P} in Γp

With the insertion of the cavity and consequently, the removal of material from the domain, the stress at points closer to this cavity increases. In order to smooth the stresses in regions close to these cavities, a weighted average was created, see Fig. 3. Thus, the removal of material from the domain is performed in such a way that the stress at the internal point (i), is given by (j) within a circle of radius (R), whose weights are the distances of element (i) to elements (j) and are compared with the stresses at the interior points. Mathematically, the mean of the stress is given by Eq. (17):

∫ dΩ − Vmax ≤ 0 Ω

Ji − RR∗Jimax ≤ 0 ∀ i = 1, 2, …, NPI

(16)

Thus, when 0 < DY < 2 , the adjacent cavities no longer have elements in common, and no regional coupling is required. If d DX = DY = 2 , we have the common nodes and elements in the domain region of the structure. In this case, it is necessary to connect the nodes and elements during the optimization process. Considering 3 DX = DY ⋅ 2 the hexagonal cavity will be regular. The flowchart with cavity creation process is shown below, see Fig. 2.

∫ dΩ − Vmax ≤ 0 vm σivm − RR∗σi,max ≤ 0 ∀ i = 1, 2, …, NPI

DX 2 ⎛ ⎞ + (DY )2 ⎝ 2 ⎠

(15)

4. Cavity creation and average weighted stress

N

∑ wj σjvm 4.1. Cavity creation

σivm =

j=1 N

∑ wj

The cavities are created around the inner points that satisfy the inequality 8. Thus, the material is removed by inserting a hexagonal cavity whose edge are boundary linear elements. Once a steady state is reached the iterative procedure is updated by Eq. (9) which controls the step size through the evolutionary ER ratio. In the present formulation, it was verified that the optimization procedure via ESO-BEM is not sensitive to these parameters as in OT applications using FEM. The

j=1

(17)

where σivm is the mean of the von Mises stress at the internal point i (center of the circle), and σjvm is the von Mises stress at the internal point j. N is the number of internal points j within the circle of radius R and wj is the inverse of the distance between the internal points i and j given by Eq. (18). 86

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 2. Cavity creation flowchart.

wj =

1 (x j − x i )2 + (yj − yi )2

(18)

Thus, the removal criterion is given by the inequality 19 and is able to determine the influence of the stress of the elements j around the central point i when a cavity is created by sensitively softening the stress.

σivm < RRk ⋅(σivm ,max )

(19)

The smooth closed polygonal with line segments of uniform length used to create the hexagonal cavities are necessarily used in the structural analysis in the next iteration. In the analysis of contour elements, if intersection points are used directly as node points of the element, two points of intersection may be very close to each other, which could cause difficulties and instabilities during the analysis. However, all the formulation presented uses the double-node technique which allows element discontinuity, and the smoothness of the polygonal geometry used does not produce high stress concentrations. For this reason, no additional technique has been used in the present formulation to treat the stresses, other than the spatial filter which in itself has the ability to do so.

Fig. 3. Weighted average of the stresses at the central point i.

87

Structures 19 (2019) 84–95

H.L. Simonetti et al.

5. Multi-objective optimization problem 5.1. Definition In contrast to single goal optimization, a solution to a multi-objective problem is more a concept than a definition. Usually, there is not a single global solution, and it is often necessary to determine a set of points that fit into a predetermined definition for an optimal one. In general, an MO problem can be expressed as follows:

Minimize: F(x) = [F1 (x ), F2 (x ), F3 (x ),…, Fk (x )]T x Subject to hi (x ) = 0, i = 1, 2, …, n g j (x ) ≤ 0, j = 1, 2, …, m xlj < x t < x uj , t = 1, 2, …, p

(20)

where k is the number of objective functions and m is the number of inequality constraints, n is the number of equality constraints, p is the number of project variables, x ∈ RN is the vector for project variables and F(x) ∈ Rk is the vector of the objective functions, where Fi(x) : RN → R. The viable space for the project is defined as X = {x/ gj(x) ≤ 0, j = 1, 2, …m} and the viable space criterion is defined as F (x) = {F(x)/x ∈ X}. In Eq. (20), the notation suggests a goal of simultaneously minimizing more than one objective function, which is impossible. One way to determine a solution that partially satisfies the equations present in MO is contained in the definition of Pareto Optimality, Pareto [37], which is defined as follows:

Fig. 5. Design domain.

mentioned: Hierarchical Optimization Method, Trade-Off Method, Goal Programming Method, Global Criteria Method (Global Criterion Method) and Weighting Objectives Method. In addition, there are heuristic techniques that are divided into two types: 1) constructive heuristics - indicated to solve multi-objective combinatorial optimization problems to quickly obtain the Pareto-optimal set; 2) refinement heuristics, also called local search methods, are optimization methods that consist of the refinement of a solution. Among the refinement heuristics, we can highlight the Variable Neighborhood Descent Method (VND), which is based on searches in different neighborhood structures. The present study addresses only weighted summing methods for the formulation of the Pareto concept to obtain optimal topologies. Fig. 4 shows the optimal Pareto front for two objective functions. Since the definition of a weight vector leads to only one point on the Pareto curve, performing several optimizations with different weight values can produce a considerable computational load. Therefore, it is necessary to choose which combinations of different weights should be considered to reproduce a representative part of the Pareto front. In addition to this computational cost, the objective weighting method has two technical shortcomings, explained as: a) The relationship between the weights of the objective function and the Pareto curve is such that a uniform distribution of the weight parameters, in general, does not produce a uniform distribution of points on the Pareto curve, and (b) the non-convex parts of the Pareto set cannot be achieved by minimizing the convex combinations of the objective functions. Despite the shortcomings with respect to representing the optimal Pareto set, the weighted sum method for MO problems continues to be widely used, not only to provide multiple solution points by consistently varying the weights, but also to provide a single point solution that reflects the preferences presumably incorporated in the selection of a single set of weights, Marler and Arora [27]. In the goal-weighting method, the multi-objective problem is converted into a mono-objective problem using the weighted sum of the original multiple objectives. That is, the von Mises stress, Eq. (7) and the strain energy, Eq. (13), are combined to form a single new criterion, given by:

Definition 1. Pareto optimality - A point, x∗ ∈ X is a Pareto optimum if there is no other point x ∈ X, such that F(x) ≤ F(x∗) and Fi(x) ≤ Fi(x∗) for at least one function. Where i is the number of functions to be optimized as exemplified in Fig. 4. Obviously the Definition 1 needs a strict inequality somewhere to be meaningful. All Pareto optimal points lie within the bounds of the feasible criterion space. Often, algorithms provide solutions that may not be optimal for Pareto, but may meet other criteria, making them meaningful for many practical applications. For example, the Pareto weak optimality is defined as follows: Definition 2. Weak Pareto optimality - A point x∗ ∈ X is a weak Pareto optimal if there is no other point x ∈ X, such that F(x) ≤ F(x∗). There are several strategies to solve MO problems, among which can

2

Fmult = w1∗F1 + w2∗F2 = e

∑ wi Fi i=1

(21)

where Femult is the multi-objective function that will be used to define the new element removal criterion, wi is the ith weighting factor with 2

0 ≤ wi ≤ 1 and ∑ wi = 1 with F1 and F2 given by:

F1 =

{σivm }

i=1

NPI − 3

⎛⎜ ∑ σ vm ⎞⎟/3 i,max ⎝ i = NPI ⎠

Fig. 4. Optimal Pareto front for the weighted sum method. 88

(22)

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 6. Optimal topologies: (a) ESO-BEM present formulation, (b) FE-ESO Chu et al. (C) Cervera and Trevelyan [7,8] and (d) boundary-LSM Ullah and Trevelyan [29].

Fig. 7. (a) Pareto front and optimal settings, (b) influence of weights on functions.

F2 =

{Jimax } NPI − 3 ⎛⎜ ∑ J max ⎞⎟/3 i ⎝ i = NPI ⎠

objectives, favoring each function according to the weight. Thus, by defining the von Mises stress and the strain energy as the control parameter of linear elastic analysis, a series of internal points with the lowest values for Femult are removed so that the maximum increase in the von Mises stress and the maximum strain energy growth are minimum.

(23)

NPI − 3

vm ⎞ where ⎛⎜ ∑ σi,max ⎟ /3 is the average of the three largest von Mises ⎝ i = NPI ⎠ NPI − 3 stresses, ⎜⎛ ∑ Jimax ⎟⎞/3, is the average of the three largest strain en⎝ i = NPI ⎠ ergies and NPI is the number of internal points of the stress and energy vectors. Since the control is performed with the von Mises (normalized) stresses and the strain energy (normalized), a series of internal points with the smallest values of Femult will be removed so that both characteristics can be improved. In this way, Femult creates a cavity assuming that this point is “inefficient” for the structure. Ideally, the Femult in each part of the structure should be close to a safe level. Thus, by gradually removing material with lower values of Femult, the new structure becomes increasingly uniform and obeys the weight of the

5.2. Exponentially weighted criterion The method of the exponentially weighted criterion proposed by Athan and Papalambros [28] arises because of the inability of the weighted sum method to capture points in nonconvex portions of the Pareto optimal surface. Thus, the method is described as follows: k

Femult =

∑ (e p ∗wi − 1) i=1

∗ e p ∗ Fi (x ) (24)

where k represents the number of objective functions, herein k = 2, wi 89

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 8. Distribution of the objective function in the solution domain.

Fig. 9. (a) Design domain, (b) domain considering symmetry.

Fig. 10. (a) Pareto front and optimal settings, (b) influence of weights on functions.

90

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 6c and Ullah and Trevelyan [29] using Level Set Methods (LSM), Fig. 6d. The optimal design obtained with the present approach shows a good accuracy with the mentioned approaches. Fig. 7a shows the Pareto boundary with the optimal configurations for the weights of w1 = 0, w1 = 0.1, w1 = 0.3, w1 = 0.5, w1 = 0.7, w1 = 0.9, and w1 = 1.0. The influence of weights on the functions F1 (weighted von Mises stresses), F2 (weighted strain energy) and Femult = w1∗F1 + w2∗F2 (weighted sum) is shown in Fig. 7b. Fig. 8 shows the surface of the objective function in the solution domain. It is emphasized that the surface is created with the Pareto points for a w1 = 0.7 weight in the optimal iteration. 6.2. Example 2 - MBB beam - weighted sum The so-called MBB beam, Messerschmitt-Bolkow-Blohm GmbH; Payten et al. [31]; Bulman et al. [32], is a classic example of topological optimization. It consists of a simply supported beam, loaded with a vertical force centered on its upper edge, see Fig. 9a. The initial parameters were as follows: Young's modulus E = 210 MPa, Poisson's coefficient ν = 0.30, RR = 1.0% and ER = 1.5%, length 20 L and height 5 L. A load of 100 N was applied in the upper central part of the beam, divided into four linear elements of 10 cm. Due to symmetry conditions, only half of the domain is discretized, see Fig. 9b, with 325 internal points and 60 linear boundary elements. Fig. 10a presents the Pareto frontier with the optimal configurations for the weights of w1 = 0, w1 = 0.2, w1 = 0.4, w1 = 0.6, w1 = 0.8 and w1 = 1.0. The influence of the weights on the functions F1 and F2 is represented in Fig. 10b. Fig. 11 illustrates the evolutionary procedure of the objective function for a weight w1 = 0.6 with the optimum configuration obtained for a prescribed volume of 33.5% of the design volume. It is noted that the removal of inefficient elements of the mesh provides a reduction in the goal weighting function; that is, the weighted sum is minimized by the reduction of material until the designated volume is reached.

Fig. 11. Volume variation in relation to the objective function.

6.3. The L-shaped beam - weighted sum

Fig. 12. Design domain (in dm).

The L-shaped structure is another reference example in the field of topology optimization. The model is bi-supported on the upper edge and a 100 N load is applied in the middle of the right edge as shown in Fig. 12. The properties of the material used in these examples are: Poisson's coefficient 0.30, Young's modulus 210 GPa and thickness of 1 mm. The optimization parameters are: RR = 1% and ER = 1% and the optimum design configuration occurs when the prescribed volume achieves a reduction of 55% of the initial volume. The cavities inserted in the domain are all equal with an area equal to 0.45 dm2 with DX = 0.40 dm and DY = 0.375 dm. To validate the capability of the proposed multi-objective optimization method to deal with stress peaks, the von Mises stress distribution in different iterations is shown in Figs. 13 and 14 for weighting factors equal to w1 = 0.7 and w1 = 0.3 with 621 internal points in the structural domain and 100 linear elements. The comparison of the stress distribution results in the solution domain shows that the proposed optimization method allows the stress peaks observed in iteration 1 to be distributed on a larger surface, smoothing them during the iterative process. This results in an optimal “Full Stress Design” (FSD) design. In addition, calculating the weighted average within a circle of radius (r) increases the convergence to an optimal one in a smoother way. The optimal design generated with the weight factor w1 = 0.7 is similar to those available in the literature of Xia et al. [33], Allaire et al. [34] and Allaire and Jouve [38]. It is possible to observe in the presentation of the topologies shown in Fig. 14 that the geometry of L-shaped with w1 = 0.3, that is to say, favoring the strain energy in 70%, has a greater smoothing in the region where stress peaks than that observed in Fig. 16, where the weight

represents the ith weighting factor, Fi(x) is the ith objective function, and p is an integer value. It should be noted that the sum argument represents a function of individual utility for Fi(x). Although large values of p can lead to numerical excess, minimizing Eq. (24) provides a necessary and sufficient condition for Pareto optimization. 6. Results 6.1. Example 1 – short cantilever beam - weighted sum A common reference problem is the short cantilever beam with a central load and a ratio of 1.6. The definition of the initial problem is shown in Fig. 5. This structure having dimensions 160 mm × 100 mm is fixed in the upper and lower corners of the left side, while being subjected to a vertical load at the center of the free edge. The load was arbitrarily chosen to have a value of 100 N. This value can be considered very small in comparison to that of failure; however, it produces a stress distribution that is suitable to study the removal of material, such as the creation and evolution of cavities. The optimization parameters are RR = 5.0% and ER = 4.90%. The properties of the material are Young's Modulus 210 GPa, and Poisson's coefficient equal to 0.30. The boundary is initially discretized with 52 linear elements and a domain with 336 internal points. Numerical integrals are evaluated with 12 Gaussian points. Fig. 6a compares the optimal design obtained with the present multi-objective ESO-BEM approach with weight w = 0.7, with the solutions obtained by Chu et al. [30] using FE-ESO, Fig. 6b, Cervera and Trevelyan [7,8] using the boundary-ESO solution, 91

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 13. Topology history for w1 = 0.7 - (a) iteration 1, (b) iteration 8 (c) iteration 11 and (d) iteration 14.

coefficient ν = 0.30. The optimization parameters used were: RR = 1.75% and ER = 2.175%, whereby the optimum design configuration occurs when the prescribed volume reaches a reduction of 71% of the initial volume. The cavities inserted in the domain are all equal, having an area equal to 0.1875 dm2 with DX = 0.25 dm and DY = 0.25 dm. The design domain and boundary conditions are represented in Fig. 16a. The optimal configuration proposed with the present formulation can be compared with that presented by Cervera [35], see Fig. 16b and c. Fig. 16d, uses the ESO method coupled with BEM to perform topology optimization. The proposed algorithm creates the internal cavities based on NURBS (Non Uniform Rational Basis Spline) curves to perform topology changes. The boundary mesh was discretized with 59 linear elements and a domain of 171 internal points. The graph of Fig. 17 illustrates the volume reduction during the optimization process and the evolution of the objective function for the exponentially weighted criterion.

w1 = 0.7 favors the von Mises stresses by 70%. To show the relationship between normalized objective functions and material removal, a history of the evolution of optimization for weight w1 = 0.7 is plotted against the volume ratio of the structure, see Fig. 15a and b. They show that the objective functions are conflictive, since the removal of material provides an increase in stress energy and a decrease of von Mises maximum stress. Based on the stress minimization criterion, the von Mises maximum stress of the structure is noticeably reduced as the material is removed. In this sense, the inverse process is verified when the strain energy criterion is evaluated; that is, the maximum strain energy grows when the material is removed. A minimum value of the maximum stress and the minimum growth of the strain energy is the goal to be achieved when the desired volume is 39.5% of the initial volume. 6.4. Two bars frame - exponentially weighted criterion

7. Conclusions

The exponentially weighted criterion was applied to optimize the problem of two bars frame in which a weight of w = 0.8 was used as a weighting factor. In this example, the value of p = 4. The properties of the material used are: modulus of elasticity E = 210 GPa, and Poisson

The coupling of the ESO optimization method with BEM was developed to incorporate multi-criteria optimization using the goal 92

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 14. Topology history for w1 = 0.3 - (a) iteration 1, (b) iteration 8 (c) iteration 11 and (d) iteration 14.

Fig. 15. (a) Evolution history for strain energy and (b) evolution history for von Mises stress.

weighting method and the exponentially weighted criterion method. It was verified that the topologies produced by the weighting method form the Pareto solution space in the analyzed model. These designs are

based exclusively on the minimization of the von Mises maximum stress and on the minimization of internal strain energy growth and are optimal for a prescribed volume fraction. Any improvement in one 93

Structures 19 (2019) 84–95

H.L. Simonetti et al.

Fig. 16. (a) Design domain, (b) optimal topology with hexagonal cavities, (c) optimal topology with the distribution of the function in the solution domain and (d) optimal topology by Cervera [35].

for Scientific and Technological Development). References [1] Bendsøe MP. Optimal shape design as a material distribution problem. Struct Optim 1989;1:193–202. [2] Rozvany GIN, Zhou M, Birker T. Generalized shape optimization without homogenization. Struct Multidiscip Optim 1992;4(3):250–2. [3] Bendsøe MP, Kikuchi N. Generating optimal topologies in structural design using a homogenization method. Comput Methods Appl Mech Eng 1988;71:197–224. [4] Ma ZD, Kikuchi N, Chang H-C. Topological design for vibrating structures. Comput Methods Appl Mech Eng 1995;121(1–4):259–80. [5] Xie YM, Steven GP. A simple evolutionary procedure for structural optimization. Comput Struct 1993;49(5):885–96. [6] Xie YM, Steven GP. Evolutionary Structural Optimization (ESO). Springer; 1997. [7] Cervera E, Trevelyan J. Evolutionary structural optimisation based on boundary representation of NURBS. Part I: 2D algorithms. Comput Struct 2005;83(23):1902–16. [8] Cervera E, Trevelyan J. Evolutionary structural optimisation based on boundary representation of NURBS. Part II: 3D algorithms. Comput Struct 2005;83(23):1917–29. [9] Cisilino AP. Topology optimization of 2D potential problems using boundary elements. Comput Model Eng Sci 2006;15(2):99–106. [10] Marczak RJ. Optimization of elastic structures using boundary elements and a topological-shape sensitivity formulation. Brazilian Society of Mechanical Sciences and Engineering978-85-85769-30-7; 2007. [11] Novotny AA, Feijoo RA, Taroco E, Padra CC. Topological sensitivity analysis. Comput Methods Appl Mech Eng 2003;192:803–29. [12] Neches LC, Cisilino AP. Topology optimization of 2D elastic structures using boundary elements. Eng Anal Bound Elem 2008;32(7):533–44. [13] Isakari H, Kuriyama K, Harada S, Yamada T, Takahashi T, Mat-sumoto T. A topology optimisation for three-dimensional acoustics with the level set method and the fast multipole boundary element method. Mech Eng J 2014;1(4):CM0039. [14] Jing G, Matsumoto T, Takahashi T, Isakari H, Yamada T. Topology optimization for 2D heat conduction problems using boundary element method and level set method. Trans JASCOME 2013;13(19):5–10. [15] Isakari H, Nakamoto K, Kitabayashi T, Takahashi T, Matsumoto T. A multi-objective topology optimisation for 2D electro-magnetic wave problems with the level set method and BEM (to appear in). Eur J Comput Mech 2016;25(1-2):165–93. [16] Abe K, Kazana S, Koro K. A boundary element approach for topology optimization problem using the level set method. Commun Numer Methods Eng 2007;23:405–16. [17] Allaire G, Jouve F. Minimum stress optimal design with the level set method. Eng Anal Bound Elem 2008;32(11):909–18. [18] Yamada T, Izui K, Nishiwaki S, Takezawa A. A topology optimization method based on the level set method incorporating a fictitious interface energy. Comput Methods Appl Mech Eng 2010;199(45–48):2876–91. [19] Yamada T, Shichi S, Matsumoto T, Takahashi T, Isakari H. A level set-based topology optimization method using 3D BEM. Advances in boundary element techniques. XIV. 2013. p. 432–7. [20] Yamasaki S, Nomura T, Kawamoto A, Sato K, Nishiwaki S. A levelset-based topology optimization method targeting metallic waveguide design problems. Int J Numer Methods Eng 2011;87(9):844–68.

Fig. 17. Objective function by iteration.

criterion requires clear compensation with the other criterion. The solution of the exponentially weighted criterion method is also part of the Pareto solution space. Although only two optimality criteria have been considered herein, the procedure can be applied to other criteria such as: stress and displacement, stress and frequency, displacement and frequency or even more than two criteria simultaneously. This is being investigated and should be reported in the near future. It is also worth noting that the technique presented is quite simple, presenting a low computational cost because it requires a small number of internal points to capture the optimal topology desired. The results obtained are as accurate as those presented in literature.

Acknowledgements The authors are grateful to the Federal Institute of Science and Technology of Minas Gerais, their Department of Structural Engineering and Foundations, FAPESP (Foundation for Research Support of the State of São Paulo), the Federal University of Ouro Preto (UFOP) and FAPEMIG (Foundation for Research Support of Minas Gerais State) for their financial supports to this research, the CAPES (Coordination of Improvement of Higher Level Personnel) and CNPq (National Council 94

Structures 19 (2019) 84–95

H.L. Simonetti et al.

[21] Ullah B, Trevelyan J. Correlation between hole insertion criteria in a boundary element and level set based topology optimization method. Eng Anal Bound Elem 2013;37(11):1457–70. [22] Ullah B, Trevelyan J, Ivrissimtzis I. A three-dimensional implementation of the boundary element and level set based structural optimisation. Eng Anal Bound Elem 2015;58:176–94. [23] Kreisselmeier G, Steinhauser R. Systematic control design by optimizing a vector performance index. 1979. p. 113–7. [24] Toledo R, Aznárez JJ, Maeso O, Greiner D. A methodology for the multi-objective shape optimization of thin noise barriers. 2016. [25] Tanskanen P. The evolutionary structural optimization method: theoretical aspects. Comput Methods Appl Mech Eng 2002;191:5485–98. [26] Querin OM. Evolutionary structural optimization stress based formulation and implementation (PhD dissertation). University of Sydney; 1997. [27] Marler RT, Arora JS. The weighted sum method for multi-objective optimization: new insights. Struct Multidiscip Optim 2010;41(6):853–62. [28] Athan TW, Papalambros PY. A note on weighted criteria methods for compromise solutions in multi-objective optimization. Eng Optim 1996;27:155–76. [29] Ullah B, Trevelyan J. A boundary element and level set based topology optimisation using sensitivity analysis. Eng Anal Bound Elem 2016;70:80–98. [30] Chu DN, Xie YM, Hira A, Steven GP. On various aspects of evolutionary structural

[31] [32]

[33] [34] [35] [36]

[37]

[38]

95

optimization for problems with stiffness constraints. Finite Elem Anal Des 1997;24(4):197–212. Payten WM, Ben-Nissan B, Mercer DJ. Optimal topology design using a global selforganizational approach. Int J Solids Struct 1998;35(3–4):219–37. Bulman S, Sienz J, Hinton E. Comparisons between algorithms for structural topology optimization using a series of benchmarks studies. Comput Struct 2001;79:1203–18. Xia Q, Shi T, Liu S, Wang M. A level set solution to the stress based structural shape and topology optimization. Comput Struct 2012;90–91:55–64. Allaire G, Jouve F, Toader A. Structural optimization using sensitivity analysis and a level-set method. J Comput Phys 2004;194(1):363–93. Cervera E. Evolutionary structural optimisation based on boundary representation of B-spline geometry (PhD Thesis). Durham, UK: University of Durham; 2003. Anflor CTM, Marczak R. Multicriteria Topology Optimization of Heat And Mass Transfer Problems Using Boundary Elements. COBEM, 2007, Brasília. 19th Internatinal Congress of Mechanical Engineering. 2007. 2007. Brasília. Pareto V. Manuale di Economica Politica, Societa Editrice Libraria. In: A.S.Schwierand AN, editors. Milan; translated into English by A.S. Schwier as Manual of Political Economy. New York: A. M. Kelley; 1906. p. 1971. Allaire G, Jouve F. Minimum stress optimal design with the level set method. Engineering analysis with boundary elements. 32.11. 2008. p. 909–18.