Optimal Design of PMSM Based on Automated Finite Element Analysis and Metamodeling

To obtain accurate optimal design results in electric machines, the finite element analysis (FEA) technique should be used; however, it is time-consuming. In addition, when the design of experiments (DOE) is conducted in the optimal design process, mechanical design, analysis, and post process must be performed for each design point, which requires a significant amount of design cost and time. This study proposes an automated DOE procedure through linkage between an FEA program and optimal design program to perform DOE easily and accurately. Parametric modeling was developed for the FEA model for automation, the files required for automation were generated using the macro function, and the interface between the FEA and optimal design program was established. Shape optimization was performed on permanent magnet synchronous motors (PMSMs) for small electric vehicles to maximize torque while maintaining efficiency, torque ripple, and total harmonic distortion of the back EMF using the built-in automation program. Fifty FEAs were performed for the experimental points selected by optimal Latin hypercube design and their results were analyzed by screening. Eleven metamodels were created for each output variable using the DOE results and root mean squared error tests were conducted to evaluate the predictive performance of the metamodels. The optimization design based on metamodels was conducted using the hybrid metaheuristic algorithm to determine the global optimum. The optimum design results showed that the average torque was improved by 2.5% in comparison to the initial model, while satisfying all constraints. Finally, the optimal design results were verified by FEA. Consequently, it was found that the proposed optimal design method can be useful for improving the performance of PMSM as well as reducing design cost and time.


Introduction
The necessity of eco-friendly vehicles has been highlighted owing to environmental pollution and depletion of fossil fuels. Global electric car stocks are growing rapidly, crossing the 3 million vehicle threshold in 2017. The estimated demand for electric vehicles by 2030 is 100 to 140 million. The core of an electric vehicle is the electric powertrain, which consists of a traction motor, a reduction drive, an inverter, and a power delivery module. Permanent magnet synchronous motors (PMSMs) have been mainly used as a traction motor for electric vehicles because they have high efficiency and high output power density characteristics.
Several studies have been conducted on the PMSMs used in electric vehicles that require various characteristics such as torque, efficiency, and harmonic distortion (THD). Optimal design is essential to satisfy the various design requirements of PMSM at once. Optimal design is a method of finding the values of design variables to obtain an optimal solution within a range of constraints. The optimal design for PMSMs is created by combining design methods such as the analytical model, magnetic equivalent circuits (MEC) model, and finite element analysis (FEA) with optimal design algorithms [1][2][3][4][5][6][7][8]. First of proposed procedure is based on commercial tools, it has a ripple effect that can easily apply optimal design in academia and industry.
Most of the previous studies have been applied to optimal design using metamodel generated in one way. There have been a lot of optimizations recently using a single metamodeling technique such as Kriging and response surface method [10][11][12][13][14][15][16][17][18]. However, since a suitable metamodel is different according to each design problem and condition, it is necessary to select the best metamodel through accuracy evaluation after generating several metamodels. This is because the accuracy of the metamodel must be high to obtain good optimal design results. In this study, metamodels of objective functions and constraints are generated in 11 ways, and the most accurate metamodels are selected through the root mean squared error (RMSE) test, respectively.
In this study, shape optimization is performed for a PMSM to maximize the torque while maintaining efficiency, torque ripple and THD in the back electromotive force (EMF). First, the design target specification of a PMSM for small electric vehicles is established, and the characteristics of the initial model are analyzed using FEA. To improve the accuracy of the design results, DOE is performed using FEA. After the creation of metamodels using the DOE results, the optimal values are obtained by the optimal algorithm. The optimal Latin hypercube design (OLHD) technique [19] is applied for the DOE, and the appropriate DOE number and test point number are selected to produce accurate metamodels. To perform DOE easily and accurately, this study proposes an automated DOE procedure through linkage between an FEA program and an optimal design program. Using the DOE results, the relationship between the design and output variables are analyzed by screening. To generate an accurate metamodel, the RMSE tests are performed on eleven metamodels for each output variable, and the best metamodels are selected for each output variable. Optimization based on metamodels is performed, and the global optimization algorithm hybrid metaheuristic algorithm (HMA) [20] is utilized as the optimal algorithm. The overall process of this study is represented in a flowchart, as shown in Figure 1. the proposed procedure is based on commercial tools, it has a ripple effect that can easily apply optimal design in academia and industry.
Most of the previous studies have been applied to optimal design using metamodel generated in one way. There have been a lot of optimizations recently using a single metamodeling technique such as Kriging and response surface method [10][11][12][13][14][15][16][17][18]. However, since a suitable metamodel is different according to each design problem and condition, it is necessary to select the best metamodel through accuracy evaluation after generating several metamodels. This is because the accuracy of the metamodel must be high to obtain good optimal design results. In this study, metamodels of objective functions and constraints are generated in 11 ways, and the most accurate metamodels are selected through the root mean squared error (RMSE) test, respectively.
In this study, shape optimization is performed for a PMSM to maximize the torque while maintaining efficiency, torque ripple and THD in the back electromotive force (EMF). First, the design target specification of a PMSM for small electric vehicles is established, and the characteristics of the initial model are analyzed using FEA. To improve the accuracy of the design results, DOE is performed using FEA. After the creation of metamodels using the DOE results, the optimal values are obtained by the optimal algorithm. The optimal Latin hypercube design (OLHD) technique [19] is applied for the DOE, and the appropriate DOE number and test point number are selected to produce accurate metamodels. To perform DOE easily and accurately, this study proposes an automated DOE procedure through linkage between an FEA program and an optimal design program. Using the DOE results, the relationship between the design and output variables are analyzed by screening. To generate an accurate metamodel, the RMSE tests are performed on eleven metamodels for each output variable, and the best metamodels are selected for each output variable. Optimization based on metamodels is performed, and the global optimization algorithm hybrid metaheuristic algorithm (HMA) [20] is utilized as the optimal algorithm. The overall process of this study is represented in a flowchart, as shown in Figure 1.

