Torque Ripple Minimizing of Uniform Slot Machines with Delta Rotor via Subdomain Analysis

: Since the slot opening is large in the uniform slot machine, the torque ripple generated by overlapping or misaligning with the rotor cavity is remarkably large in the case of interior permanent magnet (IPM) machine. In this work, it is observed that the magnitude of torque ripple depends strongly on the phase difference between air-gap ﬁeld harmonics: The ripple is minimized when the two dominant harmonic components cancel each other. Based on this fact, a condition is developed to minimize torque ripple by adjusting the q -ﬂux channel width and d -ﬂux barrier width. The torque ripple minimizing solution is found from a level chart made by subdomain time-stepping analysis. Finite element analysis (FEA) also gives a very similar minimizing solution. A prototype machine is manufactured, and its performances are validated through experiments. M.L. and K.N.; methodology, M.L.; software, M.L. and Y.H.; validation, M.L. and Y.H.; formal analysis, M.L.; investigation, Y.H.; writing—original draft preparation, M.L.; writing—review and editing, K.N.; visualization, M.L.; supervision, K.N.; funding acquisition, K.N. All authors


Introduction
Electric vehicle (EV) machines are often classified as machines of high current-low voltage. To accommodate high current, hairpin coils are widely used [1]. However, hairpin coils require so many welding points [2]. Since the copper welding in the dense area is complicated, it costs expensive automated welding equipment. As an alternative, open-slot structures are gradually adopted in the stator design for easy coil insertion from the inner air-gap side. In such a machine, preformed coils are used, so that the welding points are minimized [3].
Among various interior permanent magnet (IPM) rotors, delta types are used widely in the EV machine, since they yield high torque density and allow wide constant powerspeed range [4][5][6][7]. However, this type requires a lot of effort to minimize torque ripple due to various cavities and bridges. Torque ripple is mainly originated from non-sinusoidal PM field distribution and non-uniform reluctance of the iron core. As the rotor rotates, the reluctance fluctuates and the air-gap fields become asymmetric, resulting in cogging torque and and torque ripple [8].
The flux barrier and PM cavity must be designed to minimize the flux leakage while reducing the air-gap field harmonics. Another design consideration is to reduce torque ripple. In [9], multi-barrier cores were studied extensively, while considering torque, efficiency, back electromotive force (EMF) and iron loss. Small holes were added on the front edge sides of PM to reduce the torque ripple [10]. Hwang et al. [11] reduced torque ripple by placing notches on the rotor surface. Rectangular and elliptical notches were applied asymmetrically to the PM ends in the 'I'-type rotor in [12]. In an effort to increase reluctance torque while reducing the ripple, extra long barriers were put asymmetrically on the surface of 'V'-type rotor [13].
The overlap and separation of slot openings and rotor cavities are major sources of cogging torque and torque ripple [14]. With the use of hair pin coils, the slot opening can be made arbitrarily small, since the coils are inserted from a stack end. In [15], the slot opening width has been modulated to reduce torque ripple to 2∼8%. A fractional slot stator with an IPM rotor makes it easy to achieve low torque ripple because the least common multiple (LCM) of number of slot and pole are high [16,17]. However, its relatively small reluctance torque and subharmonic issue make it unsuitable for EV traction.
Skewing is a simple and traditional torque ripple reduction method applicable to all machine types. Skewing mitigates torque ripple, changes vibration mode and magnitude, but weakens torque density. As a result, noise and vibration are attenuated [18]. Multi-step skewing is mostly used in the PM machine. If the PM arrangements are not symmetrical in the axial direction, unbalanced force is created along the axis, thereby torsional vibration mode takes place in the motor housing [19]. The symmetrical PM skew in the axial direction was proposed to eliminate the unbalanced axial electromagnetic forces [20]. An alternative skewing method was attempted by carving wedge-shaped grooves asymmetrically on the tooth surface [21].
Analytic optimization is essential before the detailed design by FEA. In [22], both slotless PM model and conformal mapping were applied to predict the performance of surface permanent magnet (SPM) machine. Using the analytic results, particle swarm optimization (PSO) was introduced to find a good compromise between efficiency and torque density. On the other hand, the magnetic equivalent circuit (MEC) and the conformal mapping was used to predict the cogging torque [23,24]. Though it predicted air-gap fields accurately, a large error was observed in the cogging torque when the core was saturated.
The subdomain analysis was also utilized extensively to predict the performances of SPM machine [25,26]. Shin et al. [27] obtained optimal slot opening and pole arc ratio to minimize the cogging torque of SPM machine by using the subdomain analysis. The torque ripple change due to the flux modulation slot of SPM machine was also well predicted [28]. However, the cores are assumed to be infinitely permeable in most subdomain analysis. In order to take care of core saturation, additional domains have to be set with a finite permeability. Then the equations turn out to be the r-edge boundary problem [29,30]. In [31], both FEA and subdomain analysis were used in optimizing an IPM machine, where the PM field was obtained from FEA and the armature reaction field was obtained by a simplified subdomain method.
For low cost machine design, a uniform (open) slot stator with a delta PM rotor is a practically attractive option for easy coil insertion. In this study, we investigate a torque ripple reduction method for a uniform slot IPM machine with a delta PM arrangement. Based on the Maxwell stress tensor, we focused on the phase difference between the radial and circumferential air-gap fields, adjusting the channel width of the q-axis flux path and the cavity width of the d-axis. A design criterion is obtained by eliciting cancellation between the dominant high order components. The optimal solution is found using the subdomain analysis and FEA. This paper is organized as follows: In Section 2, a reference design and its subdomain model are illustrated and the necessary equations for subdomain analysis are presented with the torque equation. It is discussed in Section 3 that the torque ripple is directly related to the phase profiles of high order torque components in the time-stepping analysis. An index is developed to reduce the torque ripple. Level charts are made over the various combinations of design parameters in Section 4. The validity is also verified through FEA. In Section 5, a prototype machine is manufactured, and its performances are validated in experiment.

