A Versatile Analytical Method of Investigating MMF Harmonics of Armature Windings

Analytical methods are of great interest for the design and analysis of AC machines. This paper proposes a new analytical method to study the harmonic rotating components of the magneto motive force (MMF) produced by various types of armature windings and input currents. The windings can be symmetrical-phase or asymmetrical-phase, single-layer or multi-layer, fractional-slot or integer-slot, with arbitrary phase numbers. The analytical method is expressed as concise steps and closed-form equations, and its validity and reliability are verified by comparison to finite element analysis (FEA) simulations and experiments. As opposed to previous works, the proposed method establishes a more comprehensive quantitative relationship among the harmonics of the resultant MMF, the arrangement of armature windings, and the input currents. Practical implications of this method include enabling more accurate prediction of machine performance, as well as optimizing machine design in both quantitative and qualitative respects.


I. INTRODUCTION
Though finite element analysis (FEA) has many advantages with regard to its mathematically and practically proven accuracy in analyzing the electromagnetic characteristics of AC machines, it is usually unintuitive to determine the degree of contribution of certain current or winding distribution harmonics to the resultant field, torque or vibration. Furthermore, FEA requires a time-consuming calculation process and has a high computational cost while lacking closed-form solutions. Thus, these limitations make it difficult for FEA to reveal the universality and commonality among various AC machines.
In this paper, a versatile analytical method for analyzing the electromagnetic characteristics of AC machines is studied and proposed to address the issue. The armature winding magneto motive force (MMF) is chosen as the subject of investigation. On one hand, MMF can reflect the field distribution [1], back EMF [2], [3], unbalanced magnetic force (UMF) [4], [5], vibration and noise [6], [7]. On the other hand, as a spatio-temporal entity, MMF is directly produced by The associate editor coordinating the review of this manuscript and approving it for publication was Padmanabh Thakur . temporal alternating currents carried in spatially distributed armature windings, making it relatively straightforward to derive mathematically.
The analytical treatment of the armature winding MMF dates back to the 1900s [8], where the mathematical derivation focuses on full-pitch concentrated and full-pitch or short-pitch distributed symmetrical armature windings. A preliminary investigation of the fundamental and harmonic components of the MMF induced by these windings is also included. However, the mathematical model is limited to the previously mentioned winding types only. Independently, [9] uses a vector-based mathematical approach to derive the stator winding rotating field, named as the gliding field, in a symmetrical polyphase induction motor. This approach also applies to MMF fields, but the harmonic characteristics are not thoroughly investigated. [10] further extends the analytical treatment of MMF in the previously stated literature to a larger variety of cases, including asymmetrical armature windings, and discusses the presence and absence of harmonics in the resultant MMF. Nevertheless, the effect of non-sinusoidal currents is not analyzed in depth. Similar mathematical treatments of MMF of fractional-slot armature windings can be seen in [11] on three-phase windings, and in [12] on general symmetrical windings. The methods of analysis from former studies are also recapitulated in [2], [13], [14], [15], and [16]. As some mathematical models for analyzing armature windings are being developed and improved, methods using these new tools to analyze armature winding MMF come into view. Space vectors and the composition of complex quantities are adopted in [17] to analyze the spatial harmonics of MMF produced by symmetrical multiphase windings. Similar conclusions to those from [10] are drawn in this work, but additional mathematical transformations are also introduced. A concept called the winding function is elucidated and used in the mathematical derivation of armature winding MMF in [18], [19], [20], and [21], and is widely studied in works such as [3] and [22] to analyze the MMF of fractional-slot concentrated windings. However, the aforementioned methods have limitations and typically apply to specific types of windings, particularly symmetrical ones. Furthermore, they involve laborious handling of pole-slot combinations and winding configurations, and only offer a limited amount of quantitative analysis regarding the correlation between winding distribution, nonsinusoidal current, and the resultant MMF.
In [23], the authors propose a generalized analytical method for the MMF harmonic analysis of armature windings, which is versatile enough to be applied to both symmetrical and asymmetrical, integer-slot and fractionalslot windings. The analysis reveals the interactions between spatial and temporal harmonics, and the resultant MMF harmonics. To maintain universality, the conclusions do not take into account factors such as the stator slots or poles. However, the derivation process omits the calculation of a crucial factor, and there is no experimental validation to verify its conclusions. This paper complements the authors' work in [23] by determining the missing crucial factor, referred to as the placement multiplier in this paper, and validating the method with both FEA simulations and experiments.
The key contributions and novelty of our proposed method in this paper, compared with the previous literature, are : 1) The proposed method can be applied to armature windings with arbitrary phase numbers, whether they are symmetrical or asymmetrical. It also offers a way to perform an equivalence transformation on specific phases. 2) The proposed method is inherently applicable to integer-slot and fractional-slot armature windings, as well as non-sinusoidal currents, since it takes into account all possible combinations of odd and even orders of harmonics in both spatial and temporal domains.
3) The proposed method aggregates the contributions resulting from various pole-slot combinations into one factor, and presents a clear procedure to calculate the factor, ruling out the distinct effect of pole-slot combinations. 4) The proposed method uncovers a comprehensive yet concise quantitative correlation among the harmonics of the resultant MMF, the arrangement of armature windings, and the currents, encompassing their orders, amplitudes, and phase angles. 5) The proposed method has been validated by both FEA simulations and experimental results. 6) The proposed method can be easily integrated with and well supplements existing studies, such as [24] and [25], which place greater emphasis on the formation of static spatial MMF harmonics in different kinds of windings.
It also supplements studies like [26], which focus more on the generation of current harmonics due to PWM.