Finite Element Analysis
The target specifications are determined by referring to the Renault's Twizy with a torque of 57 N•m and an output power of 13 kW at 2100 rpm. In this study, a PMSM is selected as the design model, and the target output power is 15 kW which should satisfy 60 N•m at 2387 rpm.

Finite Element Analysis
The target specifications are determined by referring to the Renault's Twizy with a torque of 57 N·m and an output power of 13 kW at 2100 rpm. In this study, a PMSM is selected as the design model, and the target output power is 15 kW which should satisfy 60 N·m at 2387 rpm. Figure 2 and Table 1 show the structure and specifications of a 15 kW PMSM for a small electric vehicle, respectively. The PMSM has 8 poles, 36 slots, and distributed winding.  Table 1 show the structure and specifications of a 15 kW PMSM for a small electric vehicle, respectively. The PMSM has 8 poles, 36 slots, and distributed winding.  The electromagnetic, mechanical and thermal properties of the 35PN210 core material are shown in Table 2. Core loss is the sum of hysteresis loss, eddy-current loss, and excess loss, and is calculated by Equation (1). The core loss varies with frequency, but the analysis is based on 60 Hz.

Initial Model
Where is the core loss, is the hysteresis loss coefficient, is the eddy current loss coefficient, is the excess loss coefficient, is the frequency, and is the amplitude of the alternating flux component.
The electromagnetic, mechanical and thermal properties of permanent magnets are shown in Table 3. V-shaped N38UH grade NdFeB are applied to concentrate the magnetic flux. The magnetic flux density and coercivity decreased with increasing temperature, but the analysis is conducted at 20 °C.  The electromagnetic, mechanical and thermal properties of the 35PN210 core material are shown in Table 2. Core loss is the sum of hysteresis loss, eddy-current loss, and excess loss, and is calculated by Equation (1). The core loss varies with frequency, but the analysis is based on 60 Hz.
where P c is the core loss, K h is the hysteresis loss coefficient, K c is the eddy current loss coefficient, K e is the excess loss coefficient, f is the frequency, and B m is the amplitude of the alternating flux component. The electromagnetic, mechanical and thermal properties of permanent magnets are shown in Table 3. V-shaped N38UH grade NdFeB are applied to concentrate the magnetic flux. The magnetic flux density and coercivity decreased with increasing temperature, but the analysis is conducted at 20 • C.

