Abstract
In the structural topology optimization approaches, the Moving Morphable Component (MMC) is a new method to obtain the optimized structural topologies by optimizing shapes, sizes, and locations of components. However, the optimized structure boundary usually generates local nonsmooth areas due to incomplete connection between components. In the present paper, a topology optimization approach considering nonsmooth structural boundaries in the intersection areas of the components based on the MMC is proposed. The variability of components’ shape can be obtained by constructing the topology description function (TDF) with multiple thickness and length variables. The shape of components can be modified according to the structural responses during the optimization process, and the relatively smooth structural boundaries are generated in the intersection areas of the components. To reduce the impact of the initial layout on the rate of convergence, this method is implemented in a hierarchical variable calling strategy. Compared with the original MMC method, the advantage of the proposed approach is that the smoothness of the structural boundaries can be effectively improved and the geometric modeling ability can be enhanced in a concise way. The effectiveness of the proposed method is demonstrated for topology optimization of the minimum compliance problem and compliant mechanisms.
1. Introduction
Structural topology optimization aims to find the best distribution of materials within a prescribed design domain using an optimization algorithm in order to achieve some exceptional structural performance [1]. It breaks the limitations of the designer’s thinking in traditional structure design and can get a more liberal, updated, and lighter conceptual design. At the same time, structural topology optimization is also considered a more challenging design approach compared with the size optimization and shape optimization. Since Bendsøe and Kikuchi [2] developed the topology optimization method based on the homogenization theory, after decades of development and improvement, topology optimization has gradually become an important method in the structural concept design stage and has been applied successfully to a wide range of physical disciplines [3–6]. So far, the solid isotropic material with penalization (SIMP) method [7–9] and level set method [10–14] are two popular strategies. For the other methods, the reader is referred to [15–17]. In the SIMP approach, the design domain is discretized by finite element mesh. Then the pseudodensity of each element is set as a design variable to characterize the form of material distribution. It transforms the structural topology optimization problems into material distribution problems. The main advantages include simpler concepts, high computational efficiency, and easy to implement [18–20]. However, the densities of a certain number of elements are left to represent the final optimal results; the structure boundaries are easy to generate with many grey-scale elements and jagged shapes [21–23]. In the level set method, the values of the level set function at the node points in the finite element mesh are often taken as the design variables. The structural boundaries can be specified by drawing the contour of the level set function in one-higher dimensional space. The optimization process is represented by the evolution of specific faces. Therefore, it often obtains the optimization result with smooth boundaries and easily handles topology changes compared to the traditional topology optimization methods [24–26]. But whether the SIMP approach or level set method, the design variables increase sharply with the mesh refinement, especially in the three-dimension problem, which will far exceed the solution capacity of the existing optimization algorithm [27]. In SIMP and level set methods, an implicit way is used to solve topology optimization problems; that is, without explicitly geometric information, it is implanted. It is difficult to accurately control the structural feature sizes, which is an important factor in considering manufacturability [28]. At the same time, the geometric/topology representation is not quite consistent with that in modern computer-aided design/engineering (CAD/CAE) modeling systems [29, 30]. This treatment sets up a barrier in optimizing structural dimensions and direct link between the final topology results and CAD/CAE systems [31].
In order to overcome the above problems, recently, Guo et al. [27] proposed a more accurate and geometric topology optimization method based on the topology optimization framework, a so-called Moving Morphable Component (MMC). In this framework, a series of components with movable and deformable capabilities are used as building blocks of the topology optimization, and the optimal structural topology that meets specific performance can be obtained by varying the lengths, thickness, tilt angle, and center coordinates of these components. The same idea has also been applied in continuum-based topology optimization with discrete elements [32–36]. Although any curved structural parts can be approximated by precisely controlling a certain number of components with uniform thicknesses [27], so as to enhance the smoothness of the structural boundaries and the geometry modeling capability of the method, it is always highly desirable to include components with variable shapes in the MMC approach. Recently, Zhang et al. [37, 38] improved the geometry modeling capability of the MMC approach by constructing the TDF with variable thickness components. Guo et al. [39] proposed the TDF of the component with a curved skeleton to constructed curved structural parts with a smooth boundary. In Meisam Takalloozadeh et al.’s method [40], the start and end of the component with two straight lines were modified as a concave/convex pointed shape by adding a new variable. Van-Nam Hoang et al. [33] and Deng et al. [32] constructed the rectangle moving component with semicircular ends and enhanced connectivity between components by adding constraints.
Although the geometric modeling capabilities of the above methods and the connection methods of components have been studied, in [37, 40], there are still some obvious nonsmooth boundaries in the optimized structure using the straight lines or concave/convex pointed shapes at the start and end of a component. The method of [39] finds the optimal structural by constructing the components with curved skeleton and most of the boundaries of the optimized structure are smooth, but just like as described at the end of their paper, there still exist some local nonsmooth boundaries especially in the intersection areas of the components. In the optimization process of [32, 33], the thickness of the component changed uniformly, which limits the geometry modeling capability of their methods.
Therefore, based on the MMC framework, this paper mainly studies how to use the variability of components’ shape to solve the local nonsmooth problem in the intersection areas of the components. For this purpose, we control the length of components by adding new design variables in the TDF of the quadratically varying thicknesses. The variability of components’ shape and the geometry modeling capability of the original MMC method are enhanced through the coordination of length and thickness variables. In the topology optimization process, the shape of components can be changed according to the structural responses, and relatively smooth structural boundaries are generated in the intersection areas between the components. In addition, considering the influence of the initial layout with more components on the rate of the convergence, the hierarchical variables calling strategy is proposed to optimize the layout of components in the initial stage. The effectiveness of the proposed method is verified by several numerical examples.
2. Problem Formulation
In the MMC approach, the topology description function (TDF, ) of components can be expressed in the following form:where represents a set of points within the design domain; represents the entire design area; denotes a subset of occupied by components made of the physical material. If the number of components included in the design area is , the TDF of the whole solid material in the structure area can be built aswhere is the number of components and the TDF of the component can be written aswhere is the area occupied by the and . Then the geometric representation of the corresponding structure is shown in Figure 1.