II. ANALYTICAL DERIVATION
To make the analytical method versatile for various winding structures and current waveforms, the following conditions are used in the derivation: • The superposition of MMFs is used, thereby implying linearity in the calculation to compose the resultant MMF.
• End windings and fringing effects are ignored. This implies isotropy along the axial direction, thus allowing the AC machines to be modeled in 2D space.
• To cover all possible situations thoroughly and theoretically, the currents carried in the armature windings can have arbitrary orders of harmonics up to infinity. This means that abrupt changes are allowed in currents, as a step function can be decomposed into an infinite number of harmonics.
• Similarly, the orders of harmonics in the distribution of armature windings can also be arbitrary and infinite in the derivation process.
This paper mentions two specific types of harmonics: spatial harmonics and temporal harmonics. Spatial harmonics refer to the harmonics in space that are introduced by the distribution of armature windings, while temporal harmonics refer to the harmonics in time that are introduced by the alternating currents carried in the windings. The resultant MMF is a collective product of the spatial and temporal harmonics, and manifests itself as a spatio-temporal traveling wave.
The winding function, as proposed and studied in [18], [19], [20], and [21], is defined as the MMF per unit current in a single phase. It is selected in this paper to reflect the spatial harmonics. Correspondingly, the phase current, which represents the alternating current in a single phase, is used in this paper to reflect the temporal harmonics.
This analytical method focuses on the MMF generated by the armature windings. Therefore, no rotor-associated parameters like pole pairs or pole-slot ratio are considered. Hence, there is no such distinction as electrical frequency and mechanical frequency. All spatial harmonic orders are with reference to the fundamental harmonic, which is a complete cycle in space. Meanwhile, temporal harmonic orders are with reference to the designated fundamental phase current.
In symmetrical-phase windings, only windings with an odd number of phases are considered. This is because any even-number-phase windings can be equivalent to windings with an odd number of phases and a certain number of winding sets. The method of equivalence will be presented in detail in this section.
The analytical derivation process is organized progressively. It starts with single-phase winding, then moves on to symmetrical-phase windings, and finally ends with asymmetrical-phase windings.

