An Adaptive Optimization Algorithm Based on Kriging Interpolation with Spherical Model and its Application to Optimal Design of Switched Reluctance Motor

Abstract

In this paper, an adaptive optimization strategy utilizing Kriging model and genetic algorithm is proposed for the optimal design of electromagnetic devices. The ordinary Kriging assisted by the spherical covariance model is used to construct surrogate models. In order to improve the computational efficiency, the adaptive uniform sampling strategy is applied to generate sampling points in design space. Through several iterations and gradual refinement process, the global optimal point can be found by genetic algorithm. The proposed algorithm is validated by application to the optimal design of a switched reluctance motor, where the stator pole face and shape of pole shoe attached to the lateral face of the rotor pole are optimized to reduce the torque ripple.

1. Introduction

The target of inverse problem in electromagnetic field is to implement design optimization of electromagnetic device, so that the required performance or parameters will be satisfied. The stochastic optimization algorithms can search the global optimum. However, these stochastic algorithms cannot avoid drawbacks such as huge computational cost and low rate of convergence. Recently, stochastic global optimization algorithms combined with surrogate performance models have been concerned [1-3]. In these optimization strategies, the objective function is evaluated indirectly by interpolated functions. In general, the approximated function should have low computational cost, high accuracy, and very good interpolation performance. Kriging, a spatial statistical technique with the precise interpolation and prediction features, is popular for the electromagnetic design optimization.

Kriging model predicts the objective value at the unknown point by computing the weighted average of available known samples. There are some key components such as drift function, covariance function, neighbor structure, and variance of interpolation errors. According to the usage of different drift function, Kriging models can be classified into Simple Kriging (SK), Ordinary Kriging (OK), and Universal Kriging (UK) [4]. The SK with a drift function of zero constant is the most basic form. Ordinary Kriging with a nonzero drift function is compatible with a stationary model. Due to its low computational cost and accurate interpolation capability to replace an objective function to assist optimization search [5], the OK is the most popular one among three Kriging models. The UK, the general Kriging model, is a nonstationary geostatistical method, where the drift function is modeled as a general linear function of coordinates.

In this paper, a global optimization strategy employing multiple iterations and gradual refinement is developed, in which the Ordinary Kriging algorithm with spherical covariance model is used. The genetic algorithm is applied for parameter identification of spherical model and for searching the optimum of the approximated function. The proposed optimization method has been verified through analytic function and application to the optimal design of 3-phase 6/4 switched reluctance motor.

 

2. Kriging Methodology

2.1 Review of ordinary Kriging model

Suppose N sample points with corresponding observed values Z(X1),..., Z(XN), and then the function value Z*(X) at the unknown point X can be estimated through a linear combination of the observed values as follows:

where coefficients λi (i=1,..., N) are Kriging weights.

In Kriging models, the best linear unbiased predictor is used to select coefficients λ. Finally, the OK equations are formed by Lagrange multipliers (μ) as follows [6]:

Therefore, the Kriging weights λ can be obtained from (2-a) and (2-b).

2.2. Covariance function

The values of covariance function are not defined yet, but values of some discrete covariance function can be calculated based on sampling points. According to the distance h between every two points, corresponding discrete covariance value C*(h) is obtained as follows:

where N(h) is the number of pair of points with a distance h. The series of h and C*(h) will be used to fit the ideal covariance model.

In this paper, the spherical model will be compared with the thin elastic plates model, and the expressions are shown as follows:

2.2.1 Thin Elastic Plates Model (TEPM):

2.2.2 Spherical Model (SM):

where C0 is nugget, C0+C is sill and a is range [7]. They are unknown parameters, which can be identified by the genetic algorithm (GA).

2.3. Parameter identification of spherical model

Firstly, in GA, the parameters such as initial population size, crossover probability, mutation probability and terminate condition should be defined. Encode the three parameters C0, C, and a with 6 bit binary number for each parameter. Then, each group parameter is a chromosome with 18 lengths. The chromosome encoding method is shown in Fig. 1. Each parameter has 26 codes and decoding expression is given as follows:

where x is decimal value and its range is between Umin and Umax, bi is the i-th binary number between 0 and 1.

Fig. 1.Chromosome encoding

Secondly, according to statistical theory, the objective function is formulated as follows:

where C*(h) is a discrete covariance value and C(h) is a fit covariance value in spherical model.

With the help of the selection, crossover, and mutation operators, the better generation will be selected for next iteration. In this paper, the better individual is chosen by fitness-proportionate selection, and the expression is:

where is the i-th selective probability, m is the population size, fi is fitness value of the i-th individual.

Then, three points crossover algorithm is used with a given crossover probability Pc=0.8 as shown in Fig. 2. In each parameter chromosome fragment, one point is randomly inserted and crossover behind the gene of the point. In addition, the mutation is operated by a given probability Pm=0.04. Finally, if the algorithm satisfies the convergence condition as defined in (9), terminate and output the result. Otherwise, repeat fitness calculation and genetic operation.

Fig. 2.Three-point crossover

In (9), fmax and favg are maximum and average fitness values, respectively. The flowchart of parameter estimation is shown in Fig. 3.

Fig. 3.Flow chart of parameter estimation

 

3. Gradual Refinement Assisted Global Optimization Algorithm

As it is known, for a real electromagnetic problem, the objective function used for optimization is generally implicit and related with finite element performance analysis. Therefore, in order to reduce computational efficiency and improve convergence, during optimization process, the fitness value at unknown design point is approximated by the Kriging model.

In the optimization strategy, the uniform sampling technique is used to obtain sample points in the design space. Based on the optimal result obtained from previous iteration, the design space is reduced and new sampling points of current iteration are gradually inserted to approximate the objective function, so that the efficiency and simulation accuracy will be improved. The proposed optimization algorithm is summarized as follows:

Step 1: Define the initial design space and generate initial sampling points by uniform sampling strategy in the whole design space. Step 2: Carry out performance analysis of each sampling point and obtain corresponding objective values. Step 3: Construct the response surface by OK model. Step 4: Find the current optimal point by GA, and check convergence condition, stop and output the result if converged. Step 5: Reduce design space by adaptive factor of 0.618 around the current optimal point [2]. Step 6: Uniformly generate new sampling points in the reduced design space, and go to Step 2.

In the algorithm, the iteration repeats until the optimal point converges, and the converged global optimal point is considered as a true optimal point. The flow chart is shown in Fig. 4.

Fig. 4.Flow chart of the proposed global optimization algorithm

 

4. Numerical Optimization Results

4.1. Analytic function

A two-dimensional analytic function [8] is used to check the efficiency and validity of the proposed algorithm. The optimization model is formulated as follows:

where the global optimum is (x1=2.287, x2=2.287) and corresponding objective value is 11.81.

Firstly, three global optimal parameters are found by GA. C0=35.127, C=812.698 and a=4.01587. So the spherical model is confirmed.

Fig. 5 shows the constructed response surface by Kriging model and the distribution of sampling points at each iteration. At the initial iteration, 25 sample points are obtained by the uniform sampling in the whole design space, and corresponding response values are calculated. Then the Kriging response surface is constructed, the optimum (2.40476, 2.21429) with F(x1, x2) =11.4789 is found by GA, as shown in Fig. 5(a). Then, the design space is adaptively reduced and is sampled as shown in Fig. 5(b) and Fig. 5(c). After three iterations, a converged optimal point (2.27009, 2.29018) with F(x1, x2) =11.8001 is obtained.

Fig. 5.Optimization process of analytic function

Comparison with average relative error (ARE) for response surface between SM and TEPM, and the optimal solution obtained by direct searching of GA (DS-GA) is shown in Table 1. It is obvious that the SM and DS-GA show a higher accuracy than TEPM, but DS-GA needs much expensive calculating time. For further validation of the model accuracy, we examine two metrics, ARE and root mean squared error (RMSE) for different number of sampling points in SM are defined as follows:

Table 1.a-The results are average values from 100 runs. b-the number of sample points for each iteration is 25. c-DS-GA means the GA is directly applied to (10).

where M is the number of test points, Z*(Xi) and Z(Xi) are estimated and true value of the i-th test point, respectively.

Figs. 6 and 7 show ARE and RMSE for different number of sampling points. From the results, the estimation errors of both the ARE and RMSE decrease as the number of sample points increases, cross-validation becomes more reliable, resulting in a higher accuracy of Kriging interpolation. This verifies Kriging’s claim to be a very flexible for highly nonlinear functions. For the same sampling points, the spherical model with parameters estimation is better, and the optimal solution is closer to the true value.

Fig. 6.ARE of different sampling points

Fig. 7.RMSE of different sampling points

4.2 Electromagnetic application — Shape optimal design of the switched reluctance motor

As an electromagnetic application, a 3-phase 6/4 pole switched reluctance motor (SRM) is optimized for reducing the torque ripple. Each stator pole has a concentrated winding, and each phase consists of two coils wounded on opposite poles and connected in series. The rotor has neither the winding nor the permanent magnet. Because of the doubly salient pole structure of SRM, the torque ripple is produced and causes vibration and noise [9-10].

In SRM, the exciting current and the developed torque are decided by the following equations:

where Vk, Rk, Lk, and ik are separately exciting voltage, resistance, inductance, and exciting current of a phase, respectively. T is the Torque.

From the respective of optimizing the torque through changing the inductance, the air-gap between the stator and rotor poles is selected as the most sensitive parameter. Fig. 8 shows a design parameter θ to control the air-gap when the rotor moves counter-clockwise. Even if the torque ripple from a phase is sufficiently reduced by optimizing the stator pole shape, the SRM still may have a torque ripple. A pole shoe is suggested to be attached to the lateral side of the rotor pole as shown in Fig. 9. Thus, the pole face shape of stator (θ) and the pole shoe of rotor (α) are chosen as two design variables [11].

Fig. 8.Design parameter of the air-gap

Fig. 9.The pole shoe of rotor

To find two optimal parameters by the proposed algorithm, the optimization target is defined as follows:

where Tmax and Tmin are the maximum and minimum Torque, respectively.

In this optimization algorithm, the objective function values of each sampling point are obtained by using finite element analysis (FEA). The three best parameters of the spherical model are obtained by using GA as C0=1.8889E-3, C=2.11111E-3 and a=4.71429. The global optimal point is obtained after three iterations, and through multiple iterations and gradual refinement, the final 75 sample points are used to construct a Kriging response surface. After optimization, the optimal design variables are found θ = 2.90298°, α = 4.38571°, and the corresponding torque ripple is reduced to 3.8724 N·m from the initial one of 10.90 N·m while the average torque is increased from 5.1451 N·m to 5.2259 N·m. Fig. 10 shows the comparison of the torque ripples between the initial and the optimized pole shapes.

Fig. 10.Comparison of the torque between initial model and optimized one

 

5. Conclusion

In this paper, a global optimization strategy employing multiple iterations and gradual refinement is proposed. The OK model with spherical covariance model and TEPM is used as interpolation function to approximate the objective function. Then GA can successfully estimate three parameters of spherical model and find the optimal point based on the surrogate model. Through the applications to numerical examples, the proposed OK with spherical model is proven to give a better optimal solution than that with the thin elastic plate model, and speed up the optimization algorithm to search the optimal design than the direct searching of GA. Our proposal to obtain Ordinary Kriging model for reducing computational cost opens new possibilities in electromagnetic fields to investigate optimization problems under parameter unconstraint. In multidimensional problems and sensitivity analysis of optimization problems, it should be solved by the proposed optimization algorithm. Our future research will be in this direction.