Subdomain Analysis for Delta IPM Rotor
The target power and torque of machine are 15 kW and 90 Nm. The maximum operating speed is 7000 rpm. The cross-sectional area of magnet is fixed as 3.5 × 12 mm 2 due to the limitation of PM weight and price. The main parameters of the designed machine are listed in Table 1. Differently from SPM rotor, the IPM rotor has multiple cavities to minimize the PM field leakage. The cavities for PM differentiate the rotor reluctance depending on the axis, so that L d and L q are different. However, the delta-type rotor requires relatively many bridges as shown in Figure 1a. The practical design is modified schematically as shown in Figure 1b for subdomain analysis. Note that the frontal PM to the air-gap is converted a one in the groove with a steel cap on the top. The angles of V-shaped PMs are also modified, so that they are aligned to the radial direction. Finally, the cavities in the shaft side are enlarged along the circumferential direction.
In fact, the narrow bridges around cavity are neglected in the analysis since they are always saturated by the PMs. Further, if a bridge is assumed with infinite permeability, enormous flux is supposed to leak out through it. Therefore, all the narrow bridges are removed in this subdomain model. On the other hand, we assume for simplicity that the steel cores are infinitely permeable.
In practice, small amount of PM flux leaks through the bridge. This PM flux loss is considered by lowering the remnant flux density of PM such that [32]: where B m is the remnant flux density of a real PM. w rm is the outer magnet width, w θm is the V-shape magnet width, and w boo , w boi , and w bii are the bridge width.

General Solutions and Boundary Conditions
A delta-IPM with p pole pair and Q s stator slots is considered with the following assumptions: the coil end turn effects are neglected; vector potential (A) and current density (J) have only the z-component; the core permeability is infinite; the PM is isotropic and the demagnetization curve is linear; eddy current effects are neglected; stator teeth, rotor cavities and PM have radial sides.
The subdomain analysis is performed for the equivalent model shown in Figure 1b. ω r t is the rotation angle at time t, θ si is an initial angle of i-th slot, θ m1,2j is an initial angle of spoke PM in j-th pole, θ m3j is an initial angle of radial PM in j-th pole, θ c1,2j is an initial angle of outer rotor cavities in j-th pole, and θ c3j is an initial angle of inner rotor cavity in j-th pole. Note that φ s is the slot width angle, φ so is the slot-opening width angle, φ m1,2 is the spoke PM domain width angle, φ m3 is the radial PM domain width angle, φ c1,2 is the width angles of outer rotor cavities, and φ c3 is the width angle of inner rotor cavity.
The field vectors B can be expressed as the magnetic material Equation [29] where µ 0 is the vacuum permeability, µ r is the relative permeability of material, and M is the magnetization vector for PM. In this subdomain model, both radial and circumferential magnetized magnets are used. The radial magnetization has a rectangular pulse pattern, whereas the circumferential one has a constant pattern [33]: where B Sub rm and B Sub θm are, respectively, radial and circumferential remnant flux densities obtained by (1) and (2).
By applying the Ampere's law (∇ × H = J), Gauss's law (∇ · B = 0) and vector potential definition (∇ × A = B) to (3), the Partial Differential Equation (PDE) can be obtained as Using the assumptions and magnetization sources in each domain, the governing equations for each domain are derived such as: where J i is the current density in i-th slot. Laplace or Poisson equation is solved in each region (16)∼ (24).
Subdomain 3 (Air-gap), Subdomain 4 (Outer cavity 1), Subdomain 5 (Outer cavity 2), Subdomain 6 (Radial PM), Subdomain 7 (Circumferential PM 1), Subdomain 8 (Circumferential PM 2), Subdomain 9 (Inner cavity), The boundary conditions are the following: where A sub1 and A sub2 are vector potentials in two adjacent subdomains, H θ.sub1 and H θ.sub2 are circumferential magnetic field strengths in two adjacent subdomains. The boundary conditions should be satisfied at each θ-boundaries (r = r si , r t , r a , r ob , r θm ). Each calculation of the boundary conditions is described in Appendix A.