A. SINGLE-PHASE WINDING
Assuming there is a polar axis in space referred to as the original position, then at any space position angle β, the winding function N of a single-phase armature winding with respect to β and phase axis position φ s (the spatial phase shift of the winding) can be expressed as where ν s denotes a certain order of spatial harmonics,N ν s and ϕ sν s represent the amplitude and the harmonic-order-related initial phase of the ν s th order harmonic, respectively. Similarly, the phase current I in a single-phase armature winding is a function of time t and current phase shift φ t (the temporal phase shift of the current). It can be written as where ν t represents a certain order of temporal harmonics, I ν t and ϕ tν t are the amplitude and the harmonic-order-related initial phase of the ν t th order harmonic, respectively, and ω is the fundamental harmonic angular frequency of the alternating current. The MMF, denoted as in this paper, of a single-phase winding is the multiplication of the winding function and phase current of the same phase, thus a combination of spatial and temporal harmonics. By applying trigonometric transformations, the stated MMF can be broken down into forward and backward rotating harmonics ± (+: forward,

B. SYMMETRICAL-PHASE WINDINGS
To calculate the MMF composed by a set of symmetricalphase windings from that of a single-phase winding presented above, the spatial arrangement of the phases is very crucial. It is usually determined by the desired number of poles and slot-pole combination principles during the machine design stage. Even so, as mentioned above, the analytical approach proposed in this paper does not require the introduction of slot or pole-related factors to maintain universality in conclusions.
To address this issue, a placement multiplier x is introduced. Although the spatial phase shift φ s for each phase winding can vary, the displacement of every two adjacent phases is always equal. The temporal phase shift φ t is specified as the reference value. Therefore, φ s can be expressed as the product of x and φ t , as shown in (4). It is worth noting that x can be either an integer or a fraction.
As both φ t and φ s are angles that have a period of 2π, the desired value of x is non-unique for any winding configuration. For clarity, in this paper, x is constrained to its minimal positive alternative.
In the set of symmetrical-phase windings with k phases, the temporal phase shift of the (n + 1)th phase current is 2πn k , where n = 0, 1, 2, · · · , k−1. By introducing Euler's formula, forward and backward rotating MMF harmonics composed by the k symmetrical phases can be derived as In (5), the function Re returns the real part of its argument. e and i are mathematical constants representing Euler's number and the imaginary unit, respectively. By observing the items of the summation, there are two cases to be discussed concerning the relationship between xν s ∓ ν t and k: (i) If xν s ∓ ν t is an integer multiple of k, the corresponding harmonic in (5) can be simplified to The produced harmonic has an amplitude of k 2N ν sÎ ν t , a rotation angular speed of ± ν t ν s ω, and an initial phase of ϕ sν s ∓ ϕ tν t . (ii) If xν s ∓ ν t is not an integer multiple of k, the produced harmonic is 0. This means no MMF harmonic will be composed in this case. It should be noted that the above-mentioned derivation is satisfied on the premise that k is an odd number.

C. ASYMMETRICAL-PHASE WINDINGS
A set of asymmetrical-phase windings is constructed with multiple sets of symmetrical phase windings. The currents in each set get shifted by a temporal angle of π gk consecutively, where k is the number of phases in each set and g is the number of total sets. The spatial phase shifts of the windings in each set can be determined accordingly using (4). Obviously, symmetrical-phase windings can be seen as a special case of asymmetrical-phase ones where g equals 1.
If k is an even number, a method of equivalence is introduced to transform it into an odd number: 1) Let s equal the number of total regions that the space is divided into by all the phase axes: gk; 2) Express s as a product of an even number and an odd number; 3) Let the even number be 2g ′ and let the odd number be k ′ . Replace the original values of g and k with g ′ and k ′ respectively. According to the definitions in the previous subsection, the forward and backward rotating MMF harmonics in the asymmetrical g-set k-phase windings can thus be expressed as Based on the conclusions regarding the symmetrical phase windings, it is assumed that xν s ∓ ν t is an integer multiple of k. Therefore, (7) can be further simplified to By observing the items of the summation, there are three cases to be discussed concerning the relationship of xν s ∓ ν t with g and k: (i) If xν s ∓ ν t is an even multiple of gk, the corresponding harmonic in (8) can be simplified to The produced harmonic has an amplitude of gk 2N ν sÎ ν t , a rotation angular speed of ± ν t ν s ω, and an initial phase of ϕ sν s ∓ ϕ tν t . (ii) If xν s ∓ ν t is an odd multiple of k, and the odd multiple is denoted as 2m + 1 where m ∈ Z (Z stands for integral domain), the corresponding harmonic in (8) can be simplified to The produced harmonic has an amplitude of kN νsÎν t 2| sin 2m+1 2g π | , a rotation angular speed of ± ν t ν s ω, and an initial phase of ϕ sν s ∓ ϕ tν t − arctan cot 2m+1 2g π .
(iii) If the relationship of xν s ∓ ν t with g and k does not belong to the above-mentioned cases, no harmonic will be produced. Table 1 summarizes the analytical quantities of asymmetrical phase winding MMF harmonics. Since symmetricalphase windings can also be seen as asymmetrical-phase windings with only one set, the results in this section are applicable to both symmetrical-phase and asymmetricalphase windings.