No Load Analysis
Characteristic analysis of the initial model is performed by FEA under the no load condition without current excitation. When the rotor of the initial model rotates at the rated speed, back EMF is induced in the stator winding. Since the back EMF simulation is performed while rotating the rotor under the no load condition, the equivalent circuit when the PMSM operates as a generator should be considered, as shown in Figure 3. The voltage equation of the equivalent circuit is shown in Equation (2). However, since no current flows in the armature winding under the no load condition, the terminal voltage and the no load EMF are the same. The back EMF of phase A can be obtained by Equation (3), and the analysis result by FEA is illustrated in Figure 4a.
V is the terminal voltage, . E 0 is the no load EMF, X a is the armature reation reactance, X l is the leakage reactance, R a is the armature resistance and I a . is the load current.
where E a . is the back EMF of phase A, N is the number of the amarture turns, ∅ is the magnetic flux and ω is the electrical angular velocity. THD is an important factor in the electrical equipment and power systems. THD can be obtained by adding the harmonic components to the fundamental wave components of voltage or current as shown in Equation (4) [21]. A higher THD increases the core loss in electric machines, which can reduce the efficiency and generate excessive heat. The harmonic analysis result of the back EMF waveform is shown in Figure 4b, and the THD of the back EMF calculated by Equation (4) is 3.52%.
where V THD is the THD of the back EMF, V 1 is the RMS voltage of the fundamental frequency and V n is the RMS voltage of nth harmonic.
where is the terminal voltage, is the no load EMF, is the armature reation reactance, is the leakage reactance, is the armature resistance and is the load current.
where is the back EMF of phase A, is the number of the amarture turns, ∅ is the magnetic flux and is the electrical angular velocity. THD is an important factor in the electrical equipment and power systems. THD can be obtained by adding the harmonic components to the fundamental wave components of voltage or current as shown in Equation (4) [21]. A higher THD increases the core loss in electric machines, which can reduce the efficiency and generate excessive heat. The harmonic analysis result of the back EMF waveform is shown in Figure 4b, and the THD of the back EMF calculated by Equation (4) where is the THD of the back EMF, is the RMS voltage of the fundamental frequency and is the RMS voltage of nth harmonic.
where is the terminal voltage, is the no load EMF, is the armature reation reactance, is the leakage reactance, is the armature resistance and is the load current.
where is the back EMF of phase A, is the number of the amarture turns, ∅ is the magnetic flux and is the electrical angular velocity. THD is an important factor in the electrical equipment and power systems. THD can be obtained by adding the harmonic components to the fundamental wave components of voltage or current as shown in Equation (4) [21]. A higher THD increases the core loss in electric machines, which can reduce the efficiency and generate excessive heat. The harmonic analysis result of the back EMF waveform is shown in Figure 4b, and the THD of the back EMF calculated by Equation (4) where is the THD of the back EMF, is the RMS voltage of the fundamental frequency and is the RMS voltage of nth harmonic.

Load Analysis
Target specifications are 15 kW and 60 N•m at 2387 rpm as shown in Table 1. To perform load analysis, a current condition satisfying the target torque of 60 N•m at 2387 rpm should be found. Through static torque analysis, the current condition is determined as the RMS value of 146 A and phase angle of 25°. As shown in Figure 5a, the average torque is 59.95 N•m and the torque ripple is 5.09% of torque. The core loss is interpreted as shown in Figure 5b, and total losses are the sum of the