Based on the MMC topology optimization framework, a typical topology optimization problem can be described aswhere and , , denote the objective function and constraint function, respectively, and is the admissible set that belongs to.
Topology optimization for a mean compliance minimization under available volume constraint problem is considered, and the corresponding mathematical model can be written aswhere , , is the vector of the design variables of the component. is the objective function. , , , and are the body force density, the displacement field, the prescribed surface traction on Neumann boundary , and the linear strain tensor, respectively. is the displacement on Dirichlet boundary . is the corresponding test function defined on , and the symbol means the Heaviside function. denotes the topology description function set for the overall component. is an integer (in this work, is used), and is the isotropic elastic tensor. and are Young’s modulus and Poisson’s ratio of the solid material, respectively, and and denote the symmetric part of the fourth-order identity tensor and the second-order identity tensor, respectively. is the upper limit of the volume of the solid material.
3. TDF for the Variability of Components’ Shape and Hierarchical Variable Calling Strategy
In this part, we will discuss how to construct the TDF for the variability of components’ shape and implement the hierarchical variables calling strategy based on the MMC framework.
3.1. TDF for the Variability of Components’ Shape
The difference consists that, in the MMC method, it is possible to give an explicit description of the boundary and geometry features of a component. This cannot be achieved in the conventional level set method. In [27], the uniform thickness component can be described using the hyperelliptic TDF. In order to improve the geometry modeling capabilities of the MMC approach, the TDF of quadraticall varying thicknesses (QVT) can be expressed aswhere is a relatively large positive integer, which takes (6) in this paper. Here, controls the shape of componentswhere , , and denote the coordinate of the center, the half-length, and the inclined angle of component, respectively. , , and are the three thicknesses for different positions of the component, respectively, as shown in Figure 2.