Air-Gap Field and Torque
At each time-step, the air-gap fields are obtained from the vector potential such that The instantaneous torque is derived from the Maxwell stress tensor using (26) and (27): , L z is stack length, and r air is center radius in air-gap. Torque can be stated differently according to flux sources. Let the field densities be denoted by B r = B r.PM + B r.AR and B θ = B θ.PM + B θ.AR , where B r.PM is the radial field component created by the rotor PM, and B r.AR is the one driven by stator armature winding. The tangential components, B θ.PM and B θ.AR are also defined in the same manner. Then, the torque equation has four terms as follows [33]: Note that B r.PM B θ.PM is the cogging torque which is independent of the stator current, and B r.PM B θ.AR + B r.AR B θ.PM is the electromagnetic torque resulted from the interaction between PM and armature reaction fields. The last one, B r.AR B θ.AR is called reluctance torque since it is originated from the variation in the armature reaction field.

Air-Gap Field Comparison between Subdomain and FEA
To verify the accuracy of the subdomain results, their air-gap fields are compared with the FEA results under different current loading and angle conditions. Two results are compared first when only the PM is magnetized, and second when only the stator current is conducting. Note that Q s /GCD(Q s , p) = 12, which is the number of slots per pole pair. Since the air-gap field does not have even harmonics, the most dominant harmonic numbers are found at Q s /GCD(Q s , p) − 1 = 11 and Q s /GCD(Q s , p) + 1 = 13. Figure 2 shows the air-gap fields when only the PM is magnetized at ω r t = 0 where the cogging torque is zero. It also shows magnitudes and phases of the field harmonics. Both the magnitude and phase of dominant harmonics show good agreements between the subdomain and FEA results. When the cogging torque is maximized at ω r t = 5.75 • , similar results are obtained as shown in Figure 3. The second part is shown in Figure 4, where the PMs are unmagnetized. Instead, current i q = 185 A is injected to the stator coil with zero current angle. Note that the current angle (β) is defined as the angle of the current vector with respect to the q-axis. When β = 0, the current vector consists of only the q-axis component, and the d-axis current is equal to zero. Similarly to the previous results, the most dominant harmonic numbers are found at 11 and 13 orders. The subdomain and FEA results generally agree well, but there are some differences in the magnitude of the harmonics. However, the phases of the dominant harmonic components appear to be approximately the same.  Similar results are obtained in Figure 5 for ω r t = 5.75 • with the largest torque ripple. Note again that the harmonic phases agree well.

