Topology optimization incorporating external variables with metamodeling

The objective of conventional topology optimization is to optimize the material distribution for a prescribed design domain. However, solving the topology optimization problem strongly depends on the settings specified by the designer for the shape of the design domain or their specification of the boundary conditions. This contradiction indicates that the improvement of structures should be achieved by optimizing not only the material distribution but also the additional design variables that specify the above settings. We refer to the additional design variables as external variables. This paper presents our work relating to solving the design problem of topology optimization incorporating external variables. The approach we follow is to formulate the design problem as a multi-level optimization problem by focusing on the dominance-dependence relationship between external variables and material distribution. We propose a framework to solve the optimization problem utilizing the multi-level formulation and metamodeling. The metamodel approximates the relationship between the external variables and the performance of the corresponding optimized material distribution. The effectiveness of the framework is demonstrated by presenting three examples.


Introduction
Topology optimization is a promising methodology for structural design. Compared with other structural optimization methodologies such as size optimization and shape optimization, topology optimization has a high degree of freedom in structural representation, and can therefore generate new design concepts. Since the pioneering work by Bendsøe and Kikuchi (1988), various methods to achieve topology optimization have been proposed, and the scope of application has widely expanded. The basic theory and development of the approach are described in the literature (e.g., Eschenauer and Olhoff 2001;Rozvany 2009;Sigmund and Maute 2013;Deaton and Grandhi 2014;Zhu et al. 2016 Any optimization framework requires pre-processing that depends on the designer's decisions. In the case of conventional topology optimization, design problems are defined using a fixed design domain and boundary conditions. These conditions need to be specified in advance by the designer. These settings are usually referred to as design parameters in the sense that they are outside the framework of optimization.
Although the design parameters influence the solution of topology optimization in general, the specific values that would need to be assigned to the design parameters to improve the performance of optimized material distribution are not necessarily obvious in advance. For example, let us consider the problem of maximizing the stiffness of a cantilever beam. Clearly, different optimal solutions would be obtained if the design parameters, such as the shape of the design domain or the positions specified for the boundary conditions of the fixed supports, are different. From the viewpoint of the original purpose of designing the cantilever beam, these design parameters could also be regarded as design variables that specify the structure. Therefore, it is desirable to extend the conventional topology optimization problem to an optimization problem that includes these design parameters in the form of design variables. We refer to these design variables by introducing a term external variables. Our work involves topology optimization that incorporates these external variables, and aims to find the optimal material distribution and external variables. Topology optimization incorporating external variables allows us to explore a larger design space.
The introduction of the concept of external variables reveals the existence of a dominance-dependence relationship between the external variables and material distribution. That is, the specification of external variables determines the design domain and boundaries. Then, depending on these boundaries, the material distribution is given within the domain. Obviously, the two types of design variables are qualitatively different. Focusing on this characteristic, we replace the original optimization problem with a set of two problems, which we refer to as a multi-level formulation approach. The first is a sub-problem of the original problem, which optimizes only the material distribution by assigning constant values to the external variables. This is termed the lower level optimization problem. The other is an optimization problem in which the original problem is rewritten in consideration of the dominance-dependence relationship. In this problem, the optimal external variables are searched based on the result of the lower level optimization. This is termed the upper level optimization problem.
Contrary to our point of view, most related works employed a single-level formulation approach, which simultaneously updates the material distribution and the external variables. A well-studied application is the so-called integrated layout design of multicomponent system, which is a design problem that entails simultaneously optimizing the material distribution of load-bearing members and the positions of fixed-shaped parts embedded in them (Qian and Ananthasuresh 2004;Zhu et al. 2009). The design variables that specify the position can be regarded as one form of external variables. Derived from these studies, Zhu and Zhang (2010) solved the optimization problems of material distribution and the positions at which fixed support is provided, by regarding the fixed support as highly rigid parts. In recent years, it has been applied not only to rigidity maximization problems but also to thermal problems (Li et al. 2017).
However, owing to the dominance-dependence relationship, the modeling of structural representation is challenging when a single-level formulation approach is employed. This is because variation of the external variables in an optimization process could modify the domain in which the material distribution is given, in which case the function representing the material distribution would be inconsistent in the process. In practice, variation of the external variables could distort the finite element mesh in which the material distribution is defined. Consequently, techniques are required to derive the design sensitivity with respect to the external variables. Various methods have been proposed to overcome these difficulties (Qian and Ananthasuresh 2004;Zhang et al. 2012;Xia et al. 2013;Gao et al. 2015;Kuci et al. 2018). Another problem associated with the single-level formulation approach is that, for certain design problems, the optimized solution is known to strongly depend on the initial value of the external variables (Zhang et al. 2011). In the case of close interaction between the material distribution and external variables, the initial value dependency could affect the performance of the optimized solution. In such a situation, a trial and error process would be necessary to find the appropriate initial value for the external variables.
On the other hand, Li et al. (2001) employed the multi-level formulation approach by explicitly utilizing the dominance-dependence relationship. They designed a piezoelectric actuator by combining the positional and sizing optimization of a piezoelectric component and topology optimization of a flexible mechanical structure. They optimized the position and size of the piezoelectric component in the upper level with a genetic algorithm (GA), and optimized the shape and topology of the flexible mechanical structure in the lower level with sequential linear programming (SLP). However, the external variables in their method only take discrete values, because the position of the piezoelectric component is only allowed to fit into a fixed background finite element mesh. Furthermore, the computational cost is unacceptably high owing to the double loop algorithm including the GA.
The aim of our work is to devise a new optimization framework to solve the problem of topology optimization incorporating external variables using a multi-level formulation. A key aspect of the framework is to introduce the metamodeling technique at the upper level. We construct a metamodel to approximate the relationship between external variables and the performance of the corresponding optimized material distribution. Based on this idea, a singleloop optimization algorithm is proposed. The metamodel is designed to replace the objective function at the upper level. The sensitivity analysis of the objective function at the upper level is replaced with that of the metamodel. The sensitivity of the metamodel is derived analytically and the computational cost of this calculation is usually low. In this framework, the external variables can take continuous values. Furthermore, at the lower level, conventional topology optimization methods can be used without modification.
The paper is organized as follows. The general formulation of the topology optimization incorporating external variables is described in Section 2. Details of the proposed framework are provided in Section 3. The proposed framework is verified using three numerical examples, which are included in Section 4. The limitations and possibilities of the proposed framework are discussed in Section 5. Finally, conclusions are presented in Section 6.

General formulation of target problem
In this section, we describe the problem of topology optimization incorporating external variables in a general form. First, we briefly describe the mathematical representation of the material distribution, design domain, and its boundary in the conventional topology optimization. Second, we describe the relationship between them and the external variables. Finally, we formulate the optimization problem.
Topology optimization is a methodology to optimize both the shape and topology. The basic concept of topology optimization consists of replacing the structural design problem with a material distribution optimization problem within a fixed design domain. The material is represented as a domain in the design domain D, and the cavity is represented by D\ . Several methods are available to express the material distribution, such as the homogenization, density, and level set methods. In this study, we employ the density method because of its theoretical simplicity. To represent and D\ , a pseudodensity function ρ(x) is introduced at any position x in D. ρ(x) takes continuous values from 0 to 1. ρ(x) represents that x is in a material domain if ρ(x) = 1, a cavity domain D\ if ρ(x) = 0, or intermediate state if 0 < ρ(x) < 1. The performance of a represented structure is calculated based on state equations and boundary conditions in an analysis domain including D and its boundary. ρ(x) is updated by mathematical programming to obtain the optimal structure. The material property represented by ρ(x) is given by the solid isotropic material with penalization (SIMP) model (Bendsøe and Sigmund 1999).
In conventional topology optimization, the shape of the design domain D and position of a boundary for which the boundary conditions are given are prescribed. In this study, we explicitly consider the external variables s, which specifies D and . The external variables s are, for example, certain quantities that specify the shape of the design domain or the positions of the boundaries of the analysis domain. A dominance-dependence relationship exists between ρ(x) and s. Hereafter, when explicitly describing that ρ(x) is given depending on s, the material distribution is denoted by ρ s (x). When external variables s are uniquely specified as s = s (1) , a design domain D and boundaries are determined, and then a pseudo-density function ρ s (1) (x) is given within the domain. If the external variables take a different value s (2) , another design domain D and other boundaries are determined, then another pseudo-density function ρ s (2) (x) is given.
The problem of topology optimization incorporating external variables is formulated in a general form as follows.
where F (s, ρ s ) is an objective functional to be minimized, G i (s, ρ s ) is the ith constraint functional, and n is the number of constraint functionals. In common topology optimization theory, ρ(x) is defined as a pseudo-density function in the design domain. Meanwhile when implemented, ρ(x) is discretized by the finite element method (FEM), and a set of design variables ρ is used instead. Hereafter, we treat the design variables of the material distribution as a finite dimensional vector ρ s , in which each component is defined on each finite element mesh. In this study, the optimization problem in Eq. 1 is addressed under the assumption that the size of s is small, because this would enable the computational cost of metamodeling for an accurate approximation to remain low.

Overall optimization procedure
In the proposed framework, the original optimization problem (1) is transformed into a two-level optimization problem and solved using conventional topology optimization and metamodeling techniques. An overview of the proposed framework is illustrated in Fig. 1. The graphic illustration on the left in Fig. 1 shows the objective function of the original optimization problem, where the design variables corresponding to material density ρ and external variables s are assumed to be one-dimensional for explanation purposes. Note that although ρ s is defined individually depending on s, they are represented by one representative ρ for the sake of simplicity. The two-level optimization problem consists of an upper level optimization in which s is optimized, and a lower level optimization in which ρ is optimized. At the lower level, multiple topology optimizations are performed by assigning a fixed value to s. This corresponds to finding a local minimum in each F − ρ cross section of the objective function. At the upper level, s is optimized based on the objective function values of the structures optimized at the lower level. Let H (s) be a function that returns this value ( Fig. 1 middle). However, optimization at the upper level is usually impractical due to high computational cost in evaluation of the objective function and difficulty in sensitivity analysis. Therefore, the objective function is replaced with an approximate functionH (s) using metamodeling. This is Minimize H(s) s Fig. 1 Schematic illustration of the proposed framework and an example of an optimization problem in which the design variables ρ s and s are both one-dimensional a process of approximating computationally intensive functions with simple analytical models, which are known as metamodels (Wang and Shan 2007). That is, s is optimized on the metamodel, after which topology optimization is performed again under the condition that the optimized external variables are fixed ( Fig. 1 right). The abovementioned procedure enables the optimal solution ρ and s to be obtained, where the superscript denotes the optimized values of each design variable. The two-level optimization problems are formulated as follows. First, the lower level problem is At this level, conventional topology optimization is performed under the condition that the external variables are fixed. Second, the upper level problem is where H (s) is the objective function in the upper level and is defined as where ρ s and s satisfy the constraints G i (s, ρ s ) ≤ 0 (i = 1, 2, ..., n). H (s) is a function that returns a value of the objective function in the lower level with an optimized solution F (s, ρ s ). The upper level optimization procedure involves substantial cost particularly for the evaluation of H (s), because it is necessary to conduct a sequence of topology optimizations. Moreover, sensitivity analysis is usually difficult because the design domain does not always have a specific shape. We therefore introduce metamodeling at the upper level to address the above difficulties. H (s) is approximately represented as a metamodelH (s). Accordingly, the upper level optimization problem Eq. 3 is reformulated as H (s) is constructed from multiple topology optimization results with metamodeling techniques (details are provided elsewhere 3.2). These multiple topology optimizations can be performed in parallel. The sensitivity ∂H (s)/∂s is not required in upper level optimization; ∂H (s)/∂s is used instead.
In the proposed framework, optimization is performed according to the following procedure: Note that, as the approximation error is inseparable in metamodeling in some cases, the quality of the metamodel should be confirmed after step (v). Then, if necessary, return to step (i) to reconsider the sampling method; for example, add more sampling points or reduce the sampling range.

Metamodeling techniques at the upper level
Generally, metamodeling involves: (Step a) Choosing an experimental design for generating data, (Step b) Choosing a model to represent the data, and (Step c) Fitting the model to the observed data.
Several options are available for each step, and the designer decides which one to use depending on the problem. Useful summaries of the techniques for these steps can be found in the literature (see for example Simpson et al. 2001;Wang and Shan 2007).
Although the selection of an approximation model (step b) depends on the user's decision, this step greatly influences the effectiveness of the proposed framework. We therefore recommend using leave-one-out cross-validation (Meckesheimer et al. 2002) for model assessment at step (iii). This method entails extracting one data value from the sample data group as test data value and using the remaining data to construct the metamodel. The verification is repeated until all data values have been used as test data once. In this study, leave-one-out cross-validation is used to determine an appropriate parameter setting for an approximation model. Designers first select an approximate model, and then select a parameter with the smallest average root mean squared error (RMSE) from the validation results.
In the proposed framework, the metamodelH (s) is constructed from the results of topology optimizations. First, a sampling range is set. Then, a set of external variables is arranged as sample points with any design of experimental techniques, e.g., full factorial design or Latin hypercube design. The kth sample points are: where m is the number of external variables. Second, as sampling, topology optimizations are performed under the condition of the external variables s (k) being fixed for each k, and then values of H (s (k) ) are obtained. Third, an approximation model is chosen, e.g., a polynomial function, radial basis function (RBF), or Kriging model, the parameters of which are based on the result of the leave-one-out cross-validation. Finally, metamodelH (s) is constructed from sampling data; that is, the values of H (s (k) ), with any of the fitting techniques, e.g., leastsquares regression or back-propagation.

Numerical examples
The proposed framework was verified by using three numerical examples. The first example is a relatively simple design problem, in which the theoretical optimum is obvious under the mechanics of materials. In addition, in the proposed framework, it is important to appropriately set the parameters of the approximate model in the metamodeling process; therefore, an example of this selection by leaveone-out cross-validation is given. The second example is a simplified problem of the integrated layout design of a multicomponent system. In this example, we investigate the applicability of the proposed framework to design problems with multimodal objective functions at the upper level. The third example is an industrial application of a permanent magnet synchronous motor (PMSM). In each of these three examples, the external variables represent different objects. In the first example, the external variables represent the positions of the boundary conditions; in the second, these variables indirectly represent the shape of the design domain; and in the third, both.

Design problem
The proposed framework is applied to the simple problem of designing a cantilever. This example is intended to demonstrate the procedure followed by the proposed framework, in particular the process of model selection for metamodeling.
A schematic illustration of the design problem is shown in Fig. 2. The design domain is a 1 × 3 length rectangular area. The cantilever is fixed to two rigid supports on the left side of the design domain, and a downward load is applied on the center of the right side of the design domain. The design variables are the position of support L which is an external variable, i.e., s = L, and the material distribution ρ s in the design domain. The objective is the minimization of the mean compliance. The maximum volume of material is subjected to 15% of the volume of the design domain. L takes a value within the range of 0 ≤ L ≤ 1.4.
In the size optimization problem of a V-shaped truss structure, the theoretical optimum is known. This problem is similar to the discussed example. This optimal structure is a V-shaped structure with a 45 • angle to the wall. In this example, we examine whether the proposed framework can provide a solution similar to the theoretical optimum of the size optimization problem. Note that the solution is expected to be slightly different from the structure with a 45 • angle to the wall because the support is modeled to be a fixed one, and not a pinned one. Implementation for topology optimization at the lower level was accomplished by applying a filter based on a partial differential equation (Kawamoto et al. 2011) to ensure the smoothness of ρ s . The filter radius was set to half of the length of one finite element. The penalization parameter of the SIMP model was set to 3. The design domain was discretized with linear rectangular elements, of which the length is 0.025. The sensitivities were calculated using the adjoint variable method. The optimization problem was solved with SLP, of which the move limit is set to 0.050. The commercial software, COMSOL Multiphisics 5.2a, was used for finite element analysis. The optimization routine was implemented with MATLAB 2016a.
In the upper level, the algorithm for metamodeling was implemented with MATLAB 2016a. The optimization problem was solved with SLP by using the same settings as in the lower level.

Results
The optimization was performed according to the following procedure. First, 10 sample points L (k) (k = 1, 2, ..., 10), 10 equally spaced points from L (1) = 0.1 to L (10) = 1.4, were arranged. Second, 10 topology optimizations were performed by fixing L (k) to each sample point value. Representative results of the topology optimizations are illustrated in Fig. 3. Third, a metamodel H (L) is constructed. Here,H (L) returns an approximated value of mean compliance corresponding to an optimized material distribution under the condition that L is fixed.
Leave-one-out cross-validation was applied to polynomial models of the 2nd to 6th order. We selected one out of 10 sample points and excluded it from the dataset that was used to construct the metamodel. The metamodel was generated using the remaining nine points, and the RMSE between the excluded sample point and the constructed metamodel was measured. Similarly, this operation was carried out for 10 patterns in which all points were selected, and finally, the average RMSE was calculated. This series of operations was performed for polynomials of the 2nd to 6th orders. The results are shown in Fig. 4. The results show that the average RMSE of the 5th order polynomial is the smallest. This is because the representation ability of the model improves as the degree of the polynomial increases (2nd to  Fig. 4), but if it becomes overly large, it causes overfitting (6th in Fig. 4). Therefore, in this example, a 5th order polynomial is used for the metamodel.
The metamodel constructed by the 5th order polynomial with 10 samples is shown in Fig. 3. Then, the external variable L is optimized based on theH (L) using SLP. The design sensitivity ofH (L) is derived analytically. The optimal L onH (L) is 0.8873, and the value ofH (L ) is The result shows the changes in the resulting material distributions and their performance when the value of the external variable is varied. Although there is a slight difference because of the problem settings, the obtained optimal solution is close to the theoretical optimum of the size optimization problem of a V-shaped truss. The optimized structure is physically reasonable as it eases the trade-off between the beam width restricted by the volume constraint and aspect ratio of the support points and load point.

Design problem
The purpose of the second example is to demonstrate the applicability of the proposed framework to a problem of which the upper level objective function H (s) is multimodal. In this design problem, we optimized the configuration of a cantilever with respect to the material distribution of the support structure of the cantilever and the position of a rigid object embedded therein, as illustrated in Fig. 5.
The problem configuration is described below. The design domain is 1 × 1.5 rectangular. A downward load was applied to the lower right part of the design domain (of which the width is 0.2), whereas the left edge remains fixed. The embedded rigid object had a circular shape with a radius of 0.1. Young's modulus of the support structure and the rigid object were set to 1 and 100, respectively. The support structure is represented by material distribution ρ s . The position of the rigid object is represented by the  x, y], which are external variables in this problem. x, y take 0.15 ≤ x ≤ 1.34, 0.14 ≤ y ≤ 0.86, i.e., the center of the rigid object is placed inside the dotted frame in Fig. 5. In this case, x, y are normalized such that their values are in the range from 0 to 1. Let them be x, y; hereafter, an overline denotes that the variable is normalized. The objective function to be minimized is mean compliance, and the design variables are x, y, and ρ s . A volume constraint is applied to ρ s , of which the upper limit is 75% of the design domain. In the implementation of the topology optimization, linear triangular finite elements are used. The average size of the finite elements is 1.06 × 10 −4 . The other settings are the same as those in Example 1.

Results
First, the sample points are arranged to be 121 grid combinations of x = 0, 0.1, ..., 1.0, y = 0, 0.1, ..., 1.0, i.e., [x (k) , y (k) ], (k = 1, 2, ..., 121) are arranged. Second, leave-one-out cross-validation is applied to polynomial approximation metamodels to decide which degrees of the polynomial to use. As a result, a 6th order polynomial model is selected, because the mean RMSE in the validation result of the 2nd to the 7th order models reaches a minimum at the 6th order model. The validation was terminated at the 7th order, because the mean RMSE increased again. Third, a metamodel is constructed. The sample points, the metamodel, its local optima and corresponding optimized material distributions are shown in Fig. 6. The distribution of the sample points in Fig. 6a and b shows that H (x, y) is a multimodal function. Four local solutions s (i) − s (iv) are obtained by optimizations on the metamodel with 25 different initial values using SLP (red points in Fig. 6a and b). Finally, topology optimizations are performed under the condition that values of external variables are fixed to s (i) − s (iv) . The results are shown in Fig. 6c. Among them, the structure corresponding to s (i) = [x (i) , y (i) ] = [0.0000, 1.0000] showed the best performance. In the other solutions, the performance is higher, of the order of s (iv) = [0.0000, 0.0000], s (ii) = [0.7650, 0.4050], s (iii) = [0.9590, 0.0000]. The errors between H (s ) andH (s ) are 0.03747%, 1.162%, 0.2454%, and 0.07264%, for (i), (ii), (iii), and (iv), respectively.
In the optimal structures, the rigid component is located near the fixed end for (i) and (iv), and near the loading surface for (iii). The structure in (ii) is similar to a previously reported result (Qian and Ananthasuresh 2004). The ability of these structures to support the load is reasonable because all the rigid objects in these results are located where the stress is usually high, similar to the considerations described in the previous work (Qian and Ananthasuresh 2004). The results obtained by the proposed framework are similar to the ones obtained by Sample point Local optima of external variables x , y the single-level formulation approach from the perspective of physical meaning.

Design problem
The proposed framework is applied to an industrial application of an 8 pole, 48 slot PMSM rotor design, the specifications of which are based on the Toyota Prius model (Burress et al. 2011). A PMSM consists of a rotor and a stator. Among them, the rotor consists of a yoke and permanent magnets (PMs). The void space of the yoke is known as the flux barrier, which rectifies the magnetic flux. The PMSM is mounted in hybrid or electric vehicles, and therefore it is desired to improve the torque performance. The structural design of the PMSM rotor should consider the positional interaction between the PMs and yoke. Because the PMSM generates torque consisting of a combination of magnetic torque and reluctance torque, the positional relation of the PMs, yoke, and flux barriers greatly influence the torque performance. Therefore, the material distribution of the yoke and arrangement of the PM are optimized simultaneously in this example. Several researchers have also attempted to solve this design problem (Ishikawa et al. 2015;Lim et al. 2015;Kuci et al. 2018;Jung et al. 2018).
The design domain and design variables are shown in Fig. 7. The analysis model of the PMSM is set to a section 1/8th of that of a circle by making use of magnetic symmetry. The design domain of topology optimization is the rotor domain excluding the PM domain, and the outer edge is excluded to prevent the yoke from being separated from the PMs because of centrifugal force. The feasible domain of the PM arrangement is inside the rotor domain.  R, θ, ϕ]. The ranges of [R, θ, ϕ] are normalized within the range from 0 to 1. Let them be [R, θ, ϕ]. The objective function F (s, ρ s ) to be maximized is the average torque of one period. This optimization problem includes three constraints: a rigidity constraint G 1 (s, ρ s ) to prevent separation of the yoke and PMs, a constraint G 2 (ρ s ) on the volume upper limit of the yoke, a constraint on PM arrangement G 3 (s) such that the PM does not intersect or protrude from the design domain. The upper limits of G 1 (s, ρ s ) and G 2 (ρ s ) are set equal to the mean compliance and volume, respectively, of a reference rotor design (a), which is described later.
In the topology optimization of the yoke, magnetic reluctivity ν is given by the following interpolation model (Choi and Yoo 2008) where ν 0 , ν 1 (B), B are the magnetic reluctivity of the void and that of the material and magnetic flux density, respectively. ν 1 (B) is a function of B, of which the behavior depends on the type of material. B is computed by the FEM. p is set to 1 in this study following the setting in Choi and Yoo (2008).
In the implementation, quadratic quadrilateral and triangular finite elements are used to discretize the analysis domain. The design domain for topology optimization is discretized using only the triangular finite elements. The other settings are the same as those in Example 1.

