Structural Topology Optimization Based on the Smoothed Finite Element Method

In this paper, the smoothed finite element method, incorporated with the level set method, is employed to carry out the topology optimization of continuum structures. The structural compliance is minimized subject to a constraint on the weight of material used. The cell-based smoothed finite element method is employed to improve the accuracy and stability of the standard finite element method. Several numerical examples are presented to prove the validity and utility of the proposed method. The obtained results are compared with those obtained by several standard finite elementbased examples in order to access the applicability and effectiveness of the proposed method. The common numerical instabilities of the structural topology optimization problems such as checkerboard pattern and mesh dependency are studied in the examples.


Latin American Journal of Solids and
To date, the predominant numerical method used for topology optimization is the finite element method.The finite element method encounters some difficulties when dealing with problems such as large deformation or moving boundary problems.The standard finite element method often suffers from numerical instabilities, and its solving is sensitive to element distortion because of overestimation of the stiffness matrix.To overcome these difficulties, various numerical methods have been developed and achieved remarkable progress, such as meshfree methods (Monaghan (1992); Belytschko et al. (1994); Liu et al. (1995); Atluri and Zhu (1998)) and smoothed finite element method (Liu et al. (2007)).The meshfree methods do not require maintaining the integrity and desired shape of elements due to their meshfree nature.Therefore, large deformation and crack propagation problems can be effectively modelled with meshfree methods (Shobeiri (2015a(Shobeiri ( , 2015b))).Though meshfree methods generally exhibit good numerical stability and accuracy, the complex field approximation considerably increases the computational cost.
It is clear that methods which combine finite element method with meshfree methods can exhibit advantages of computational efficiency and simplicity.The smoothed finite element method is such a typical method.This new numerical method is rooted in meshless stabilized conforming nodal integration and exhibits a number of attractive properties such as good numerical stability and accuracy, excellent convergence rate, and insensitivity to volumetric locking and mesh distortion (Liu et al. (2007)).The smoothed finite element method has been successfully applied to large variety of problems including 2D and 3D linear and nonlinear problems (Nguyen et al. (2009)), dynamic analysis (Luong- Van et al. (2014)), plate and shell structures (Nguyen-Xuan et al. (2008); Nguyen-Thanh et al. ( 2008)).
In this paper, the smoothed finite element method is proposed to carry out the topology optimization of continuum structures using the level set method.The feasibility and efficiency of the proposed method are illustrated with several 2D examples that are widely used in topology optimization problems.The optimized topologies are compared with those obtained by the standard finite element -based method in order to access the applicability and effectiveness of the proposed method.The common numerical instabilities of the structural topology optimization problems such as mesh dependency and checkerboard patterns are studied in the examples.

REVIEW OF CELL-BASED SMOOTHED FINITE ELEMENT METHOD (CS-FEM)
In this section, a brief review of the cell-based smoothed finite element method (as a branch of the smoothed finite element method) is presented.Full details can be found in Liu et al. (2007).In the cell-based smoothed finite element method, the total design domain W is first divided into e N ele- ments as in the finite element method.Depending on the necessity of stability, each element is then subdivided into a number of smoothing domains such that Here, ( or )   i j W is the domain of th i or th j smoothing domain and nc is the total number of cells inside the design domain.Fig. 1 shows the smoothing domains relating to various number of cells in the cell-based smoothed finite element method.For each smoothing domain ( ) c W associated with cell c, the smoothing strain S  can be written as: Latin American Journal of Solids and Structures 13 (2016) 378-390  where F is the smoothing function written as: ( ) where nSC is the number of cells for each element.Note that this kind of smoothing is also employed in the smoothed particle hydrodynamics method (Monaghan (1992)).Substituting F into Eq.(2.1), the smoothing strain can be written as: (2.5) Latin American Journal of Solids and Structures 13 (2016) 378-390 In the above equations, nb is the total number of boundary sections of G .It can be pointed out from Eq. (2.6) that unlike the finite element method, there is no derivative of shape function in the smoothing strain matrix.By employing the current formulation of the smoothed finite element method, the discrete equation for the cell-based smoothed finite element method is given as: where  U is the vector of nodal displacements,  f is the vector of nodal forces, and CS-FEM  K is the global smoothed stiffness matrix: where D is the stress-strain relationship matrix.It should be noted that nodal shape functions constructed in cell-based smoothed finite element method have the Delta function property.Therefore, essential boundary conditions are imposed as in the finite element method.

Optimization Problem
The objective of topology optimization is the optimal distribution of material in the design domain to minimize the cost functionals under various constraints such as stress, displacement or weight constraints.In this study, the aim is to find the stiffest structure subject to a given structural weight.Therefore, the optimization problem can be formulated as: where e c is the center position of element e .Here,  is initialized as a signed distance function and an upwind finite difference scheme is employed to accurately solve the evolution equation.The discrete level set function is updated to find new structure using the following equation: where t indicates time, u is scalar field on the design domain, w is a positive parameter that specifies the influence of g , and g is a forcing term that determines the nucleation of new holes within the structure.Note that to satisfy the weight constraint, two parameters u and g are obtained based on the shape and topological sensitivities of the Lagrangian framework.

NUMERICAL EXAMPLES
In

Example 1
Fig. 2 shows the design domain of a cantilever beam with a length to height ratio of 2:1.The objective function is to minimize the compliance and the objective weight is 45% of the total weight of the design domain.To determine the optimal structural layout, the design domain is discretized using 36 18, 48 24, 60 30 and 72 36 quadrilateral elements.And to study the effects of the number of the smoothing cells, each element is subdivided into different number of smoothing cells using 2 nSC = , 4 nSC = , 8 nSC = and 16 nSC = . The optimization results obtained from different mesh sizes with different number of smoothing cells are shown in Fig. 3, from which it can be seen that the final solutions obtained from 4 nSC = , 8 nSC = and 16 nSC = with different mesh sizes are almost identical, and are different from those obtained from 2 nSC = .The optimization results obtained using 2 nSC = show the so-called mesh dependency effect, for which different optimal topologies are generated from different mesh sizes.It can be found that the use of four ( 4 nSC = ) or more than four smoothing domains can be good choices for topology optimization problems to overcome numerical instabilities such as mesh dependency phenomenon.
Figs. 4 and 5 illustrate the history of objective function and objective weight using 4 nSC = based on different mesh sizes over iterations, respectively.It can be seen from Figs. 4 and 5 that the number of iterations for mesh sizes of 36 18, 48 24, 60 30 and 72 36 are 65, 57, 53 and 47, respec- tively and their corresponding compliances are calculated as 70.71, 71.35, 68.60 and 68.96, respectively.It can be also seen from these results that for different mesh sizes, their convergence characteristics are very similar.To verify the present method, the above problem using the same mesh sizes is solved by the Solid Isotropic Microstructure with Penalization (SIMP) method (Sigmund ( 2001)) (standard finite element-based method).The optimal structural layouts without and with sensitivity filtering are shown in Fig. 6, from which it can be seen that with and without using a filter, the topologies obtained from the SIMP method are quite different.Note that the SIMP method using a filter generates similar topologies to the designs obtained by the present method using 4 nSC = , 8 nSC = and 16 nSC = .It can be also observed that numerical instabilities such as checkerboard patterns and mesh-dependency exist in the results of the SIMP method if no filtering is employed, while for the present method no such problem can be seen.

Example 2
Fig. 7 shows the design domain of a beam with fixed supports.The beam length to height size ratio is 2:1.The objective function is to minimize the compliance, and the target weight is 30% of the total weight of the design domain.The division of the element into four smoothing domains ( 4 nSC = ) is used as default in this example.In order to show that the optimum structural layout is mesh independent and checkerboard free, the design domain is discretized using 40 20, 80 40 and 100 50 quadrilateral elements.Fig. 8 shows the evolution history of the objective function over iterations based on different mesh sizes.It can be observed from Fig. 8 that the number of iterations for mesh sizes of 40 20, 80 40 and 100 50 are 50, 53 and 46, respectively and their corresponding compliances are calculated as 13.75, 14.32 and 14.92, respectively.Fig. 9 shows the curves of convergence of the objective weight based on different mesh sizes.The almost monotonic and uniform convergence can be seen from this figure.
The optimization results obtained by the present method are shown in Fig. 10, the topology optimization history at various iterations is shown in Fig. 11, and the optimization results obtained by the Solid Isotropic Microstructure with Penalization (SIMP) method (Sigmund (2001)) (standard finite element-based method) without and with sensitivity filtering are shown in Fig. 12. From these results it can be seen that the SIMP method using a filter produces similar topologies to the present method, and the present method can effectively remove numerical instabilities such checkerboard pattern and mesh dependency phenomena.The number of iterations of the SIMP method for mesh sizes of 40 20, 80 40 and 100 50 are 41, 96 and 123 respectively and their corresponding compli- ances are respectively calculated as 15.14, 15.12 and 14.53.It is confirmed that the number of iterations of the present method is less than the SIMP method, and a smoother optimization result can be obtained by the present method.

Example 3
The design domain of a simply supported Michelle type structure with a length to height ratio of 6:1 is shown in Fig. 13.The objective function of this example is to minimize the compliance and the objective weight is 50% of the total weight of the design domain.In order to show that the solution is mesh independent and checkerboard free, the design domain is discretized using 84 14, 120 20 and 180 30 quadrilateral elements.Each element is subdivided into four smoothing domains ( 4 nSC = ).Fig. 14 shows the evolution history of the objective function over iterations, and Fig. 15 gives the curves of convergence of the objective weight based on different mesh sizes.The number of iterations for mesh sizes of 84 14, 120 20 and 180 30 are 65, 63 and 55, respectively and their corresponding compliances are calculated as 96.15, 94.82 and 95.82, respectively.It should be noted that the occasional jumps in Fig. 14 may be attributed to a remarkable alteration of topology due to the elimination of one or more bars in a single iteration.
The optimum structural layouts and the topology optimization history at various mesh sizes are shown in Figs.16 and 17, respectively.The optimization results obtained by the Solid Isotropic Microstructure with Penalization (SIMP) method (Sigmund (2001)) (standard finite element-based method) are given in Fig. 18.A comparison between the final solutions shown in Figs.16 and 18 shows that the two different optimization methods generate similar topologies and the present method can avoid numerical instabilities such as checkerboard pattern and mesh dependency phenomena.
The compliances of the solutions of the SIMP method for mesh sizes of 84 14, 120 20 and 180 30 are respectively calculated as 101.92,103.81 and 100.22 which are higher than those of the present method.These differences may be due to the over-estimated strain energy of elements in the solutions of the SIMP method.

CONCLUSIONS
In this paper, the smoothed finite element method is combined with the level set method to develop an efficient approach for topology optimization of continuum structures.The cell-based smoothed finite element method is employed to improve the accuracy and stability of the standard finite element method.Several numerical examples were presented to show the validity and feasibility of the proposed method.The examples have shown the effectiveness of the proposed method to overcome numerical instabilities such as checkerboard patterns and mesh dependency phenomena.As a future research work, the proposed approach can be efficiently used to solve structural topology optimization problems with other objective functionals and constraints such as stress or displacement constraints.

Figure 2 :
Figure 2: Design domain of the cantilever beam.

Figure 3 :
Figure 3: Optimization results obtained from different mesh discretizations with different number of smoothing cells.

Figure 4 :
Figure 4: Evolution histories of the objective function over iterations, example 1.

Figure 5 :
Figure 5: Evolution histories of the objective weight over iterations, example 1.
the structural compliance, ( ) W x is the weight of the current topology, req W of element densities.The design variable e x indi- cates the presence (1) or absence (0) of an element, where e is the element index.e u is the element displacement vector and CS-FEM e K  is the element stiffness matrix for element e .In this study, the discrete level set method is employed to find the solution for the optimization problem.The discretized level set function y can be defined as:

Figure 6 :
Figure 6: Results obtained by the SIMP method from different mesh sizes without and with sensitivity filtering.
this section, three widely studied examples in the field of topology optimization are presented to show the effectiveness of the proposed method.Poisson's ratio 0.3 m = and Young's modulus of 1 E = are used for all examples.

Figure 7 :
Figure 7: Design domain of beam with fixed support.

Figure 8 :
Figure 8: Evolution histories of the objective function over iterations, example 2.

Figure 9 :
Figure 9: Evolution histories of the objective weight over iterations, example 2.

Figure 10 :
Figure 10: Results obtained by the present method with different mesh discretizations.

Figure 12 :
Figure 12: Results obtained by the SIMP method from different mesh sizes without and with sensitivity filtering.

Figure 13 :
Figure 13: Design domain for Michelle type structures.

Figure 14 :
Figure 14: Evolution histories of the objective function over iterations, example 3.

Figure 15 :Figure 16 :
Figure 15: Evolution histories of the objective weight over iterations, example 3.