Load Analysis
Target specifications are 15 kW and 60 N·m at 2387 rpm as shown in Table 1. To perform load analysis, a current condition satisfying the target torque of 60 N·m at 2387 rpm should be found. Through static torque analysis, the current condition is determined as the RMS value of 146 A and phase angle of 25 • . As shown in Figure 5a, the average torque is 59.95 N·m and the torque ripple is 5.09% of torque. The core loss is interpreted as shown in Figure 5b, and total losses are the sum of the core loss and copper loss. The output power of the motor is calculated as the product of torque and angular velocity, and the efficiency can be calculated from the output power and total losses of the motor. The efficiency of the initial model is calculated to be 91.42%. core loss and copper loss. The output power of the motor is calculated as the product of torque and angular velocity, and the efficiency can be calculated from the output power and total losses of the motor. The efficiency of the initial model is calculated to be 91.42%.

Design Process
The design optimization process for maximizing the average torque, while maintaining THD of the back EMF, efficiency, and torque ripple, is shown in Figure 6. The objective function, constraints, thermal condition, and design variables are established, as described in Equations (5)-(9) and Figure  7. To improve the average torque, the average torque is set as both an objective function and a constraint, with a goal of 2% improvement over the initial model. THD of the back EMF and the efficiency are set as constraints to maintain the same level as the initial model. The torque ripple is set below 10%, which is an acceptable level as a traction motor for electric vehicles [22]. To improve the accuracy of the design results, DOE is performed using FEA. Because DOE using FEA requires a significant amount of time and effort, interworking is conducted between the FEA and optimal design programs, which are ANSYS Maxwell and PIAnO, respectively, to automatically perform DOE. When the analysis and extraction of results for one experiment is finished, the values of the design variables are automatically changed to perform the FEA at the next DOE point. From the DOE results, sensitivity analysis between design variables and output variables is conducted using screening, and each metamodel for output variables is generated. RMSE test was conducted to evaluate the predictive performance of the metamodels, and the best metamodel is selected for each output variable. Based on the selected metamodels, the optimal values are obtained using the HMA.

Objective function
Maximize the average torque

Design Process
The design optimization process for maximizing the average torque, while maintaining THD of the back EMF, efficiency, and torque ripple, is shown in Figure 6. The objective function, constraints, thermal condition, and design variables are established, as described in Equations (5)-(9) and Figure 7. To improve the average torque, the average torque is set as both an objective function and a constraint, with a goal of 2% improvement over the initial model. THD of the back EMF and the efficiency are set as constraints to maintain the same level as the initial model. The torque ripple is set below 10%, which is an acceptable level as a traction motor for electric vehicles [22]. To improve the accuracy of the design results, DOE is performed using FEA. Because DOE using FEA requires a significant amount of time and effort, interworking is conducted between the FEA and optimal design programs, which are ANSYS Maxwell and PIAnO, respectively, to automatically perform DOE. When the analysis and extraction of results for one experiment is finished, the values of the design variables are automatically changed to perform the FEA at the next DOE point. From the DOE results, sensitivity analysis between design variables and output variables is conducted using screening, and each metamodel for output variables is generated. RMSE test was conducted to evaluate the predictive performance of the metamodels, and the best metamodel is selected for each output variable. Based on the selected metamodels, the optimal values are obtained using the HMA.