Results
The optimization is performed according to the following procedure. First, 125 sample points s (k) (k = 1, 2, ..., 125) are prepared by filling up the design space under the condition that the PM neither intersects with nor protrudes from the design domain. Second, the topology of the yoke is optimized under the condition that the PM arrangement s (k) is fixed. Third, a metamodelH (s) illustrated in Fig. 8 is constructed. From the result of leave-one-out cross-validation, 4th order polynomial is selected for the metamodel. The areas in Fig. 8  Finally, the topology optimization is performed under the condition of s , and an optimized structure and external variables are obtained. The error between H (s ) andH (s ) is 0.09911%. A motor design constructed based on this result is illustrated in Fig. 9c, hereafter referred to as design (c). Design (c) is constructed by binarizing the threshold value of ρ = 0.5 for comparison with design (a) described later.
We compare the optimized motor design (c) with two designs. Design (a), shown in Fig. 9a, is a reference model built on the basis of a PMSM product (Burress et al. 2011). Design (b), shown in Fig. 9b, is a design in which only the material distribution is optimized and the position of the PM is the same as in design (a). Design (b) is constructed by binarizing as same way in (c). The torque performance of each of the designs is 198.8[N · m], 229.8[N · m], and 232.6[N · m], respectively. The error of performance before and after binarization of designs (b) and (c) is 0.7848% and 0.5465%, respectively. Design (b) is an improvement of 15.60% relative to design (a), the improvement of which Fig. 9 Comparison of the rotor structures. a Reference rotor design, b optimized rotor design where the position of the PM is not considered, and c optimized rotor design where both the material distribution and PM position are considered. In each figure, the figure on the left shows the overall PMSM, and that on the right shows 1/8 of the section of PMSM is achieved only by topology optimization. Design (c) is an improvement of 17.06% compared with design (a), the improvement of which is achieved by incorporating external variables in the topology optimization. In other words, considering the external variables (PM position), the improvement in performance from designs (a) to (c) is 9.39% greater than that from (a) to (b).
Here, we discuss the geometric differences. First, we compare the shape of the yoke of designs (a) and (b). A feature common to these two designs is that flux barriers, which inhibit magnetic flux leakage, are placed between the PMs and the outer narrow sides of them. However, the flux barrier between the PMs is larger in design (b); therefore, the leakage is more effectively suppressed than in design (a). The shape of the yoke around the rotating axis is considered to contribute to satisfying the rigidity constraint G 1 (s, ρ s ). Second, we compare the position of the PM and the shape of the yoke of designs (a), (b), and (c). The position of the PM in design (c) results in R and θ increasing to the maximum within the feasible region as in Fig. 8f. The angle between the PMs becomes steeper by increasing ϕ. This PM arrangement contributes to widen the magnetic path on the outer side of the PMs and suppress the magnetic flux leakage at the outer narrow sides of PMs more effectively. Variation in the position of the PM causes the shape of the yoke to change, thereby placing a smaller flux barrier between PMs.