It is worth noting that, in this method, the start and end of the component were two straight lines perpendicular to the direction of the component. When two or more components are connected at an angle, as shown in Figure 3 which plots the components with dashed boundary, the local nonsmooth region of the structural boundary is often caused. If we assume that the straight line of the start and end of the component can be bent and deformed, the variability of the shape of the component can be achieved by combining the thickness variations. Then, the local smooth structural boundaries can be formed in the intersection parts of the components according to the structural response in the optimization process, as shown in Figure 3 which plots the components with solid boundary. The geometry description of the components’ shape with variability is shown in Figure 4. The deformation of components’ shape can be achieved by adjusting the length and thickness variables. The corresponding TDF can be constructed aswhere is the same as equation (7) and can be constructed aswhere , , and are the three lengths for different positions of a component as shown in Figure 4, respectively. and are the center coordinates and radius of the arc of the start and end of a component, respectively, and we have


By inputting different parameters in equation (9), the variability of components’ shape can be easily achieved. Some components are shown in Figure 5, and the values for the lengths and thicknesses variables of components are taken as (a) [0.6 0.6 0.9 0.3 0.3 0.3], (b) [1.0 1.0 0.7 0.3 0.3 0.3], (c) [1.0 0.4 0.7 0.3 0.3 0.3], (d) [1.0 0.1 0.8 0.3 0.3 0.3], (e) [1.0 0.05 0.7 0.08 0.5 0.15], and (f) [1.0 0.1 0.7 0.7 0.7 0.03], respectively. It is worth noting that although the topology optimization process using the variability of components’ shape based on the MMC framework is somewhat similar to the level set method (TDF has been used to represent the geometry of components), the proposed method can also possibly give an explicit description of the boundary and geometric features (such as length and thickness) of components, which cannot be achieved in the level set method that uses the free transformation of the level set functions to represent the boundary of a structure.

(a)

(b)

(c)

(d)

(e)

(f)
3.2. Hierarchical Variables Calling Strategy
In this section, we will discuss how to optimize the initial layout of components to improve the rate of convergence. In MMC framework, although structural topology can be constructed by moving, deforming, disappearing, and overlapping of a set of components, the different initial layouts of components have different effects on the rate of convergence. In particular, when there are more components in the design domain, the effective force transfer path may be affected by improper positioning or deformation of some components, which will delay the rate of convergence.
In order to improve the rate of convergence of the optimization process, a hierarchical variables calling strategy is proposed. In this strategy, the optimization process is divided into two levels: the level 1 is the layout variables calling and the level 2 is all variables calling. The average variation of components’ positioning (the central coordinates and the rotation angles) within the adjacent iterations is set as a constraint. In the initial stage of optimization, only the layout variables (the central coordinates and the rotation angles) are called, and the initial layout of the components can be modified according to the structural responses in this process. When the constraint is satisfied, the full variables are called to construct the topology optimization structure. It is worth noting that the output of the level 1 is used as input for the level 2.
The constraint of the hierarchical variables calling strategy can be expressed aswhere is the iteration, and are the set maximum upper limits of change values, respectively, and generally, the smaller positive numbers are selected. and are the average change values of the center coordinates and rotation angle of components:withwhere and are the center coordinates and the rotation angles change value of the component in twice adjacent iterations, respectively, as shown in Figure 6.