Objective function
Maximize the average torque Thermal condition Reference temperature = 20 C • (8)  The mechanical constraints of the barrier length are set from a range that facilitates the flow of magnetic flux to a range that inhibits the flow of magnetic flux very much, as shown in Figure 8a. Rib thickness is set to be at least 1 mm in consideration of workability at manufacture and mechanical rigidity at high speed, as shown in Figure 8b. Figure 8c shows the mechanical constraints of the teeth width, and the range is set so that the slot is smaller than the initial model and maintains the proper fill factor. The range of teeth thickness is set up to reduce the saturation of the magnetic flux at the tip of the teeth and maintain the proper fill factor, as shown in Figure 8d. Barrier gap size affects the  The mechanical constraints of the barrier length are set from a range that facilitates the flow of magnetic flux to a range that inhibits the flow of magnetic flux very much, as shown in Figure 8a. Rib thickness is set to be at least 1 mm in consideration of workability at manufacture and mechanical rigidity at high speed, as shown in Figure 8b. Figure 8c shows the mechanical constraints of the teeth width, and the range is set so that the slot is smaller than the initial model and maintains the proper fill factor. The range of teeth thickness is set up to reduce the saturation of the magnetic flux at the tip of the teeth and maintain the proper fill factor, as shown in Figure 8d. Barrier gap size affects the The mechanical constraints of the barrier length are set from a range that facilitates the flow of magnetic flux to a range that inhibits the flow of magnetic flux very much, as shown in Figure 8a. Rib thickness is set to be at least 1 mm in consideration of workability at manufacture and mechanical rigidity at high speed, as shown in Figure 8b. Figure 8c shows the mechanical constraints of the teeth width, and the range is set so that the slot is smaller than the initial model and maintains the proper fill factor. The range of teeth thickness is set up to reduce the saturation of the magnetic flux at the tip of the teeth and maintain the proper fill factor, as shown in Figure 8d. Barrier gap size affects the formation of magnetic flux and the motor performance since the permanent magnet position also changes. Therefore, its mechanical constraints are set as shown in Figure 8e.

Automated DOE Procedure
To perform DOE, the shape of design variable should be changed. However, when DOE is processed manually, the shape of each model is drawn using the CAD tool. Next, the designed shape should be imported into the FEA program and the FEA should be performed for each model. After FEA, the post process is required to calculate the desired result. As manual DOE requires a lot of effort and time, this study suggests the automation of the DOE process. First, the Maxwell's parametric sweep setup function is used to change the shape of design variables without using the CAD tool. The use of this function can change the shape of the FEA model by inputting numerical values in the Maxwell program. Next, Maxwell's Macro function is used to perform DOE using PIAnO, an optimal design program. As the design variables change, the vbscript and batch files are created to automatically change the shape of the FEA model. In addition, vbscript and batch files are generated for FEA under the no load and load conditions for each experiment. Vbscript and batch files also are generated to output and quantify the torque, core loss and FFT analysis results of the back EMF obtained through FEA. The files created through the Macro function are shown in Figure  9. Next, the interface configuration between Maxwell and PIAnO for automation is shown in Figure  10. To process DOE in PIAnO, the files created in Maxwell are imported as shown in Figure 10a, and the script for calculating the output variables is shown in Figure 10b.

Automated DOE Procedure
To perform DOE, the shape of design variable should be changed. However, when DOE is processed manually, the shape of each model is drawn using the CAD tool. Next, the designed shape should be imported into the FEA program and the FEA should be performed for each model. After FEA, the post process is required to calculate the desired result. As manual DOE requires a lot of effort and time, this study suggests the automation of the DOE process. First, the Maxwell's parametric sweep setup function is used to change the shape of design variables without using the CAD tool. The use of this function can change the shape of the FEA model by inputting numerical values in the Maxwell program. Next, Maxwell's Macro function is used to perform DOE using PIAnO, an optimal design program. As the design variables change, the vbscript and batch files are created to automatically change the shape of the FEA model. In addition, vbscript and batch files are generated for FEA under the no load and load conditions for each experiment. Vbscript and batch files also are generated to output and quantify the torque, core loss and FFT analysis results of the back EMF obtained through FEA. The files created through the Macro function are shown in Figure 9. Next, the interface configuration between Maxwell and PIAnO for automation is shown in Figure 10. To process DOE in PIAnO, the files created in Maxwell are imported as shown in Figure 10a, and the script for calculating the output variables is shown in Figure 10b.

