On-Load Analytical Modeling of Slotted Interior Magnet Synchronous Machines Using Magnetic Islands Method

In general, the Magnetic Circuit (MC) is an effective and valuable tool for calculating the PM flux and predicting machine behavior during the no-load operation. Since the rotor PMs are the sole source of the flux, the flux path and the respective MC model is independent of the rotor position. However, this is not the case for the loaded operating condition. As the rotor spins, the path of the armature flux necessarily changes to pass through the rotor with minimum reluctance. Therefore, the MC model would vary at every rotor position and is no longer a feasible solution particularly, for cases with complex rotor structure. To resolve this, an analytical solution procedure based on the concept of the Magnetic Islands (MI) is presented to predict the on-load response of the motor such as armature flux and the winding inductance. The air-gap is defined as a function along the perimeter of the stator to include the stator slotting effect. In addition, the presented method accounts for the magnetic saturation. To assess the accuracy of the model, analytical results are compared against the finite element results.


I. INTRODUCTION
The performance characteristic of Permanent Magnet (PM) machines is closely related to the quality of the air-gap flux distribution. In particular, the spatial harmonics of magnetic flux can generate a pulsating torque, which could potentially vibrate the motor and cause damage on the load and the shaft [1]- [3]. A significant portion of the low order spatial harmonics are related to the stator slots and the rotor flux barriers. Therefore, a deep understanding of the relation between the rotor/stator layout and the flux distribution could serve as the foundation of a successful design. There are primarily two approaches (i.e. Finite Element (FE) analysis and analytical model) to achieve this objective. Despite the accuracy of the FE method, it could take excessive computational time, and is not justified at the early design stage. By contrast, analytical solution could converge in a timely manner with a good precision. In fact, a coherent design framework involves The associate editor coordinating the review of this manuscript and approving it for publication was Pierluigi Siano . analytical modeling at the early design stage and FE modeling at the final stage for the design refinement and validation.
One of the most accurate analytical techniques available in the literature is the subdomain method, where the air gap is divided into different subdomains [4], [5]. In this technique, the magnetic flux density in each subdomain is calculated by solving Poisson's equation and applying boundary conditions. However, the calculation of magnetic flux in the stator slot subdomain is redundant and increases the computational volume.
An analytical model based on the surface magnetization current is presented in [6] to calculate d-q axes inductance values and the electromagnetic torque in synchronous PMinset machines. However, the solution procedure requires iterative steps to calculate equivalent surface currents. In [7]- [10], the motor performance is modeled by developing a MC model and a matrix based on the permeance of the flux paths. In [8], the stator slotting effect is also included by using conformal mapping method. However, the on-load operation cannot be modeled via a single MC model. As the motor rotates, the armature flux necessarily changes to find a new path with minimum reluctance. Therefore, the MC method could be complex at the loaded operating condition as it requires different MC models at every rotor position. In addition, the use of conformal mapping to account for stator slots would further increase the complexity of the process.
In [11] the distribution of the magnetic field is obtained in a V-shaped IPM synchronous machine by transferring the rotor into the polar coordinate system and utilizing the magnetic scalar potential. In [12] the motor torque ripple, and in [13] the iron loss are calculated by computing the magnetic voltage drop across flux barriers, and the magnetic potential across the rotor. Further, in [12], the rotor flux barriers are designed with favorable shape and position to reduce the torque ripple.
In this paper, to overcome the shortcomings associated with the MC, the loaded analysis of the motor is carried out via the concept of Magnetic Islands (MI). With this method, the relative position between the rotor and the stator is automatically accounted during the modeling process. For this purpose, the rotor is classified into a number of segments (islands) based on the magnetic potential of each segment. The air-gap armature flux and the winding inductance are then obtained by calculating the respective magnetic potential of each segment at every positions of the rotor. In addition, the air-gap is defined as a function along the perimeter of the stationary stator to account for the stator slotting effect. Finally, FE analysis is performed to verify the validity of the analytical results.

II. DESCRIPTION OF THE MOTOR
The multi-layer Interior Permanent Magnet (IPM) structure can partially bring the flux distribution into uniformity and reduce the pulsating torque [8]. Another key factor that affects the quality of the flux distribution, is the structure of the PMs, flux barriers and the overall rotor topology [14]. In this paper, a double-layer, U-shape PM rotor structure is adopted.
The schematic of the selected motor with the related parameters is illustrated in Fig. 1. The stator has 36 slots with three coils/pole/phase. The motor data/specifications are presented in Table 1.
In Fig. 1, 2απ/N p is the angular length of the inner (layer I) magnet arc on the rotor circumference. Similarly, 2α 1 π/N p is the angular length of the outer (layer II) magnet arc on the rotor circumference.
Iron bridges are utilized in the rotor to enhance the mechanical stress. However, to avoid PM flux getting wasted by forming short loops inside the rotor core, the thickness of these bridges is very small (i.e. 0.5 mm). Considering the nonlinear magnetic characteristic of the iron core (Fig 2), these narrow bridges tend to get saturated fairly quickly. As a result, the flux density inside and in the vicinity of these bridges is assumed to be constant (≈2 T) during the motor operation.