4. Numerical Solution Aspects
In this section, numerical solution aspects in the proposed approach will be discussed in detail.
4.1. Finite Element Analysis
In this study, in order to analyze the structure constructed using the MMC method, the Eulerian mesh is utilized and TFD is used to calculate each node. In order to enhance the computational efficiency, the ersatz material model is used for FEM analysis. Young’s modulus of every element can be calculated by equation (15), which was proposed in [37]:where is an integer number, we take , is Young’s modulus for the material of every element, and is a regularized Heaviside function:where is the regularization parameter and to prevent the occurrence of the singular of the global stiffness matrix, and the stiffness matrix of the element can be obtained by equation (17):
In the above formula, is the stiffness matrix of the element, regardless of the number of materials occupied by components.
4.2. Sensitivity Analysis
Only design sensitivity analysis for the compliance minimization problem is discussed here. Other objectives/constraints can be achieved by resorting to adjoint sensitivity analysis. When we know the stiffness matrix of each element, the global stiffness matrix can be assembled. The structural compliance and the displacement vector will be obtained by using equation (18):
When is assumed to be any one of the variables, then the sensitivity of the structural compliance of can be expressed as
Since the force is independent of the design variables, equation (20) can be transformed into
By equations (15) and (17), we can get equation (22):
By substituting equations (21) and (22) into equation (19), which can be written as follows:where is the element stiffness matrix corresponding to . Because of is a function of , the derivative of the structural compliance can be calculated using equation (23).
For the derivative of volume constraint, we also havewhere can be calculated easily since is an explicit function of .
It is worth noting that, in the above derivations, as mentioned in [37], the nondifferentiable issue arising from the max operation when is constructed by multiple components overlapping is neglected. The numerical case also proves that this has no effect on the optimization process.
5. Numerical Examples
In this section, the effectiveness of the proposed method is verified by several representative numerical examples. Unless otherwise stated, all examples use a four-node bilinear square unit discrete design area. Young’s modulus and Poisson’s ratio of the solid material are taken as and , respectively. The two-dimensional topology optimization problem is considered. The Method of Moving Asymptotes (MMA) is used in this method, and the hyperparameters of MMA are , , , , , and . In addition, all the example codes analysis is performed on the LENOVO-TY5060 workstation with an Intel® CoreTm i5-6400 CPU 2.70 GHz, 4G RAM, and Windows10 and ran with MATLAB R2017a.
5.1. The Compliant Mechanism Problem
To demonstrate the effectiveness of the proposed method in dealing with local nonsmoothness problems in the intersection areas of the components, the compliant mechanism design problem is considered. The design objective is to minimize the geometric advantage of the compliant mechanism. The design domain, geometry data, and boundary conditions are shown in Figure 7. In this problem, the objective is to minimize the geometric advantage (GA) of the compliant mechanism , where and , meaning the displacement in the horizontal direction under a unit force applied at point . The elastic constant of the springs is . The design domain with length and width of is discretized by a FE mesh. The upper bound of volume constraint is set to 0.3. Since the problem under consideration is symmetric in nature, only half of the structure is used to optimize in this example. Sixteen components are placed in the initial design domain, and the overall number of design variables is . The initial values of the design variables for each component of the QVT method and this paper method are set as and , respectively. The initial layout of the components is shown in Figure 8.


(a)

(b)
The intermediate iteration steps and the final topology optimization results for the two methods are shown in Figures 9 and 10. It can be seen that the final structure obtained by the QVT method forms an obvious incomplete connection phenomenon in the intersection area between the components, which easily cause the local structural boundaries nonsmoothness (see Figure 9(d)). From the contour plot of the final structures shown in Figure 10, the components exhibit some flexibility in the form of the connection, and a relatively smooth structural boundary transition is obtained. In the case of the proposed method, the hinge-like connection can no longer be produced in the final optimized structure as illustrated in Figure 10(d).

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)
Figure 11 provides the convergence history of the objective function and the constraint function. It can be observed from Figure 11(b) that, for this example, the value of the objective function has dropped sharply in the earlier iteration step, which indicates that an effective loading transmission path can be faster constructed using the method. The comparison results are shown in Table 1. The values of the objective functions obtained by the two methods are and , respectively. Compared with the QVT method, the proposed method has a faster convergence rate and less time consumption.

(a)