Design of Experiment
The number of experiments and the number of test points are determined in three steps [23]. First, the number of experiments should be selected according to the number of design variables. When the number of design variables is ten or less, the number of experimental points is determined by Equation (10), nEXP 1.5 nSAT 1.5 nDV 1 nDV 2 2 (10) where nEXP is the number of DOE, nSAT is the number of saturation, and nDV is the number of design variables.
Next, the number of DOEs that can be used as test points for evaluating the accuracy of the metamodel should be secured by Equation (11): nEXP min nDV 1 nDV 2 2 , 10 nDV 5 nDV (11) Because five design variables are used in this study, the number of DOE should be more than 46 by Equations (10) and (11). Therefore, the number of DOE is determined to be 50, which is a multiple of the design variables. If fifty experiments are manually operated, a significant amount of effort and time would be required. However, in this study, automation is implemented so that DOE can be easily developed and design cost can be reduced. Finally, the number of test points for evaluating the accuracy of the metamodel is determined to be five by Equation (12), nEXP_ts min nEXP 10%, 10 nDV (12)

Design of Experiment
The number of experiments and the number of test points are determined in three steps [23]. First, the number of experiments should be selected according to the number of design variables. When the number of design variables is ten or less, the number of experimental points is determined by Equation (10), where nEXP is the number of DOE, nSAT is the number of saturation, and nDV is the number of design variables. Next, the number of DOEs that can be used as test points for evaluating the accuracy of the metamodel should be secured by Equation (11): Because five design variables are used in this study, the number of DOE should be more than 46 by Equations (10) and (11). Therefore, the number of DOE is determined to be 50, which is a multiple of the design variables. If fifty experiments are manually operated, a significant amount of effort and time would be required. However, in this study, automation is implemented so that DOE can be easily developed and design cost can be reduced. Finally, the number of test points for evaluating the accuracy of the metamodel is determined to be five by Equation (12), nEXP_ts > min[nEXP × 10%, 10 × nDV] (12) where nEXP_ts is the number of test points. The OLHD technique is applied to determine the sampling point of the DOE. OLHD is a type of DACE sampling technique developed for computational experiments. In the computer experiment, because there are no random errors, only the bias error should be considered and the test point should be spread evenly inside the design area. OLHD improves the space-filling property by using the optimum conditions and spreads the test points evenly; thus, even if there are several test points, they can be selected efficiently. DOE for fifty test points selected by OLHD is easily performed using the automated program. Sensitivity analysis is conducted to analyze the correlation between the design variables and design results. Figure 11 shows that the barrier length has the highest impact on the output variables among the five design variables. However, as shown in Table 4, even the most optimal experimental point among the 50 experiments does not satisfy the constraints. Therefore, metamodeling based on DOE results is conducted. where nEXP_ts is the number of test points. The OLHD technique is applied to determine the sampling point of the DOE. OLHD is a type of DACE sampling technique developed for computational experiments. In the computer experiment, because there are no random errors, only the bias error should be considered and the test point should be spread evenly inside the design area. OLHD improves the space-filling property by using the optimum conditions and spreads the test points evenly; thus, even if there are several test points, they can be selected efficiently. DOE for fifty test points selected by OLHD is easily performed using the automated program. Sensitivity analysis is conducted to analyze the correlation between the design variables and design results. Figure 11 shows that the barrier length has the highest impact on the output variables among the five design variables. However, as shown in Table 4, even the most optimal experimental point among the 50 experiments does not satisfy the constraints. Therefore, metamodeling based on DOE results is conducted. Figure 11. Sensitivity analysis using screening.

Metamodeling
Five test points are selected to evaluate the metamodel, and eleven metamodels are generated for each output variable. The metamodel can be classified into a regression model and an interpolation model. The regression model, i.e., polynomial regression (PR), radial basis function regression (RBFr), ensemble of decision trees (EDT), and multi-layer perceptron (MLP), smoothens the noise data because they do not pass through the test points exactly. Therefore, this model is useful for real experiments with random errors. PR allows free choice of regression terms [24]. RBF is easy