Discussion
Although the effectiveness of the proposed framework is demonstrated by using three examples, some points need to be taken into consideration when applying the framework to other design problems. The proposed framework inherits the limitations of the underlying methods, namely topology optimization and metamodeling. The requirements for the effective use of the proposed framework are described below.
First, it is necessary to assess whether the configurations optimized in the topology optimization are reasonable from an engineering point of view. This is because of the risk of constructing a metamodel that is meaningless in engineering. A metamodel constructed based on meaningless solutions of topology optimization would not accurately estimate the value of H (s); consequently, the optimized s would be equally meaningless.
Second, the computational costs for metamodeling should be considered. The increase in computational costs associated with an increase in the number of dimensions is known to be an important issue in metamodeling. This issue is often referred to as the "curse of dimensionality" (Simpson et al. 2008;Shan and Wang 2010). Computational cost may be enormous when the proposed framework is applied to a problem with many external variables. Note that, in the proposed framework, the limit in the number of dimensions is considered fewer than that of the metamodels that approximate responses only by conducting forward analysis. This is because the metamodel in the proposed framework approximates responses by solving multiple topology optimization problems.

Conclusion
In this study, we introduced the concept of an external variable to address the design problem of topology optimization by incorporating external variables, where the material distribution and external variables are concurrently optimized. Then, a novel structural optimization framework is presented. The framework employs multi-level formulation by focusing on the dominance-dependence relationship between external variables and material distribution. In addition, metamodeling was utilized to reduce the computational cost and perform sensitivity analysis in the upper level optimization. Leave-one-out cross-validation was used to determine the parameter of the metamodel. The existing density-based topology optimization method was used without modification in the lower level optimization. The effectiveness of the framework was demonstrated by using three examples. The results showed that the performance of the optimized material distribution derived by conventional topology optimization depends on the values of the external variables. Furthermore, we showed that the proposed framework can obtain the optimal solution of the external variables and material distribution. The proposed framework can be applied to optimization problems for which the upper level objective function is expected to be multimodal if metamodels can be constructed with good accuracy.

Compliance with ethical standards
Conflict of interest The authors declare that they have no conflict of interest.

Replication of results
The necessary information for replication of the results are present in the manuscript. The interested reader may contact the corresponding author for further implementation details.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creativecommons. org/licenses/by/4.0/.