(b)
5.2. The MMB Problem
In the literature on structural optimization studies, this example is often used to verify the validity of a topology optimization method. The design domain and boundary conditions are shown in Figure 12. The unit vertical load is imposed on the middle point of the upper surface of the beam.

The design domain with length and width of is discretized by a FE mesh. The upper bound of volume constraint is set to 0.4. Because of the symmetric property of the problem, only half of the structure is used to be optimized in this example.
Twenty-four components are placed in the initial design domain, and the overall number of design variables is only . The initial values of the design variables for each component of the QVT method and this paper method are set as and , respectively. In equation (12), the maximum upper limits of change values are and . The initial layout of the components is shown in Figure 13.

(a)

(b)
Some steps of the optimization and the convergence histories of the MBB example are shown in Figures 14–16, respectively. From Figure 14(f) and the 810th step in Figure 16, we can observe that the final topology optimization structures obtained by the two methods are very similar.

(a)

(b)

(c)

(d)

(e)

(f)


However, it is worth noting that the QVT method shows a strong instability in the optimization process. In contrast, the proposed method has a more stable convergence process. It can be clearly seen that when equation (12) is not satisfied, the shape and size of the components do not change, as shown the 20th step contour plot in Figure 16. The entire process only optimizes the center coordinates and the rotation angles of components before the 40th iteration. Once equation (12) is satisfied, the components can be moved and deformed according to the structural response, as shown in the contour plot of the 41th and 42th step in Figure 16.
The comparison results by using the two methods are presented in Table 2. The values of the objective function are and , respectively. Compared with the QVT method, the value of the objective function obtained is slightly lower and the time consumption has been reduced with this paper method. The optimization process has been converged at the 810th iteration. It is worth mentioning that the proposed method enhances the stability of the solution process and accelerates the convergence rate.
5.3. The Bridge Problem
In order to demonstrate that this paper method cannot only solve the structural boundaries nonsmoothness in the intersection areas of components but also realize the variability of components’ shape to some extent, the bridge problem is considered in this section. The loading and boundary conditions are shown in Figure 17. The unit vertical load is evenly imposed on the upper surface of the design field. The design domain with length and width of is discretized by a FE mesh. The upper bound of volume constraint is set to 0.4. A rectangle domain with length and width of on the top of the design domain is a nondesign area.

The initial layout is shown in Figure 18. The initial values of the design variables for each component of the QVT method and this paper method are set as and , respectively.

(a)

(b)
From the intermediate steps of the optimization process shown in Figures 19 and 20, the components of this paper method have more flexible deformation ability in the optimization process, and the smooth optimization structural boundaries with circular arc degree are obtained. Although the number of components in the design domain is small, the effective load transmission paths can still be obtained by proper bending deformation of the components. The objective function values are and , respectively. It can be observed that the curved bridge structures have a smaller objective function value.

(a)

(b)

(c)

(d)

(a)

(b)

(c)

(d)
6. Conclusions
In the present paper, a topology optimization method based on the MMC to solve the local nonsmooth problem of components connection areas is proposed. It enhances the variability of components’ shape by coordinating multiple length and thickness variables, and the smooth structural boundaries can be formed in the connection areas between the components. This method preserves the advantages of the original MMC approach and enhances its geometry modeling capabilities. The effect of undesirable positioning or deformation of some components in a layout with more components on the rate of convergence can be optimized by the hierarchical variables calling strategy. Several numerical examples are used to verify the effectiveness of the proposed method. It is also observed that the initial layout of the components is very important for the rate of convergence. In principle, the hierarchical variables calling strategy can be applied to the structural topology optimization that needs to optimize the initial layout of components to improve the rate of convergence. But how to achieve a reasonable initial components layout way is a more challenging task, and the related research will be reported in separate works.
Data Availability
The data used to support the findings of this study are available from the corresponding author upon request.
Conflicts of Interest
The authors declare there are no conflicts of interest regarding the publication of this paper.
Acknowledgments
This project was supported by the National Key Research and Development Program of China [Grant no. 2017YFB1102804].