Metamodeling
Five test points are selected to evaluate the metamodel, and eleven metamodels are generated for each output variable. The metamodel can be classified into a regression model and an interpolation model. The regression model, i.e., polynomial regression (PR), radial basis function regression (RBFr), ensemble of decision trees (EDT), and multi-layer perceptron (MLP), smoothens the noise data because they do not pass through the test points exactly. Therefore, this model is useful for real experiments with random errors. PR allows free choice of regression terms [24]. RBF is easy to design and generalize, and has strong tolerance to input noise [25]. EDT is advantageous for expressing nonlinearity in large amounts of data. MLP is type of deep learning algorithm and has the advantage of being able to represent the nonlinear relationships between input and output variables [26]. In contrast, the interpolation model, i.e., Kriging and radial basis function interpolation (RBFi), is well suited for function approximation using analytical results without random errors because it passes through the test points exactly. The estimated equation of the Kriging model was defined to eliminate bias and thereby minimize error variance [27]. Thus, a numerically robust model is provided. RBFi was first popularized in the machine learning community and has been used in computer graphics [28].
The accuracy of the metamodel is a very important factor in the optimal design using metamodel [12]. This is because the predictive performance of the metamodel affects the reliability of the optimal design. Most of the existing studies have been metamodeled by a single method such as Kriging and RSM, and the accuracy evaluation has not been performed [10][11][12][13][14][15][16][17][18]. In this study, however, metamodels for the objective function and constraints are generated in 11 ways provided by PIAnO, and the best metamodels are selected, respectively, by comparing the RMSE test results to evaluate the metamodel accuracy. The predictive performance of the metamodel is evaluated by the RMSE test and is calculated by Equation (13) [23], where y(X i ) is the value of the real function andŷ(X i ) is the value of the metamodel. Through the RMSE test, the predictive performances of the metamodels are evaluated for the output variables. The RMSE test showed the best predictive performance of RBFr as a metamodel of the average torque as shown in Table 5. Similarly, the RMSE tests are conducted on the metamodel for efficiency, torque ripple, and THD of the back EMF. Based on the test results, the metamodels with the best predictive performance for each output variable are selected for use in the optimal design, as shown in Table 6.

Design Optimization Based on Metamodel
The HMA, a global optimization algorithm, is used for the optimal design based on the metamodel. The HMA was proposed in 2016 by Park [20]. HMA can determine the global optimum faster than other global optimizers owing to the combined advantages of improved constrained differential evolution and modified cuckoo search.
The optimum design results predicted from the metamodel based HMA are shown in Table 7 and verified through FEA. The predicted results showed that the average torque, THD of the back EMF, efficiency and torque ripple results are similar to the FEA results. Therefore, the automated DOE procedure and the generation and evaluation of the metamodel were verified. The average torque of the optimal model was 2.5% better than the initial model, and the torque ripple increased slightly, as shown in Figure 12. THD of the back EMF and efficiency set by the constraints were slightly improved. Although the torque ripple of the optimal model is 7.822%, it is very acceptable as a traction motor for electric vehicles [22].

