Average-torque-maximization and cogging-torque-minimization of permanent-magnet-assisted synchronous reluctance motor using topology optimization

Abstract

Although rare earth materials are actively utilized in permanent magnets for a wide variety of applications including electric vehicles and trains, the cost of these materials is on the rise. Therefore, the synchronous reluctance motor (SynRM) which does not require permanent magnets is remarkable from the perspective of high efficiency and low cost. However, the torque characteristic of the SynRM is much lower than that of the interior permanent magnet synchronous reluctance motor (IPMSM). Then, a permanent magnet assisted synchronous reluctance motor (PMASynRM) which bonded magnet was attached to the rotor structure of SynRM has been proposed. The manufacturing cost of the PMASynRM is lower than that of the IPMSM because a bounded magnet is utilized former. In this work, the structure of the PMASynRM has been drastically changed by topology optimization to average torque maximization and cogging torque minimization.

Keywords

Average torque cogging torque Heaviside function permanent magnet assisted synchronous reluctance motor topology optimization

1. Introduction

The interior permanent magnet synchronous reluctance motor (IPMSM) using a neodymium magnet made from rare earth materials is currently used in electric vehicles, trains and other applications. However, because the cost of rare earth material is continuously increasing, the synchronous reluctance motor (SynRM), which does not require a permanent magnet, is of great interest from the perspective of high efficiency. Specifically, the efficiency of SynRM has been improved beyond the induction motor [1]. However, the torque characteristic of SynRMs is much lower than that of the IPMSM. To improve the torque characteristic, the permanent magnet assisted synchronous reluctance motor (PMASynRM), in which a bonded magnet [2,3] is attached to the rotor structure of the SynRM, is proposed in [4].

In conventional SynRM and PMASynRM research, the rotor structure design has been attempted such as raise the torque [5] and torque ripple reduction structure [6]. However, topology optimization of PMASynRM has not been reported to date. In this investigation, the aim is to design the rotor structure of the PMASynRM to improve its torque characteristic using topology optimization. Because topology optimization has a high flexibility of design, the rotor structure is drastically changed.

In this report, topology optimization based on the Heaviside function aimed at average torque maximization and cogging torque minimization was applied to the rotor structure of PMASynRM. Suppression of motor noise and vibration is expected with a reduction in the cogging torque. Coefficients are provided for each objective functions, the optimization which increasing the coefficient for the cogging torque minimization is investigated. In the optimized structure, the torque average and cogging torque are improved.

2. Method of topology optimization

2.1. Finite element method

The governing equation is formulated in magnetostatic field, which can be written as follows: $\begin{eqnarray}\frac{\partial }{\partial x}\left({\nu}\frac{\partial A}{\partial x}\right)+\frac{\partial }{\partial y}\left({\nu}\frac{\partial A}{\partial y}\right)+J_{0}(t)=0,\end{eqnarray}$ (1) where A is a magnetic vector potential of the z-axis direction, J ₀(t) is the magnetizing current density at time t, 𝜈 is the magnetic resistivity. Equation ((1)) is changed to the weak form G _i by using the Galerkin method as follows: $\begin{eqnarray}G_{i}=\iint _{S}\left\{\frac{\partial N_{i}}{\partial x}\left({\nu}\frac{\partial A}{\partial x}\right)+\frac{\partial N_{i}}{\partial y}\left({\nu}\frac{\partial A}{\partial y}\right)\right\}dxdy-\iint N_{i}J_{0}(t)=0,\end{eqnarray}$ (2) where N _i is the shape function. The linear equation is defined as follows: $\begin{eqnarray}K\boldsymbol{A}=\boldsymbol{F},\end{eqnarray}$ (3) where K is the stiffness matrix, F is the right side vector.

Fig. 1.

Topology modelling based on Heaviside function.

2.2. Heaviside function