D. COMPUTATION OF PLACEMENT MULTIPLIER
Even though the exact value of x will not affect the conclusions, a deterministic way to compute it is essential for the practical application of this analytical method.
In the g-set k-phase windings, let q represent the slot number and p represent the pole pair number. Define h as follows: The minimum mechanical phase shift span τ ∈ Z + (Z + stands for positive integral domain), expressed in the number of slots, satisfies the following linear congruence equation: This means that the left-hand side ( 2π q pτ ) and right-hand side ( 2π h ) of ≡ have the same remainder on division by the modulus (2π ). By eliminating all fractions, (12) can be rearranged as hpτ ≡ q (mod hq) (13) According to modular arithmetic, the necessary and sufficient condition for the existence of solutions to (13) is expressed by the following condition. In this condition, | denotes exactly divides and GCD is the function that calculates the greatest common divisor of its arguments: Assuming that condition (14) holds, the linear congruence equation (13) can be rewritten into its equivalent binarylinear-equation form, with respect to the integers τ and y: To find a special solution of τ , the extended Euclidean algorithm is first used to find a special pair of integer solutions a and b that satisfies another binary linear equation: Equation (16) can then be further transformed into By comparing (15) with (17), a special solution for τ , denoted as τ 0 , can be directly obtained: The minimal positive integer solution for τ , denoted as τ m , can be expressed as: In (19), mod is the modulo operation that results in a value with the same sign as the divisor q GCD(p,q) . This ensures that τ m is positive. The value of x can then be easily determined, as expressed in (20): The whole procedure for determining x is illustrated in Figure 1.
Compared to calculating x, verifying it requires less computational cost. From equation (20), it is known that τ m must equal x q h . However, as shown in (19), τ m is the remainder of a modulo operation and must be smaller than the divisor q GCD(p,q) . This indicates that the following condition must be satisfied: Considering that τ m is a special solution of τ , all possible values of τ can be expressed as shown in equation (22), where m ∈ Z: By substituting (22) into (13), a simple form of linear congruence equation can be derived: Therefore, equations (14), (21), and (23) are the necessary and sufficient conditions to verify the existence and correctness of x.

E. USAGE OF ANALYTICAL METHOD
The proposed analytical method can be of help to rapidly predict the characteristics of the resultant MMF harmonics produced by various armature windings, or, to reverse-reason and trace the source of a specific MMF harmonic. A general procedure to use this analytical method is summarized as following: 1) Acquire the set number g, phase number per set k, polepair number p and slot number q of the target machine; 2) Calculate the placement multiplier x using g, k, p and q according to Figure 1; 3) If qualitative analysis is sufficient, acquire the single-phase winding configuration of the target machine by referring to its winding connection diagram or by dismantling it. Then, qualitatively draft the MMF distribution curve of one phase with a unit current flowing through it. If quantitative analysis is desired, a relatively accurate single-phase MMF distribution can be calculated using the magnetic circuit method or FEA. The curve can also be acquired through experiments. The obtained single-phase MMF distribution is denoted as N ; 4) Use the given current, or measure and save the current data of a single phase at the desired working point. The single-phase current is referred to as I ; 5) Apply the Fourier transformation to N and I . The harmonic orders, amplitudes, and initial phases of N are denoted as ν s ,N ν s and ϕ sν s , respectively. Similarly, the corresponding quantities of I are denoted as ν t ,Î ν t and ϕ tν t , respectively; 6) Observe the significant harmonic components of N and I , and check whether the orders of them satisfy the conditions presented in Table 1 with the previously introduced g, k and x. The amplitudes, initial phases and angular speeds of the resultant MMF harmonics can also be accordingly computed or decided. The entire process can be implemented into computer-aided design programs to provide support for machine design, prototyping and defect diagnosis.

III. VALIDATION WITH FEA SIMULATIONS AND EXPERIMENTS
FEA simulations and experiments are conducted to verify the proposed analytical method. Practical factors must be taken into account when designing the experiments.

A. DESIGN OF EXPERIMENTS
As this study aims to calculate the harmonics in the resultant MMF, rich harmonics, including even and odd ones in both space and time, are desirable to thoroughly test the conclusions demonstrated in Table 1 and cover all the possible cases. Therefore, fractional-slot windings that generate multiple pole pairs and low-order harmonics are preferred for constructing spatial harmonics. As for temporal harmonics, although step currents rarely appear in armature windings in applications, they possess rich harmonics with orders that can reach infinity along the spectrum, and can produce greater low-order harmonics with relatively smaller current. Hence, step currents are selected to produce a sufficient number of temporal harmonics. Additionally, step currents consist of only several discrete constant values, thereby facilitating the design and implementation of experiments.
Since the proposed analytical method is applicable to machines with arbitrary phases, without losing practicality, two types of frequently studied and used multiphase windings, namely dual-three-phase windings (denoted as machine A) and five-phase windings (denoted as machine B), are selected as the experimental units. Three-phase windings are not used as an example, since these types of windings have been well studied and their study is relatively simple. However, it should be noted that the proposed analytical method works for three-phase windings, and even other winding structures, though no related details are discussed here.  Since the investigated MMF is composed only of armature windings, and to keep the magnetic path intact, the rotor is designed as a solid iron cylinder without magnets.
Since MMF cannot be measured directly, the magnetic field, which can also be viewed as MMF spatially modulated by geometric variations and unevenly distributed air gap permeability, is chosen as the metric for measurement and calculation instead. However, to mitigate the potential effect of permeability saturation on the field, the maximum value of currents carried in the windings is carefully restricted.
The machine parameters used in the experiments are shown in Table 2. The models and their parameters used in the FEA simulations are identical to those used in the experiments.
In order to measure the magnetic field in the air gap, an SS49E Hall effect sensor is placed on the surface of the  rotor. Using an STM32H7 MCU, the analog output of the sensor is converted into a digital signal and then collected. When the rotor is manually rotated slowly, the air gap field measured by the sensor, along with the rotor position representing the space position β measured by a quadrature encoder, are simultaneously collected and transmitted to a PC using the UART communication protocol. The data is then decoded and processed with MATLAB. The complete experiment diagram and experiment set are illustrated in Figure 4 and Figure 5, respectively.
The basic experiment procedure can be summarized as follows: 1) Set up the experimental apparatus and ensure that all electrical and mechanical connections are in place; 2) Choose the first of the states as shown in Figure 8a and Figure 9a; 3) Adjust the current outputs of the DC/DC regulators (i.e., current supply) to match the chosen state; 4) Manually rotate the rotor slowly, collect measured air gap field data along with rotor position data, and transmit the data to the PC; 5) Move on to the next state, return to step 3), and continue from there, until all states have been traversed; 6) Process all the received data.