Torque Ripple and Harmonic Phase Condition
According to (28), the phase difference α rl − α θl of the field harmonics determines the lth order torque component along with the magnitude, T l . In the machine under consideration, the 11th and 13th harmonic components are dominant. Thus, we consider just the two harmonic components and let ∆T(t) = T 11 (t) cos(σ 11 (t)) + T 13 (t) cos(σ 13 (t)), (29) where σ 11 (t) = α r11 (t) − α θ11 (t) and σ 13 (t) = α r13 (t) − α θ13 (t). Note that ∆T represents approximately the reluctance torque ripple since the PM is unmagnetized. When the rotor rotates, ∆T varies as the flux paths change microscopically around the rotor cavities and bridges. Specifically, the magnitudes and phases of the space harmonics vary as the time-stepping progresses. Figure 6 shows the time-stepping analysis of two different designs: Design A and Design B have the same arc φ d = 20 • for the d-axis flux barrier, but different arcs φ q = 5 • and 8 • for the q-axis channel. Figure 6b,c shows the magnitudes T 11 (t) and T 13 (t), and the phase differences σ 11 (t) and σ 11 (t). Comparing T 11 (t) and T 13 (t) of Design A and B, they, though different, are bounded in both cases. On the other hand, the phases difference, σ 11 (t) and σ 13 (t) of the space harmonics show completely different profiles in time. With Design A, the 11th and 13th have almost the same phase in time. Thereby, T 11 (t) cos(σ 11 (t)) and T 13 (t) cos(σ 13 (t)) are hardly canceled out, resulting in a large torque ripple ∆T. In contrast, σ 11 (t) and σ 13 (t) of Design B are almost 180 • out of phase, so that T 11 (t) cos(σ 11 (t)) and T 13 (t) cos(σ 13 (t)) tend to have opposite signs. Correspondingly, ∆T of Machine B is small. Based on this observation, we are motivated to design the machine such that the signs of cos(σ 11 (t)) and cos(σ 13 (t)) are opposite. It is proposed in this paper as a torque ripple minimizing method that Thereby, the ripple minimizing design can be preceded by reducing the index |π − (σ 11 (ω r t) + σ 13 (ω r t))|∆ω r t.
The closer χ to zero, the smaller the ripple under the assumption that T 11 (t) = T 13 (t).

Rotor Cavity Design Minimizing Torque Ripple
The slot area is regarded as the magnitude of ampere-turn. However, a wide slot results in a narrow tooth width, increasing the teeth reluctance. Thereby, the stator magnetomotive force (MMF) must be balanced with the teeth reluctance. In most uniform slot machines, the tooth width is about the same as the slot opening, i.e., the slot angle is generally φ s = φ so ≈ 0.5 × 2π/Q s . We let the angles of cavity bridges are fixed such that φ c1,2 = 2 • and φ m1,2 = 3 • . As the design variables, we select the arcs of d-axis flux barrier and q-axis channel. Specifically, φ d and φ q shown in Figure 6a are considered the design variables. They are normalized by the slot pitch such that γ q ≡ φ q 2π/Q s and γ d ≡ φ d 2π/Q s . The proposed index χ is the error phase difference sum to the ideal phase difference when the rotor rotates half-slot. That is, it is an error accumulation value and its value is changed by the number of samples. In this work, χ on each design is calculated using 3 time-steps when rotation angle is 1, 2, 3 • . The maximum harmonics for each subdomain are set as: K = 50, G = 50, L = 100, M = 50, N = 50, U = 50, V = 50, W = 50, and H = 50. So, the computation time for subdomain analysis is 5 s under following computer system: Intel(R) Core(TM) i5-6600 CPU @ 3.30 GHz, installed memory 16.0 GB, and MATLAB program. Figure 7a shows a level chart of χ over (γ d , γ q ) obtained by subdomain analysis. The minimum χ is found from the numerous combinations of (γ d , γ q ). Note that the contour lines represent equilevel curves, i.e., sets of the same χ values. The optimal design is obtained when γ d = 2.67 and γ q = 1.06. Figure 7b shows a level chart of the torque ripple. Note that the minimum value is found at an almost same point predicted by the subdomain analysis. It shows clearly a strong correlation between the torque ripple and the phase criterion χ. FEA for torque ripple analysis is set as following: 19,158 meshes, the rotation angle per step is 0.25 • , the rotation angle is from 0 to 15 • , and 61 time-steps. So, the computation time for FEA is 270 s using 2D magnetic field transient analysis of JMAG Designer.
The same optimization study was done via FEA with the practical model shown in Figure 1a under the same condition for PM and current. The equi-level contours are shown in Figure 8a,b. The general trend look similar as the one shown in Figure 7. Though a little differences are shown where γ d and γ q are small, the optimum point is found at the same location as the subdomain analysis. Furthermore, the time to find the optimal point is greatly reduced.  The above results are made with PM unmagnetized. Now, Figure 9a shows an optimization result by FEA when the PM is magnetized and current angle is β ≡ tan −1 (−i d /i q ) = 44 • . It also yields the same optimal point. Figure 9b shows the torque wave plots. When it is compared with an initial design, torque ripple is shown to be reduced significantly.  Figure 10a,b shows a delta-IPM rotor and uniform slot stator with a stranded winding. Experiments were done using a dynamometer as shown in Figure 10c. The test inverter was set up as follows: pulse width modulation switching frequency was 5 kHz, current control period was 200 µs, TMS320F28377 processor was utilized, and switch component was Infineon FF450R12KT4. Due to the low resolution of the torque sensor, the torque ripple were not measured. Line to line back electromotive force (EMF) is shown in Figure 11a,b.    Figure 13 shows the comparisons between FEA results and experiment results. The back EMF magnitudes agree well in Figure 13a. The torque-speed performances are compared in Figure 13b. The same results are obtained below 5000 rpm, but less torque is generated in practice over 5000 rpm due to stray and iron losses. The FEA-based efficiency is calculated as follows: Efficiency = Tω r Tω r + P cu + P f e + P stray + P mech ,