Figure 1 shows the topology optimization based on the Heaviside function. When the smoothed Heaviside function H (𝜓) [7] is applied to the approximated characteristic function, the reluctivity 𝜈( B ², 𝜓) in the design domain is formulated as follows: $\begin{eqnarray}{\nu}(\boldsymbol{B}^{2},{\psi})=\{1-H({\psi})\}{\nu}_{0}+H({\psi}){\nu}_{e}(\boldsymbol{B}^{2}),\end{eqnarray}$ (4) where 𝜈₀ is the reluctivity in air, 𝜈_e ( B ²) is the nonlinear magnetic reluctivity which is a function of the magnetic flux density B . 𝜓 are the design variables, and H(𝜓) determines distribution of the material density in the design domain. The area H(𝜓) = 0 represents the air region, the H(𝜓) = 1 represents the magnetic body. H(𝜓) is formulated as follows: $\begin{eqnarray}\displaystyle H({\psi})=\frac{3}{16}\left(\frac{{\psi}}{h}\right)^{5}-\frac{5}{8}\left(\frac{{\psi}}{h}\right)^{3}+\frac{15}{16}\frac{{\psi}}{h}+\frac{1}{2}\quad (-h\leq {\psi}\leq h), & & \displaystyle\end{eqnarray}$ (5) where h is the transition width between H (−h) = 0 and H (h) = 1.

2.3. Adjoint variable method

Sensitivity analysis based on the adjoint variable method (AVM) [8] allows for the calculation of the sensitivity for all design variables. In finite element analysis considering magnetic nonlinearity, the update formula of the Newton-Raphson (NR) method is defined as follows: $\begin{eqnarray}\left(K_{n-1}+\left.\frac{\partial K}{\partial \boldsymbol{A}}\right|_{n-1}\boldsymbol{A}_{n-1}\right){\delta}\boldsymbol{A}_{n-1}=-(K_{n-1}\boldsymbol{A}_{n-1}-\boldsymbol{F}).\end{eqnarray}$ (6) The converged solution A _n is obtained from the iteration of the NR method (n −1)-th iteration. The extended objective function is defined as follows: $\begin{eqnarray}\overline{W}=W+\boldsymbol{{\lambda}}^{T}(K_{n}\boldsymbol{A}_{n}-\boldsymbol{F}),\end{eqnarray}$ (7) where W is the objective function, and 𝝀 is the adjoint variables. When Eq. ((7)) is differentiated by 𝜓 and transformed, the adjoint equation is obtained as follows: $\begin{eqnarray}\left(K_{n}+\left.\frac{\partial K}{\partial \boldsymbol{A}}\right|_{n}\boldsymbol{A}_{n}\right)^{T}\boldsymbol{{\lambda}}=-\left.\frac{\partial W}{\partial \boldsymbol{A}}\right|_{n}^{T}.\end{eqnarray}$ (8) The objective function W is the average torque maximization as follows: $\begin{eqnarray}W=-t_{a}.\end{eqnarray}$ (9) ∂W∕∂ A is defined as follows: $\begin{eqnarray}\frac{\partial W}{\partial \boldsymbol{A}}=-\frac{1}{n}\mathop{\sum }_{l}^{N_{p}}\frac{\partial \boldsymbol{t}_{l}}{\partial \boldsymbol{A}},\end{eqnarray}$ (10) where t _l is the torque obtained by the nodal force method [9] at the rotor position l, and N _p is the total motor position. ∂ t _l∕∂A is defined as follows: $\begin{eqnarray}\frac{\partial \boldsymbol{t}_{l}}{\partial \boldsymbol{A}}=\mathop{\sum }_{i}\left(\boldsymbol{r}_{i}\times \frac{\partial \boldsymbol{f}_{i}}{\partial \boldsymbol{A}}\right)_{z},\end{eqnarray}$ (11) where f _i is the nodal force of node number i as follows: $\begin{eqnarray}\boldsymbol{f}_{i}=-\int _{{\varGamma}_{f}}T∼\text{grad}N_{i}dS,\end{eqnarray}$ (12)where T is the Maxwell stress tensor as follows: $\begin{eqnarray}T={\nu}_{0}\left(\boldsymbol{B}\boldsymbol{B}^{T}-\frac{1}{2}\boldsymbol{B}^{2}I\right).\end{eqnarray}$ (13)

By solving the adjoint equation of (8), the adjoint variable 𝝀 is derived, and the sensitivity is obtained as follows: $\begin{eqnarray}\frac{\partial \overline{W}}{\partial {\psi}}=\boldsymbol{\boldsymbol{\it\lambda}}^{T}\left(\left.\frac{\partial K}{\partial {\psi}}\right|_{n}\boldsymbol{A}_{n}-\frac{\partial \boldsymbol{F}}{\partial {\psi}}\right).\end{eqnarray}$ (14)