B. INPUTS OF ANALYTICAL METHOD
First, when a unit DC current flows through the first phase of the windings, the resulting magnetic fields represent the       Moreover, it is observed that in addition to odd harmonics, a considerable amount of even harmonics (4th, 2nd, 12th . . . ) are also introduced by the five-phase windings, but not by the dual-three-phase windings.
The selected step phase currents, which depict the temporal characteristics, are presented in Figure 8 and Figure 9. The entire temporal period consists of 5 or 6 stable states marked with circled numbers, as shown in Figure 8a and Figure 9a. In each state, the current remains constant. The selected 60188 VOLUME 11, 2023 Authorized licensed use limited to the terms of the applicable license agreement with IEEE. Restrictions apply.   Figure 8b and Figure 9b show that both even and odd harmonics are injected into the phase current in these two machines.

C. VALIDATION OF ANALYTICAL METHOD
The resultant air gap fields obtained from FEA simulations and experiments are shown in Figure 10 and 11, for machine A and machine B, respectively. The circled numbers in both Figure 10a and Figure 11a correspond to those states in Figure 8a and Figure 9a. In each state, the air gap field varies along the space position β, reflecting spatial harmonics. At the same time, for any given space position β, the field changes from state 1 ⃝ to 2 ⃝, 3 ⃝, and so on, reflecting temporal harmonics. Therefore, Figure 10a and Figure 11a represent both spatial and temporal harmonics of the resultant MMFs produced by multiphase currents in the two types of windings, respectively. Again, the simulation results agree very well with test results, as the curves overlap each other.
A 2D Fourier Transform is used to extract amplitude and phase angle information from the resultant magnetic field harmonics for clearer comparisons in the frequency domain. Since the fields acquired from simulations are nearly identical to those obtained from experiments, only single-phase air gap fields per unit current obtained from experiments are used as inputs for the analytical method. Correspondingly, outputs of the analytical method are compared with experiment results. Figure 10b and Figure 11b compare the amplitudes of the resultant air gap fields in machine A and machine B, respectively. The primary harmonic components have a maximum relative error of 1.9% and a maximum absolute error of 0.3 mT. Additionally, Figure 10c and Figure 11c show the absolute errors of phase angles of resultant air gap fields in machine A and machine B, respectively. And the absolute errors of the primary harmonic components are below π 120 . It can be observed that the results calculated by the analytical method closely match those obtained from experiments, providing support for the quantitative relationship summarized in Table 1.

IV. CONCLUSION
This paper presents an analytical method for investigating MMF harmonics of armature windings in AC machines. The complex but comprehensive relationships between spatial harmonics of winding functions, temporal harmonics of currents, and the spatio-temporal harmonics of MMFs which are generated by them in either symmetrical-phase or asymmetrical-phase windings with an arbitrary number of phases are mathematically derived and generalized. The complete procedure for deriving and using this analytical method is presented in detail, and its correctness is fully verified with FEA simulations and experiments. The analytical method is versatile for analyzing MMF characteristics of AC machines with any type of winding structure, and has practical implications to help design, test, or diagnose AC machines. In the future, this analytical method has the potential to facilitate related studies involving the use of MMF for further in-depth analysis of AC machines.