III. MAGNETIC ISLANDS METHOD
Since the magnetic flux tends to follow a path with the least magnetic resistance, the path of the armature flux path  continuously changes during the operation of the motor. As a result, a single MC model cannot adequately represent the motor response, but multiple MC models are required. This means that the MC model should be updated per rotor position to accurately reflect the new flux path. Fig. 3 illustrates the magnetic flux lines of phase a of the armature winding at four different rotor positions of θ To avoid the necessity of using multiple MC models, an analytical modeling method based on the concept of MI is presented based on the following procedure. The magnets are removed and the machine is split into a number of islands, each with a different magnetic potential. The whole stator is considered as a single segment with magnetic potential U s , and the rotor core is classified into 10 segments with magnetic potentials U r 1 to U r 5 in the upper half and U r1 to U r5 in the lower half ( Fig. 4). Rotor segments are separated from adjacent rotor segment(s) via a flux barrier [15], [16]. The first step is to calculate the magnetic potential of these segments (i.e., stator and rotor segments). In Fig. 4, P 1 to P 32 pinpoint the angular location of barriers along the air-gap. These points define the borders per rotor island. Due to the symmetrical structure of the motor, the solution procedure is presented only for the upper half of the machine but it can be easily extended to cover the entire machine.
It is worth mentioning that in the overall analysis of electric motors, the local magnetic saturation inside the iron bridges is modeled only once, (normally during the no-load analysis) as a constant current source [8]. Therefore, at the loaded condition the saturation points are modeled only as reluctance elements [15].
The stator winding is excited by an arbitrary DC current to calculate the winding inductance (since the inductance value is independent from the current amplitude). Therefore, the stator magnetic potential is time invariant and is constant per coil per phase (independent from the rotor position) as expressed in (1). The solution procedure remains identical for AC excitation current, but the time factor will be added. The MMF of phase a of the stator winding is [15]: where n a m (ϑ) is the turn function, N a m (ϑ) is the winding function, and I is the current in phase a. The variable m is the number of coils per phase.
Relation (1) can be expressed in the following form (2) based on the position of the coils in phase a The resultant MMF of the three phase stator winding (U s ) is simply calculated by the sum of the MMF of each phase. Next, the magnetic potential of the nth rotor island is obtained by integrating the armature MMF along the angular length shared between stator and the same (nth) rotor island.  Finally, the profile of the armature magnetic flux in the airgap (along 360 mechanical angle) is obtained by connecting the magnetic flux profile across every rotor islands.
The magnetic potential of each of the rotor segments U r1 to U r5 can be expressed as the following [15]: where R b1 to R b4 are the reluctance of the flux barriers in between the rotor islands as given in (5), and ϕ b1 to ϕ b4 are the flux passing through these barriers.
where R o1 represents the reluctance of the radial iron bridge located at the center of the inner layer I, R o2 is the reluctance of tangential iron bridges located at both ends of the inner layer I, R F1 is the barrier reluctance of inner layer I and R F2 is the barrier reluctance of outer layer II [15]. Since the thickness and the length of the saturated area in both inner and outer flux barriers are identical, the respective magnetic resistances are equal (R o1 = R o11 and R o2 = R o12 ). Due to the symmetry of the rotor, relation (4-a) is the negative of (4-d), and the relation (4-b) is the negative of the (4-c). The symmetrical structure of the motor reduces the number of unknown variables but it is still higher than the number of equations. In fact, parameters ϕ b1 , ϕ b2 ,R b1 , R b2 and U r3 are required to solve (4). These five parameters can be calculated with the help of five additional relations of (6) through (10). Because the length of the motor (L stk ) is fixed, the surface integral in (6), (7) and (8) can then be converted into a line integral along the stator circumference: where the parameter R is defined as the sum of the rotor radius and half of the length of the air-gap (g), L stk denotes the motor length, and B gar(1:3) represent the air-gap armature flux density associated with each of the rotor islands one, two and three.
where U s (1:5) are the stator MMF along the angular interval of rotor islands one to five and obtained from the relation (3). U r (1:5) are the scalar magnetic potential of the rotor islands one to five, and obtained by solving the relation (4). By inserting relation (9) into (8), the scalar magnetic potential of island three, U r3 can be expressed as relation (10).
As the rotor spins, the relative position between rotor islands and the stationary stator, and therefore the stator MMF across the rotor islands vary. This means the amount of the armature flux passing through rotor islands and their respective magnetic potential will vary at every rotor position. Therefore, the magnetic potential per island is recalculated at every rotor position.
The overall air-gap flux density along 360 angle are given by: The air-gap function over the circumference of the stationary stator [15].
where n represents the number of rotor islands. U s(n) denotes the stator scalar magnetic potential right across the nth rotor segment. Likewise, B gar(n) denotes the air-gap flux density across the nth rotor segment. The term in (11) is not a summation sign but only an indicator that connects the magnetic flux profile in the air-gap, across each of the five rotor segments to complete the 360-degrees profile of the armature flux in the air-gap.
Once the air-gap flux density is determined at every rotor position, the inductance values can be calculated from (12). The overall per flux linkage per phase equals to the summation of the magnetic flux passing through each coil (The winding is arranged by 3 coils/pole/phase).
where N a (1:3) represents the turn number for each of the three coils in phase a. Similarly, N b(1:3) represents the turn number in phases b and c. In this paper, the number of turns per coil per phase is identical and equal to N in Table 1.
The length of the air-gap is 0.5 mm and the effective airgap along the stator slots is 2 mm. Therefore, the area where the flux passes through the air-gap is proportional to the ratio of the stator teeth and the stator inner surface as indicated in the air-gap function in (13) [15].
where ϑ sq and γ s are the location and the slot opening angle. Fig. 6 illustrates the air-gap function over the stator reference frame, based on (13).
Once the motor inductance matrix is found, the reluctance torque can be calculated via (14-a). Also, the magnet torque can be calculated via (14-b) using magnetic circuit analysis [8], [15]. Finally, the electromagnetic torque is obtained by adding the reluctance torque and the magnet torque as expressed in (14-c).
4 L stk N s I s B gpm (14-b) T el = T rel + T pm (14-c) where B gpm is the air-gap flux density of the rotor magnets. Fig. 7 illustrates radial distribution of the armature current flux density (per-phase) in the air-gap over the circumference of the stator based on (11). The flux density is calculated and plotted at rotor positions of θ r = 20 and θ r = 60 degrees in regards to the stationary stator. The resultant three phase armature flux is plotted in Fig. 8. To verify the validity of the proposed technique, analytical results are plotted and compared against the FEA results. As observed, the stator slots create higher order harmonic contents in the flux profile. The armature flux profile in Fig. 7 and Fig. 8 is a function of ϑ(mechanical angle along the circumference of the stationary stator), and θ r (rotor angular position) and are plotted at a particular time instant t = t 1 .  The self-inductance of phase a and the mutual inductances are calculated based on (12) and illustrated versus rotor angular position (θ r ) in Fig 9. The analysis are all carried for both slottless and slotted motor to highlight the slotting effect.
Finally, the profiles of the reluctance torque, magnet torque and the total electromagnetic torque are calculated based on (14) and plotted in Fig. 10. These results are obtained when the stator winding is excited by a 3-phase, balanced sinusoidal current with the amplitude of 2 A.  The presented analytical solution is general and may be applied to other IPM machine topologies. For cases with different stator winding arrangement, the MMF waveform is updated through relation (1).
On the other hand, the number of rotor islands and therefore, the number of equations in relation (4) are a function of layers of flux barriers. For rotor structures with different layers of flux barrier, the number of rotor islands, and the number of equations in (4) would change as listed in Table 2, but the solution procedure remains intact.
For rotor structures with different barrier arc angle, the number of rotor islands remains constant, but the angular position of barrier ends (P points in Fig. 4) would change. The P points define the integration limits in relations (6) to (8).

IV. CONCLUSION
An analytical modeling technique was presented to model the steady-state performance of the slotted IPM synchronous motor at loaded operating condition. For this purpose, the distribution of the armature magnetic flux density was modeled through the concept of the MI. Unlike the equivalent magnetic circuit model, the presented technique can operate without requiring any information on the path of the magnetic flux. Also, the relative position between the rotor and the stator is automatically accounted for, throughout the process. Moreover, the magnetic saturation is accounted through the process by modeling saturated region at the iron bridge with reluctance elements. The presented methodology is general and maybe extended and applied to other classes of IPM machines with a different number of magnetic layers and different magnet shape. In addition, the effect of the stator slotting were included by applying air-gap function over the perimeter of the stationary stator. The air-gap function accounts for the stator teeth with a good precision and is not as complex as the conformal mapping method. To verify the accuracy of the presented technique in predicting the motors performance, analytical results were all validated through the FE analysis.