Prototype Motor and Experimental Results
where P cu = 3 2 R ph i 2 d + i 2 q is the copper loss and P f e is the iron loss obtained by FEA. The stray loss is estimated as 10% of iron and copper losses, P stary = 0.1 × (P cu + P f e ). The mechanical loss is assumed to be linear in speed and approximated as P mech = cω r with c = 0.05 W/rpm [34]. Meanwhile, the experimental efficiency is calculated as Tω r P in , where P in is the input power obtained by the power analyzer. Table 2 shows the comparison of the efficiency between the experimental and calculated results. Below the base speed, similar results are obtained. However, the difference is more pronounced in the field-weakening region, and the difference becomes larger as the speed increases. It seems to be caused by unmodeled AC loss of the coil and inaccurate stray losses.

Conclusions
An IPM machine is developed with uniform stator slots for small EVs. This type machine has the advantage of reducing the winding cost, since a rectangular coil can be inserted from the inner wall side without an expensive coil inserter or copper welding machine. However, in general, there is a problem of large torque ripple. According to the Maxwell stress tensor, the torque is obtained as an integral sum of B r × B θ in the air-gap. Therefore, the phase differences of the radial and circumferential harmonic fields remain as the ripple components. When the number of slots per pole pair is 12, the phase differences of 11th and 13th are dominant. If the two harmonic components cancel each other, the torque ripple is minimized. Specifically, they should have 180 • difference in order to minimize the torque ripple. The desired cancellation condition is formulated as an index, and the minimizing solution is found by using subdomain analysis. The similar optimal solution was also verified by FEA. However, the subdomain analysis took 5 s to compute for each design, whereas FEA took 270 s. A prototype motor was manufactured and its performance was verified through dynamometer test.

Conflicts of Interest:
The authors declare no conflict of interest.

Abbreviations
The following abbreviations are used in this manuscript: where E k = ( r si r so ) F k , F k = kπ φ s , E g = ( r si r t ) F g , and F g = gπ φ so , 2 )dθ, In the same way, continuity of tangential magnetic field strength (H θik = H θig ) gives At r = r t and θ ∈ [θ si − φ so 2 , θ ci + φ so 2 ], continuity of vector potential(A ig = A l ) yields where E l = ( r a r t ) l p , η il0 = 1 φ so θ si + φso 2 ))dθ, 2 ))dθ. In the same way, continuity of tangential magnetic field strength (∑ i H θig = H θl ) gives At r = r a and θ ∈ [0, 2π p ], continuity of vector potential (A l = A jm ) yields L ∑ l (a l E l + b l )η jc1 + L ∑ l (c l E l + d l )ζ jc1 = a jm0 + b jm0 ln r a , (A9) L ∑ l (a l E l + b l )η jlm + L ∑ l (c l E l + d l )ζ jlm = a jm + b jm E m , where E m = ( r a r ob ) F m , F m = mπ φ c1 . Continuity of vector potential (A l = A jn ) yields L ∑ l (a l E l + b l )η jc1 + L ∑ l (c l E l + d l )ζ jc1 = a jm0 + b jm0 ln r a , (A11) where E n = ( r a r ob ) F n , F n = nπ φ c2 . Continuity of vector potential (A l = A jv ) yields L ∑ l (a l E l + b l )η jm1 + L ∑ l (c l E l + d l )ζ jm1 = a jv0 + b jv0 ln r a − µ 0 r a M jθ1 , (A13) L ∑ l (a l E l + b l )η jlv + L ∑ l (c l E l + d l )ζ jlv = a jv + b jv E v , where E v = ( r a r ob ) F v , F v = vπ φ m1 . Continuity of vector potential (A l = A jw ) yields L ∑ l (a l E l + b l )η jm2 + L ∑ l (c l E l + d l )ζ jm2 = a jw0 + b jw0 ln r a − µ 0 r a M jθ2 , (A15) where E w = ( r a r ob ) F w , F w = wπ φ m2 , η jc1 = 1