Abstract
A robust topology optimization (RTO) approach with consideration of loading uncertainties is developed in this paper. The stochastic collocation method combined with full tensor product grid and Smolyak sparse grid transforms the robust formulation into a weighted multiple loading deterministic problem at the collocation points. The proposed approach is amenable to implementation in existing commercial topology optimization software package and thus feasible to practical engineering problems. Numerical examples of two- and three-dimensional topology optimization problems are provided to demonstrate the proposed RTO approach and its applications. The optimal topologies obtained from deterministic and robust topology optimization designs under tensor product grid and sparse grid with different levels are compared with one another to investigate the pros and cons of optimization algorithm on final topologies, and an extensive Monte Carlo simulation is also performed to verify the proposed approach.
1. Introduction
Structural optimization is an intrinsic element of engineering design for the purpose of improving the structural performance and reducing costs. The major ingredients of structural optimization are, in general, divided into three levels, namely, sizing, shape, and topology optimization. The aim of classic topology optimization is to obtain an optimal material or design parameter distribution at a given structural domain under nominal (i.e., deterministic) material properties, geometry, and loading conditions. Until now, a range of topology optimization strategies has emerged, including the homogenization method, introduced by Bendsøe and Kikuchi [1], the density-based method [2], bilinear evolutionary structural optimization (BESO) [3], and level set method [4]. Details of various representative schemes can be found in comprehensive reviews and text books [5–7].
Traditionally, structural topology optimizations are conducted in a deterministic manner, known as deterministic topology optimization (DTO), where the design is determined without taking into account various sources of uncertainties [6]. However, uncertainties are unavoidably observed in real-world applications due to insufficient knowledge, manufacturing errors, changeable environment, and so forth. This may lead to the vulnerable optimum structure or infeasible topologies due to the fluctuation of the structure performance. Therefore, there is a strongly increasing requirement to take the effect of uncertainty into consideration for optimal topologies in structural design.
Recent years have witnessed various statements accounting for uncertainty in topology optimization process. In reliability-based topology optimization (RBTO), optimal designs with acceptable failure probabilities have been achieved, where uncertainty is manifested in the probability constraint functions. Many works dealing with applications of RBTO in structural design have appeared over the last decades. One of pioneering works can be traced back to Kharmanda et al. [8], where reliability analysis is decoupled from optimization procedures and carried out at the beginning of the optimization loops, followed by equivalent DTO. Different optimization strategies have been increasingly applied for reliable designs in various mechanic and multiphysical fields. To cite a few, optimization algorithms via the density-based method have been demonstrated for microelectromechanical systems (MEMS) [9], geometrically nonlinear structures [10], and an automotive control arm [11]. Optimization problems via the BESO have been reported for a vehicle’s hood reinforcement [12] and electrothermal-compliant mechanisms [13]. Topology optimization of compliant mechanisms with loads, material properties, and member geometries uncertainties is demonstrated for the level set method [14].
Distinct from the abovementioned formulations considering the uncertainties in structural topology optimization problems, robust topology optimization (RTO) has recently emerged to optimize the objective performance while simultaneously minimizing its sensitivity with respect to uncertainties. When the uncertainty is assigned unknown-but-bounded convex models, the worst case approach is followed in general. De Gournay et al. [15] propose the so-called worst case or robust optimal design problem for minimal compliance based on the level set method. Amir et al. [16] demonstrate a fair comparison of a worst case formulation and a stochastic formulation in robust topology optimization procedure. Guo et al. [17] present structural topology optimization considering the uncertainty of boundary variations via level set approach, where the objective function is modeled by the compliance and fundamental frequency of structure enduring the worst case perturbation.
In the context of uncertainty parameters defining probability distributions, the probabilistic description of the performance function is usually characterized by its statistical quantities such as mean value and standard deviation. In this scenario, a fundamental issue associated with these RTO procedures is to calculate the statistical moments accurately and efficiently. The Monte Carlo simulations (MCS) are the most widely used sampling-based methods in robust design due to their accuracy and easy implementation. However, it is time-consuming. Schevenels et al. [18] present a robust topology optimization approach for the design of macro-, micro-, or nanostructures, accounting for spatially varying manufacturing errors, where the statistical moments and their sensitivities are estimated by means of the MCS. Jansen et al. [19] propose a robust approach to topology optimization taking into account geometric imperfections, where a sampling method (i.e., 100 Monte Carlo samples) is used to estimate these statistics that are defined as a weighted sum objective function. A new efficient and accurate approach to robust structural topology optimization under loading uncertainty is developed by Zhao and Wang [20]. The Monte Carlo method and matrix decomposition are employed for concentrated loads, and the Karhunen-Loeve expansion and its orthogonal properties of random variables are used for distributed loads. Alternatively, the stochastic collocation method has been gaining more attention in robust design. The idea is to calculate the statistical moments by direct numerical integration. Chen et al. [21] employ the univariate dimension-reduction (UDR) method combined with Gauss type quadrature sampling for calculating statistical moments with consideration of random field uncertainty in loading and material properties. Lazarov et al. [22] introduce the stochastic collocation methods in topology optimization for mechanical systems with material and geometric uncertainties. The robust designs are obtained by utilizing well-developed deterministic solvers.
Efficient and accurate generation of statistic moments is a central problem for robust design, especially for structural topology optimization with infinite-dimensional property. In this paper, we present a methodology for robust topology optimization problems with the stochastic collocation method to integrate the statistic moments considering loading uncertainty. Then, the robust topology optimization formulation is transformed into a weighted multiple loading deterministic problem at the collocation points. For the construction of the stochastic collocation method, the collocation points with corresponding weights are generated by full tensor product grid and Smolyak sparse grid, respectively. Integrations based on full tensor product suffers from the curse of dimensionality, while the sparse grid method originating from the Smolyak algorithm will dramatically reduce the number of collocation points without losing much accuracy. In addition, existing commercial software packages are capable of application in robust topology optimization in these settings.
The paper is organized as follows: A brief description of the deterministic and robust topology optimization formulation is presented in Section 2. In Section 3, the stochastic collocation methods including tensor product grid and sparse grid are demonstrated. The computational results for deterministic and robust topology optimization together with several numerical examples are illustrated in Section 4. Conclusion and future works are discussed in Section 5.
2. Topology Optimization Formulation
2.1. Deterministic Topology Optimization Formulation
The generalized topology optimization problem as shown in Figure 1 is that the predetermined design domain , with respect to the given loads and boundary conditions, is represented by the material distribution function so as to optimize given objective function (e.g., compliance) under limited volume constraint.
The structural topology optimization considering the general elasticity problem is reformulated in the uniform context of the principle of virtual displacement. The optimization problem for deterministic design can be written in a weak form aswhere is the space of kinematically admissible displacement fields. is the equilibrium displacement field, and is the arbitrary virtual displacement in the space . is the prescribed volume constraint. The energy bilinear form and the load linear form are defined as follows:where is the strain tensor. and are the body forces and the boundary tractions on the boundary , respectively. Considering linear isotropic materials, Young’s modulus is a function of design variables given by the density-based method [6] as follows: where is Young’s modulus of solid material and is the penalization power to ensure solid-void solutions.
By using FEM discretization for the linear elasticity, the optimization problem can be rewritten as follows:where is the compliance of the continuum structure, is the displacement vector, and is a set of finite elements. and are the element displacement and the element stiffness matrix with unit Young’s modulus, respectively, is the element volume, and the density of each element is represented by one density variable . Therefore, the standard FEM global stiffness matrix is expressed aswhere is the element stiffness matrix, is the strain-displacement matrix of shape function derivatives, is element domain, and is the element constitutive matrix, defined in 2D aswhere denotes Poisson’s ratio. In general, the optimal solution obtained by direct implementation of the method is prone to problems with checkerboards and mesh dependency. A large number of restriction schemes have been proposed to ensure manufacturability and mesh independency, such as perimeter control method, mesh-independent filtering methods, and Heaviside projection. A comprehensive comparison of the existing restriction schemes is provided in the paper of Sigmund [23].
2.2. Robust Topology Optimization Formulation
The topology optimization of linear elastic structures subject to loading uncertainties is considered in this paper. Due to the variation of uncertain parameters , the performance function has probability distribution description. Thus, the robustness of the design objective can be achieved by optimizing the mean performance and minimizing the standard deviation of the performance simultaneously, as shown in Figure 2. In Figure 2, point denotes the deterministic optimum not considering the uncertainty, while point denotes the robust optimum. The performance of deterministic optimum is better than that of the robust optimum. However, its variation is wider than the robust optimum.
Therefore, a main challenge in performing RTO for realistic problems is the implementation of an algorithm, which inherently represents a multiobjective optimization problem. One efficient way is to combine all objective functions into a single weighted function. In mathematical terms, a general robust topology optimization task can be stated aswhere is the weight to be determined by the designer that satisfies . Note that the weights have to be set a priori, which make a trade-off decision between an optimal mean performance and minimization of standard deviation of the performance to obtain a compromise solution. and are the mean value and the standard deviation of the performance function, respectively. is the material resource constraint. The design variables are bounded by their lower and upper bounds and , respectively. In a probabilistic setting, the mean value and the standard deviation can be expressed aswhere is the expectation operator and is the joint probability density function (PDF) of the random variables. In practice, it is very difficult and even impossible to calculate the mean value and the variation of the performance function through the multidimensional integration in (8). Such difficulties have motivated the development of various numerical methods such as simulation-based methods, local expansion-based methods, and the stochastic collocation methods. The reader is referred to Lee and Chen [24] for a more detailed review.
3. Stochastic Collocation Method
The main idea of stochastic collocation method is to sample the quantity of interest at specific collocation points in the random space, where the sample points are deterministic and are associated with quadrature formula for the evaluation of integral statistics such as mean and standard deviation. Examples of such quadrature rules include the Newton-Cotes (midpoint, rectangle, and trapezoidal) rule, the Gauss quadrature (Legendre, Chebyshev, Laguerre, Hermite, Jacobi, Kronrod, and Patterson) rule, and the Clenshaw-Curtis rule [25–27]. In the case of a one-dimensional random space, these points are quadrature abscissas selected according to the probability measure of the random variable. Tensor product grid and sparse grid constructions are then utilized to generate the multidimensional abscissas.
3.1. Tensor Product Grid
Without loss of generality, consider the class of function over the -dimensional unit hypercube with numerical integration representationIn the one-dimensional case, a sequence of univariate quadrature operators by -point quadrature formulas with level iswhere and are the collocation points and weights, respectively. To obtain a quadrature formula for the multivariate case , the full tensor product formula is defined asThen, the mean and standard deviation of the objective performance are derived from the following equations:Clearly, the tensor product grid requires collocation points, which are sampled on the full grid. Assuming the same level accuracy in each dimension—that is, —the total number of points is . Thus, the computational burden increases significantly with the dimension . This is often referred to as the curse of dimensionality [28].
3.2. Sparse Grid
The sparse grid method is a special discretization technique, which can be traced back to the Smolyak algorithm [29]. It is based on hierarchical basis, a representation of a discrete function space which is equivalent to the conventional nodal basis, and a sparse tensor product construction [30, 31]. The readers are referred to the review paper by Xiong et al. [32] and the references therein for more information. The one-dimensional difference quadrature formula is defined as Then the Smolyak quadrature formula for -dimensional functions with level is given by Note that denotes the summation of the multi-indices . Alternatively, the above formula can be written asTo compute , a specific tensor product rule of sparse grid is needed, defined as where denotes the set of collocation points for one-dimensional -level accuracy. With the Smolyak algorithm, the weight corresponding to the th collocation point is defined as Then, the mean and standard deviation of the objective performance by sparse grid method can be computed bywhere is the total number of points in sparse grid that is determined asSuppose that the same quadrature rule is chosen under the dimension with same level in each dimension, and represents the collocation points in each dimension and the tensor-product grid scales as , whereas the sparse grid scales as , significantly mitigating the curse of dimensionality. Accurate quadrature can be created by using higher values of the level ; however, this also leads to more sampling points.
The sparse grid method discretizes the stochastic space in hierarchical structure based on nested collocation nodes, such as the Chebyshev nodes or Gauss-Patterson nodes, leading to Clenshaw-Curtis rule or Gauss-Patterson rule, respectively. Here, the Clenshaw-Curtis type sparse grid design , with nonequidistant nodes, is constructed [33–35]. The collocation points are defined asThe sequence is chosen as and the corresponding weights are given by
Examples of two-dimensional collocation points of full tensor product grid and Smolyak sparse grid based on Clenshaw-Curtis rule with level are shown in Figure 3. It is clearly observed that, compared with full tensor product grid, the sparse grid dramatically reduces the number of collocation points, while preserving a high level of accuracy.
(a) Tensor product grid
(b) Sparse grid
4. RTO Algorithm and Numerical Examples
4.1. Optimization Algorithm
The implementation of optimization algorithm for RTO is based on the commercial software package Optistruct 11.0 of Altair Hyperworks, where the solid isotropic material with penalization is employed as material interpolation model and convex programming method as optimization solver. In addition, the coordinates and corresponding weights of the collocation points for random variables are evaluated in an external environment of Optistruct for auxiliary conditions. The optimization algorithm can be summarized as follows:(1)State problem discretization and initialization, including initial design domain and boundary conditions .(2)Identify the random loading uncertainty as well as the related loading angle and magnitude probability distributions.(3)Compute collocation points and the corresponding weights by tensor product gird and sparse grid, according to the dimension and level information.(4)Define multiple loading cases through mapping the loading uncertainty with the collocation points.(5)Establish the mean and standard deviation of the performance function through specifying the user-defined equations, according to the prescribed format.(6)Implement the weighted function by combining the mean and standard deviation functions with predefined weight , and choose the function as the objective function to be minimized.(7)Define the manufacturing restriction schemes, optimization strategy, and volume fraction constraint.(8)Run the optimization loop until the convergence criterion is satisfied.
4.2. Numerical Examples
In this section, we will illustrate the proposed RTO approach using two numerical example problems—a 2D simple supported beam and a 3D L-shaped structure—to demonstrate the effectiveness and efficiency for tensor product grid and sparse grid. A simple Monte-Carlo simulation is also performed to verify the accuracy of the statistical moment calculated by the proposed approach. The procedure of topology optimization is conducted on a Windows 7 workstation (Intel Xeon (R) E5440, 2.83 GHz, 4.00 GB RAM, 8 cores).
4.2.1. A 2D Simply Supported Beam
As the first example, the design domain, the boundary, and loading conditions of the simply supported beam are illustrated in Figure 4. The dimensions of the beam are mm and mm and the thickness is mm. The isotropic material is assumed to have Young’s modulus of 1 MPa and Poisson’s ratio of 0.30. A single random load case is considered with which the angle and magnitude are characterized by two independent random variables. The angle is assumed to follow a uniform distribution with the interval of . The magnitude follows a Gaussian distribution with mean of 1 and standard deviation of 0.3. In robust design, the optimization problem is to minimize the weighted summation of the mean and standard deviation of the compliance of the structure under the volume constraint, where the weight is set to be 1. The design domain is discretized by 2700 (90 × 30) four-node quadrilateral elements.
The optimal topologies obtained by deterministic and the proposed robust topology optimization approaches are presented in Figures 5 and 6, and the corresponding results are compared in Table 1, including the required number of points, statistical information of the performance function, volume fraction constraint, iterations, and computing time. Here, full tensor product grid and sparse grid under level , and 3 are investigated in this section. Extensive 1000 Monte Carlo simulations are also performed to verify the accuracy of the proposed robust approach.
(a) Deterministic design
(b) Monte Carlo simulation
(a) Tensor product formula
(b) Sparse grid method
(c) Tensor product formula
(d) Sparse grid method
(e) Tensor product formula
(f) Sparse grid method
For comparison purposes, the deterministic design is performed using the mean values of the load and the robust design is conducted with multiple load cases due to randomness of load angle and amplitude. The DTO and RTO resulted in significantly different topologies, as shown in Figure 5. The remarkable difference lies in the robust design exhibiting an asymmetric layout compared with the deterministic design characterized by symmetric configuration. It can be attributed to the fact that the deterministic design suffers from symmetric vertical constraints at the lower left and the lower right corner. However, in robust design, the structure possesses an asymmetrical degree of freedom constraint in horizontal direction under the condition of horizontal component of load.
As the results listed in Table 1 show, the robust designs present more robust topologies than the deterministic approach, according to the results summary of the statistical information. This emphasizes the importance of considering the uncertainties in load for structure design. It is also noted that the corresponding computation cost for robust designs is larger than for the deterministic design since multiple load cases are required for calculating the mean and deviation values of the performance function. The results of Monte Carlo simulations in Table 1 show that the proposed robust approach possesses improved efficiency without losing much accuracy.
From the comparison of the results, with the improvement of level accuracy, the number of collocation points in tensor product grid and sparse grid corresponding increases. As the results listed in Table 1 show, the mean and standard deviation of the compliance for sparse grid and tensor product grid are gradually close to those for Monte Carlo simulation under the increase of level . This indicates that more robust designs are clearly obtained under higher level . Meanwhile, the higher computation cost is required due to the increase in the number of evaluations of performance function responses. This means robustness enhancement is inevitably accompanied with the sacrifice of convergence time. It can also be seen from the results listed in Table 1 that, under the same level , sparse grid uses fewer collocation points than those of tensor product grid to produce similar even more robust results. At the same time, compared with sparse grid, the computation burden for the exponential dependent calculation of tensor product grid can be overcome to a certain extent. This means sparse grid, in a manner, performs better than the tensor product grid in compromise between robustness and cost in robust design.
4.2.2. A 3D L-Shaped Structure
The design domain and loading conditions of this example are illustrated in Figure 7. The geometrical dimensions of the structure are mm and mm. Young’s modulus and Poisson’s ratio of the fully solid material are 2.1 × 105 MPa and 0.30, respectively. Two concentrated forces, and , are applied at central point of right bottom surface simultaneously. The magnitudes of the forces are characterized by two uncorrelated random variables, which follow Gaussian distributions with mean of 1000 N and 30% variance from the mean values. The optimization problem is still to minimize the weighted summation of the mean and standard deviation of the compliance of the structure under the volume constraint. The design domain is divided into 1.2 × 104 eight-node hexahedron elements.
The optimal topologies obtained by deterministic, Monte Carlo simulation and the robust approaches with density value above 0.3 thresholds are presented in Figures 8 and 9, respectively. The corresponding optimization results are summarized in Table 2, including the required number of points, statistical moments, volume constraint, iterations, and computing time. Meanwhile, full tensor product grid and sparse grid are investigated under level , and 3.
(a) Deterministic design
(b) Monte Carlo simulation design
(a) Tensor product grid
(b) Sparse grid
(c) Tensor product grid
(d) Sparse grid
(e) Tensor product grid
(f) Sparse grid
Comparing the deterministic design with the proposed robust designs, obvious distinctions can be clearly observed in optimal topologies as shown in Figures 8 and 9. From the viewpoint of geometry, the robust designs represent more rational material layouts, which lead to more robust structural configuration. According to Table 2, the final statistical moments of the robust designs are superior to that of deterministic design. This indicates that the robust design is less sensitive to the variations of load magnitude. It can also be seen from the results listed in Table 2 that as the level is increasing, the mean and standard deviation of the structure compliance for tensor product grid and sparse grid are closer to those for Monte Carlo simulation in which more robust designs are obtained, but it is also revealed that the corresponding computing time is increasing. This means robustness and cost are conflicting criteria in robust topology optimization design. For comparison of the results of tensor product grid and sparse grid, these indicate that sparse grid, in a manner, performs better than the tensor product grid comprehensively considering the accuracy and efficiency.
5. Conclusion
This paper integrates the stochastic collocation methods into structural topology optimization for solving the RTO problem with random loading uncertainty. Specially, tensor product grid and sparse grid are employed to transform the RTO into a simplified weighted multiple loading deterministic problem at the collocation points. The proposed robust approaches are implemented, discussed, and tested on two numerical examples by implementation in existing commercial topology optimization software package. An extensive Monte Carlo simulation is also performed to verify the accuracy and validation of the proposed approach. From the comparison of the results, the following conclusions are summarized:
The proposed RTO approach possesses more robust topologies than DTO and may exhibit manifest distinct topology configuration. It is revealed that robust topology optimization is essential when involving the uncertainty parameters in optimization design. The stochastic collocation methods are compatible with RTO for calculating statistical moments of the objective performance response and are readily implemented in existing commercial topology optimization software package.
For tensor product grid and sparse grid, more robust results are clearly obtained under higher level setting. Meanwhile, the required computation cost is more expensive. This indicates that robustness enhancement is inevitably accompanied with the sacrifice of computational expense. Sparse grid demands fewer collocation points than tensor product grid to generate similar even more robust results under the same level. At the same time, the computation burden for tensor product grid can be overcome to a certain extent. It is demonstrated that sparse grid, in a manner, performs better than the tensor product grid in compromise between robustness and computation cost.
It is important to notice that, compared to deterministic topology optimization, the application of robust topology optimization to practical structural engineering is still a challenge task because of the high RTO computation cost. Future research is still recommended to implement the proposed RTO approach for optimization problems subject to non-Gaussian type random loading uncertainty.
Conflict of Interests
The authors declare that there is no conflict of interests regarding the publication of this paper.
Acknowledgments
This research is supported by National Natural Science Foundation of China (Grants no. 51275040 and no. 50905017) and the Programme of Introducing Talents of Discipline to Universities (no. B12022).