Design Optimization Based on Metamodel
The HMA, a global optimization algorithm, is used for the optimal design based on the metamodel. The HMA was proposed in 2016 by Park [20]. HMA can determine the global optimum faster than other global optimizers owing to the combined advantages of improved constrained differential evolution and modified cuckoo search.
The optimum design results predicted from the metamodel based HMA are shown in Table 7 and verified through FEA. The predicted results showed that the average torque, THD of the back EMF, efficiency and torque ripple results are similar to the FEA results. Therefore, the automated DOE procedure and the generation and evaluation of the metamodel were verified. The average torque of the optimal model was 2.5% better than the initial model, and the torque ripple increased slightly, as shown in Figure 12. THD of the back EMF and efficiency set by the constraints were slightly improved. Although the torque ripple of the optimal model is 7.822%, it is very acceptable as a traction motor for electric vehicles [22].   Figures 13 and 14 show the flux distributfion and flux density of the initial and optimal models, respectively. In comparison to the initial model, the rib thickness and barrier length of the optimal model were reduced, and the barrier gap was increased. As the rib thickness decreased, the unnecessary flux flow between the north pole and south pole through the rotor rib was reduced. In addition, the flux flow was smoothly improved owing to the reduction in barrier length, and consequently, more flux passed through the stator core. The improvements in the magnetic flux flow  Figures 13 and 14 show the flux distributfion and flux density of the initial and optimal models, respectively. In comparison to the initial model, the rib thickness and barrier length of the optimal model were reduced, and the barrier gap was increased. As the rib thickness decreased, the unnecessary flux flow between the north pole and south pole through the rotor rib was reduced. In addition, the flux flow was smoothly improved owing to the reduction in barrier length, and consequently, more flux passed through the stator core. The improvements in the magnetic flux flow and the change in reluctance can be considered to be the cause of the increase in the back EMF and torque [29]. Owing to the improvements of in the flux flow, the back EMF of the optimum model was 35.0 V, which showed an improvement of 7.4% in comparison to 32.6 V of the initial model, as shown in Figure 15. In addition, owing to the sinusoidal improvement in the waveform of the back EMF of the optimal model, the THD was slightly improved from 3.414% to 3.065%. torque [29]. Owing to the improvements of in the flux flow, the back EMF of the optimum model was 35.0 V, which showed an improvement of 7.4% in comparison to 32.6 V of the initial model, as shown in Figure 15. In addition, owing to the sinusoidal improvement in the waveform of the back EMF of the optimal model, the THD was slightly improved from 3.414% to 3.065%.  torque [29]. Owing to the improvements of in the flux flow, the back EMF of the optimum model was 35.0 V, which showed an improvement of 7.4% in comparison to 32.6 V of the initial model, as shown in Figure 15. In addition, owing to the sinusoidal improvement in the waveform of the back EMF of the optimal model, the THD was slightly improved from 3.414% to 3.065%.

Conclusions
This paper presented shape optimization of a PMSM for small electric vehicles to maximize torque while maintaining efficiency, torque ripple and THD of the back EMF. To improve the accuracy of the optimal design results, DOE was performed using FEA. This study proposed an automated DOE procedure through linkage between an FEA and optimal design programs to perform DOE easily and accurately. Parametric modeling was performed for the FEA model to change the shape variables automatically, and automation-related files were created using Maxwell's Macro function. In addition, an interface was established to link the FEA program with PIAnO, an optimal design program. Using the built-in automation program, 50 FEAs for the experimental points selected by OLHD were easily performed. From the DOE results, the relationship between the design and output variables was analyzed by screening. Among the five design variables, the barrier length was found to have the greatest effect on the output variables. Eleven metamodels were created for each output variable and RMSE test was conducted to evaluate the predictive performance of the metamodels. Consequently, the metamodels with the best predictive performance for each output variable were selected. Finally, the optimization design based on the metamodel was determined using the HMA to find the global optimum. The objective average torque improved by 2.5% over the initial model while satisfying all the constraints. The optimal design results were finally verified by FEA.
The proposed automated design method is expected to reduce design cost and time. Moreover, it can be used to find the optimal solution for various design problems as well as PMSMs. By following the procedure given below, the proposed optimal design method can be applied to any type of motor without any special constraints. First, in order to change the shape of the optimum design variable automatically, the dimension of the optimal design variable should be set using Maxwell's parametric sweep setup function. Next, determine the values that you want to extract from Maxwell and create vbscript and batch files to extract them. Finally, an interface setting must be performed to accommodate Maxwell's output values in an optimization program called PIAnO.
Optimization of multi-physics systems by simulation takes significant computing time for each simulation run, and its process depends on numerous runs, making it difficult and expensive [30]. However, using the automated DOE procedure suggested in this study can reduce design cost and time, so I think multi-physics analysis is possible in the near future. In the next project, I will consider multi-physics analysis that takes into account the mechanical and thermal properties.