2.4. Sequential linear programming with fixed move-limit

The optimization problem is defined as follows: $\begin{eqnarray}\begin{array}{@{}l@{}}\text{min. }g_{0}(\boldsymbol{\boldsymbol{\it\psi}})\\ \text{s. t.}∼∼∼g_{i}(\boldsymbol{\boldsymbol{\it\psi}})\leq g_{i0},\quad ∼∼(i=1,2,\ldots ,m),\\ \quad ∼∼∼∼-h\leq {\psi}_{j}\leq h,\quad (j=1,∼2,\ldots ,n)\end{array}\end{eqnarray}$ (15) where $g_0$ (𝜓), g _i(𝜓) and $g_{i0}$ are the objective function to minimize, the i-th constraint condition, and the upper bound of the constraint condition. The expanded function on the k-th iteration of the topology optimization in sequential linear programming with relaxed move limit (R-SLP) [10] is obtained as follows: $\begin{eqnarray}G_{i}^{(k)}(\boldsymbol{\boldsymbol{\it\psi}})=g_{i}(\boldsymbol{\boldsymbol{\it\psi}}^{(k)})-\mathop{\sum }_{j=1}^{n}\left\{\left.\frac{\partial g_{i}}{\partial {\psi}_{j}}\right|_{\boldsymbol{\boldsymbol{\it\psi}}^{(k)}}({\psi}_{j}^{(k)}-{\psi}_{j})\right\},\end{eqnarray}$ (16) where ∂g _i∕∂𝜓_j is evaluated by the AVM. The sub-problem is generated at 𝜓^(k): $\begin{eqnarray}\begin{array}{@{}l@{}}\text{min. }G_{0}^{(k)}(\boldsymbol{\boldsymbol{\it\psi}})\\[3pt] \text{s. t.}∼∼∼G_{i}^{(k)}(\boldsymbol{\boldsymbol{\it\psi}})\leq g_{i0},\quad (i=1,2,\ldots ,m),\\[3pt] \quad ∼∼∼∼{\alpha}_{j}^{(k)}\leq {\psi}_{j}\leq {\beta}_{j}^{(k)},\quad (j=1,2,\ldots ,n)\end{array}\end{eqnarray}$ (17) where ${\alpha}_{j}^{(k)}$ and ${\beta}_{j}^{(k)}$ are the variation constraints. To solve this problem, the Lagrange function can be defined as follows: $\begin{eqnarray}L(\boldsymbol{\boldsymbol{\it\psi}},\boldsymbol{\boldsymbol{\it\lambda}})=G_{0}^{(k)}(\boldsymbol{\boldsymbol{\it\psi}})+\mathop{\sum }_{i=1}^{n_{c}}{\lambda}_{i}(G_{i}^{(k)}(\boldsymbol{\boldsymbol{\it\psi}})-f_{i0}).\end{eqnarray}$ (18)

In the R-SLP, the move limit ζ is introduced to prevent over-modification, and it is defined as the step size of the design variable update. This limit is reduced only when the objective function vibrates, and the move limit at ζ^(k) the k-th iteration steps decelerates as follows: $\begin{eqnarray}\displaystyle {\xi}^{(k)}=\left\{\begin{array}{@{}ll@{}}{\xi}^{(k-1)} & ({\psi}_{j}^{(k)}-{\psi}_{j}^{(k-1)})({\psi}_{j}^{(k-1)}-{\psi}_{j}^{(k-2)})\geq 0\\ {\tau}{\xi}^{(k-1)} & ({\psi}_{j}^{(k)}-{\psi}_{j}^{(k-1)})({\psi}_{j}^{(k-1)}-{\psi}_{j}^{(k-2)})\lt 0,\end{array}\right. & & \displaystyle\end{eqnarray}$ (19) where τ is a constant less than 1.0.

Fig. 2.

PMASynRM model.

3. Optimization model

Figure 2(a) shows the reference model of PMASynRM. The magnetic material for the rotor and stator core is T50A1300, and the material for the rotor shaft is S45C. Figure 2(b) shows the optimization model, in which the rotor area except the region of the magnet and rotor bridge are set to the design domain 𝛺_d. To facilitate manufacturing, the mirror symmetry condition of the design variable is applied, and half of the domain for 1 pitch is set to the design domain. The type of permanent magnet is bounded magnet (B _r = 0.67 T). The design variables are the values for the material density in the design domain. PMASynRM is operated using a three-phase current as i _u(t), i _v(t) and i _w(t) as follows: $\begin{eqnarray}\displaystyle & \displaystyle i_{u}(t)=I_{m}\sin ({\omega}t+{\beta}), & \displaystyle\end{eqnarray}$ (20) $\begin{eqnarray}\displaystyle & \displaystyle i_{v}(t)=I_{m}\sin ({\omega}t-2{\pi}/3+{\beta}), & \displaystyle\end{eqnarray}$ (21) $\begin{eqnarray}\displaystyle & \displaystyle i_{w}(t)=I_{m}\sin ({\omega}t-4{\pi}/3+{\beta}), & \displaystyle\end{eqnarray}$ (22) where I _m is the amplitude of the input current (3 Arms × 35 turns), and 𝛽 is the current phase angle (30°). The design goal is average torque maximization and cogging torque minimization under the area S _i(𝜓) of the magnetic body generated in the design domain, which is less than S ₀. The optimization problem is formulated as follows: $\begin{eqnarray}\begin{array}{@{}l@{}}\text{min}.\quad W=k_{1}W_{a}+k_{2}W_{c}\\ \displaystyle \text{s.t}.\quad ∼S_{i}({\psi})=\int _{{\Omega}_{d}}H({\psi})dS\leq S_{0},\end{array}\end{eqnarray}$ (23) where k ₁, k ₂ are constants. The objective function W _a, W _c are defined: $\begin{eqnarray}\displaystyle & \displaystyle W_{a}=-t_{a}, & \displaystyle\end{eqnarray}$ (24) $\begin{eqnarray}\displaystyle & \displaystyle W_{c}=\mathop{\sum }_{l}t_{cl}^{2}, & \displaystyle\end{eqnarray}$ (25) where t _a is the average torque and t _cl is the cogging torque of the rotor position l for the optimized structure. Table 1 lists the optimization parameter. The convergence criterion ϵ_opt is defined as follows: $\begin{eqnarray}{\varepsilon}_{\text{opt}}=|W^{(k)}-W^{(k-1)}|.\end{eqnarray}$ (26)

4. Optimization result

In this section, average torque maximization and cogging torque minimization of PMASynRM is investigated.

Table 1
Optimization parameter

Rotating speed [rpm] h [mm] 𝜓⁽⁰⁾ N _p ζ τ ϵ_opt

1,500 0.2 0 30 h∕10 0.99 1.0 × 10⁻⁵

Rotating speed [rpm]	h [mm]	𝜓⁽⁰⁾	N _p	ζ	τ	ϵ_opt
1,500	0.2	0	30	h∕10	0.99	1.0 × 10⁻⁵

Figure 3 shows the result of sensitivity analysis at node ID:350. In the case of $\Delta$ 𝜓 = 10⁻³, the relative error between AVM and FDM became 0.45%. From the above, only the discretization error is included in the design sensitivity. Therefore, the accuracy deterioration of AVM does not affect the convergence result.

Fig. 3.

Sensitivity analysis result of node ID:350.

Fig. 4.

Optimized structure.

Figure 4 shows the optimized structure. With the increase in k ₂, many gray-scales remained at the edge of the rotor. The optimized structure Fig. 4(d) is similar to Fig. 4(e). Because the case of k ₁ = 1, k ₂ = 0 and k ₁ = 1, k ₂ = 1 is no gray-scale, a satisfactory rotor structure was obtained. From the convergence characteristics of the objective function, vibration is suppressed by decelerating the move limit. It is necessary to reduce the initial ζ because oscillations of the objective function are strongly generated in the region less than the 50-th step optimization iteration.

Fig. 5.

Convergence characteristics of the objective function.

Figure 5 shows the results for the torque characteristics. The torque characteristic for the case of k ₁ = 1, k ₂ = 100 is similar to that of k ₁ = 1, k ₂ = 100, h∕4. Figure 6 shows the results of the cogging torque characteristics. When k ₂ is larger than k ₁, the cogging torque can be reduced significantly. Table 2 presents the average torque and torque ripple t _r of each structure. In the case of k ₁ = 1, k ₂ = 1, the average torque is higher, and the cogging torque is lowered than the reference model. Although the cogging torque is lower by increasing k ₂, the average torque does not exceed that of the reference model. The elapsed time in each optimization case became less than half a day. The source code is written in C language.

Fig. 6.

Torque characteristics.

Fig. 7.

Cogging torque characteristics.

5. Conclusions

The average torque and the cogging torque of the PMASynRM were improved by topology optimization. When the coefficient of the cogging torque minimization is larger than the coefficient of the average torque maximization, the cogging torque can be reduced significantly. However, the average torque is lower than that of the reference model. Even though the transition width of the Heaviside function was reduced, many gray-scales remained. When the coefficients were the same, a higher average torque and lower cogging torque was obtained for the rotor structure.

Table 2

Optimization result

	Reference	k ₁ = 1,	k ₁ = 1,	k ₁ = 1,	k ₁ = 1,	k ₁ = 1,
	model	k ₂ = 0	k ₂ = 1	k ₂ = 10	k ₂ = 100	k ₂ = 100, h∕4
t _a [N ⋅ m]	1.72	1.76	1.79	1.68	1.48	1.49
t _r × 10⁻¹	4.62	3.39	4.96	2.34	3.62	3.53
W _c	5.96 × 10⁻¹	1.04	1.29 × 10⁻³	9.15 × 10⁻⁵	6.21 × 10⁻⁵	6.55 × 10⁻⁵
Elapsed time [h]		6.6	7.3	7.4	5.9	5.9

CPU:Intel Core i7-6850K 3.6 GHz & 128 GB.

References

Takeuchi

Mathushita

and Hashiba

, High-efficiency reluctance motor, Toshiba Review 70(5) (2015), 20–23, (in Japanese).

Kaidou

Ishibashi

and Higashimura

, Advance and diversification of rotating electrical machines III. Progress in applied materials, T. IEE Japan 115-D(7) (1995), 835–837, (in Japanese).

Kitamura

Ishihara

Todaka

and Inoue

, Magnetic field analysis of polar anisotropic plastic magnet and application for DC brushless motor, T. IEE Japan 117-D(8) (1997), 1001–1007, (in Japanese).

Nishiura

Morimoto

Sanada

and Inoue

, Structure and characteristics of interior permanent magnet synchronous motor with bonded rare-earth magnets, IEEJ Transactions on Industry Application 134(10) (2014), 863–869, (in Japanese).

Saitou

Morimoto

Sanada

and Inoue

, A study on the rotor design of high-torque synchronous reluctance motor for vehicle propulsion, The Annual Meeting record IEE Japan (5-040) (2014), 71–72, (in Japanese).

Obata

Morimoto

Sanada

and Inoue

, Rotor structure of ferrite magnet assisted synchronous reluctance motor for torque ripple reduction, IEEJ Transactions on Industry application 134(10) (2014), 876–883, (in Japanese).

Wang

M.-Y.

Wang

and Guo

, A level set method for structure topology optimization, Comput Methods, Appl. Mech. Engrg. 192 (2003), 227–246.

Park

I.-H.

Lee

B.-T.

and Hahn

S.-Y.

, Design sensitivity analysis for nonlinear magnetostatic problems using finite element method, IEEE Trans. Magn. 28(2) (1992), 1533–1536.

Kameari

, Local force calculation in 3D FEM with edge elements, Int. J. Applied Electromagnetics in Materials 3 (1993), 231–240.

10.

Masuda

Kanda

Okamoto

Hirono

Hoshino

and Wakao

, Topology optimization of induction heating model using sequential linear programming based on move limit with adaptive relaxation, Open Phys. 15 (2017).

Average-torque-maximization and cogging-torque-minimization of permanent-magnet-assisted synchronous reluctance motor using topology optimization

Abstract

Keywords

1. Introduction

2. Method of topology optimization

2.1. Finite element method

Table 1 Optimization parameter Rotating speed [rpm] h [mm] 𝜓(0) N p ζ τ ϵ opt 1,500 0.2 0 30 h∕10 0.99 1.0 × 10−5

References

Table 1
Optimization parameter

Rotating speed [rpm] h [mm] 𝜓⁽⁰⁾ N _p ζ τ ϵ_opt

1,500 0.2 0 30 h∕10 0.99 1.0 × 10⁻⁵