Analytical model to calculate radial forces in permanent-magnet synchronous machines

In permanent-magnet synchronous machines (PMSMs) electromagnetic forces are identified as the principal source of vibrations and noises. In this article, first of all the sources of magnetic forces in PMSMs are analyzed. Then, an analytical method for the calculation of radial forces on the base of the stator tooth created by electromagnetic sources in PMSMs is presented. Finally finite element method (FEM) is used to corroborate analytical results.


Introduction
Today, electrical machines are considered a mature technology. A good indicator is that in developed cities, more than 65% of the electrical consumption is due to electrical machines [1]. Nevertheless, the need to improve products due to increasing competitiveness and quality standards has pushed forward the study of aspects like vibration and noise [2]. These concerns are especially significant if the application environment is close to human beings, such as elevators, electric traction vehicles, and so on [3].
In permanent-magnet synchronous machines (PMSMs) not under severe working conditions, the fan can be removed. In addition, the aerodynamic noise is to a great extent white noise, which is not perceived as uncomfortable. On the other hand, the noises and vibrations generated by bearings and shaft misalignment can be reduced considerably with the improvement of the manufacturing process. Hence, as shown in Figure 1 and as described in [2][3][4][5], the electromagnetic sources are the dominant ones. In electrical machines that do not use the reluctance principle, the electromechanical energy conversion is based on an interaction between two magnetic fields in the air gap [4]. Through the analysis of these magnetic fields, many of the characteristics of the machine, such as the power, torque, and electromotive forces (EMF), can be obtained [6][7][8].
In [9,10], researchers employed FEM to design electrical machines, and in [9], FEM was also used to calculate magnetic forces. The main advantage of FEM simulations is the high accuracy of the results. However, the computational load is high, and this problem is accentuated when the mesh is not properly defined. In order to avoid this problem, analytical models can be used instead of FEM. The low computational load of the analytical models makes them a very interesting tool, especially when a substantial number of iterations need to be performed, as is usual in design processes for electrical machines. In addition, analytical models give insight into the behavior of electrical machines which is more difficult to obtain via FEM simulations.
Several methods for developing analytical models can be found in the literature. Probably the most extended methods are the lumped parameter models [11,12] and models developed by Fourier series [13,14]. The lumped parameter models can be described as very simplified FEM calculations, since the problem is divided into big elements and several parameters, such as flux sources and reluctances, are associated with each element. Once all elements are defined and the circuit is assembled, the problem is solved by applying numerical methods, as in FEM. Thus, the higher the number of elements the higher the computational load. On the other hand, when using Fourier series, spectra are considered and all harmonics can be described by analytical expressions. These expressions depend on design parameters such as dimensions, number of turns, and remnant fields of magnets, so they give insight into the problem, which is particularly important for designers.
It is well known that the performance of electrical machines, for instance magnetic forces or torque ripple, depends strongly on the flux density in the air gap. Thus, in order to obtain accurate results, it is essential to properly define the spatial distribution of this variable. One possible approximation is to suppose a rectangular distribution. Few researchers employ this method because the distribution is not exactly rectangular and this assumption may lead to significant errors [15]. In order to improve accuracy, some authors have used trapezoidal approximations for direct current (DC) motors [16] or included fringing functions [17,18]. Others have proposed empirical approximations, as shown in [19,20], through some exponential functions. Others have based their models on In electrical machines that do not use the reluctance principle, the electromechanical energy conversion is based on an interaction between two magnetic fields in the air gap [4]. Through the analysis of these magnetic fields, many of the characteristics of the machine, such as the power, torque, and electromotive forces (EMF), can be obtained [6][7][8]. In [9,10], researchers employed FEM to design electrical machines, and in [9], FEM was also used to calculate magnetic forces. The main advantage of FEM simulations is the high accuracy of the results. However, the computational load is high, and this problem is accentuated when the mesh is not properly defined. In order to avoid this problem, analytical models can be used instead of FEM. The low computational load of the analytical models makes them a very interesting tool, especially when a substantial number of iterations need to be performed, as is usual in design processes for electrical machines. In addition, analytical models give insight into the behavior of electrical machines which is more difficult to obtain via FEM simulations.
Several methods for developing analytical models can be found in the literature. Probably the most extended methods are the lumped parameter models [11,12] and models developed by Fourier series [13,14]. The lumped parameter models can be described as very simplified FEM calculations, since the problem is divided into big elements and several parameters, such as flux sources and reluctances, are associated with each element. Once all elements are defined and the circuit is assembled, the problem is solved by applying numerical methods, as in FEM. Thus, the higher the number of elements the higher the computational load. On the other hand, when using Fourier series, spectra are considered and all harmonics can be described by analytical expressions. These expressions depend on design parameters such as dimensions, number of turns, and remnant fields of magnets, so they give insight into the problem, which is particularly important for designers.
It is well known that the performance of electrical machines, for instance magnetic forces or torque ripple, depends strongly on the flux density in the air gap. Thus, in order to obtain accurate results, it is essential to properly define the spatial distribution of this variable. One possible approximation is to suppose a rectangular distribution. Few researchers employ this method because the distribution is not exactly rectangular and this assumption may lead to significant errors [15]. In order to improve accuracy, some authors have used trapezoidal approximations for direct current (DC) motors [16] or included fringing functions [17,18]. Others have proposed empirical approximations, as shown in [19,20], through some exponential functions. Others have based their models on Fourier spatial series. In [21][22][23][24], Zhu created a complete model of a DC machine, including the effect of stator slots. The method was based on a Fourier spatial series resolved in polar coordinates, and took into account the effect of the slots by applying a conformal transformation. It is now possible to find other authors who use similar models [25].
To date, there have been different studies about magnetic forces in which the calculi were based on application of the law of Maxwell's stress tensor and the use of the Fourier series [26][27][28]. In these studies, the influence of the shape of the slots, pole/slot ratio, acoustic behavior, and the unbalance of the system were also analyzed.
Nevertheless, a significant weakness can be detected in the literature, namely the complexity of the developed models. We were unable to find a simple model where the analytical expressions were related to design parameters. Therefore, in these models, the traceability of the problematic harmonics in the magnetic forces that need to be eliminated or reduced becomes a problem. An even greater difficulty is the identification of the design parameters that must be changed to improve the vibroacoustic performances of the machine. In addition, it was noticed that most previous studies were not validated with different machine topologies (integral and fractional), and even fewer in industrially implemented motors. With the aim of addressing this gap, in this study, the Fourier spatial series were resolved in Cartesian coordinates. In addition, the effect of the slots was modeled by applying a similar method to that proposed by Zhu in [23], but in a simpler way. Once the flux density in the air gap was defined, the magnetic forces were computed by applying Maxwell's tensor law. Due to the fact that the developed analytical expressions depended on design parameters, the tracing of problematic harmonics in the magnetic forces that needed to be eliminated or reduced was easy. Therefore, it was possible to identify the design parameters that should be modified to upgrade the vibroacoustic performance of the design. The proposed model is based on previous works presented in [13,14,17]. In addition, it should be pointed that this article is the continuation of [29]. However, in this article, the effect of the eccentricity is modeled, demonstrating the spatial and temporal harmonics originated by each eccentricity type. It is also important to mention that the proposed model was validated using both FEM simulations and experimental measurements performed over two types of machines: a fractional machine comprising 36 slots in the stator and 30 poles in the rotor (known as Qs36p15), and an integer machine comprising 48 slots in the stator and 16 poles in the rotor (known as Qs48p8). Both are commercialized in the elevation sector and their main characterized are summarized in Table 1.

Sources of Magnetic Forces in PMSMs
The functioning principle of a PMSM is based on the interaction of two magnetic fields in the air gap. These fields generate magnetic pressures and forces around the surfaces of the stator and rotor that are adjacent to the air gap. The magnetic forces generated in PMSMs can be created by three different sources: magnetic pressures originating on the stator, magnetic pressures originating on the rotor, and forces created by eccentricities. In this section, all these sources are analyzed.

Magnetic Pressures on the Stator
The flux density on the periphery of the stator is formed by a vector field that can be decomposed into radial and tangential components. These components generate a magnetic radial pressure that can be computed using the law of Maxwell's stress tensor. where µ 0 is the permeability of the vacuum, and B m sr (t, θ) and B m st (t, θ) are the radial and tangential components of the flux density around the periphery of the tooth of the stator.
Due to the high permeability of the material of the stator, the radial component of the flux density is much bigger than the tangential component, so it can be disregarded to obtain the next simplified equation.
The radial component of the magnetic flux density on the periphery of the stator can be computed by Fourier spatial series as shown in (3).
where B m sn are the coefficients of the Fourier series, p is the number of pole pairs, n is the order of the harmonic, ϕ is the initial position of the magnet, Ω m is the mechanical speed, t is the time, and θ is the position around the path.
From (5), Table 2 is obtained, representing the relationship between the harmonics of the flux density and the harmonics of the magnetic pressure. Table 2. Relationship between the harmonics of the flux density and the harmonics of the magnetic pressure on the stator.

Magnetic Pressure Harmonics
It is important to note that all the harmonics of the magnetic flux density are odd. Therefore, all the harmonics of the magnetic pressure are even for any combination of pole pairs and stator slots. In the same way, all the harmonics of the flux density are odd. Table 2 also provides useful information for designers regarding the flux density harmonics that could be minimized or eliminated in order to modify the spectrum of the magnetic force. It is well known that the harmonics of flux density can be modulated by acting upon design parameters such as the magnet pitch or the winding distribution.

Magnetic Pressures on the Rotor
The magnetic pressure on the rotor surface oscillates due to the variation of the flux density that crosses this surface. This variation can be caused by both the magnetic flux created by the winding of the stator and the stator slots.

Oscillating Magnetic Pressure Due to the Stator Slots
In a slotless machine, the air gap reluctance seen in the magnets that are located in the rotor is constant. In this case, when the machine works in open-circuit mode, the flux density in the air gap is created only by the magnets; therefore, the flux density that crosses the surface of the rotor is constant. On the other hand, the flux density crossing the surface of the stator is not constant; it varies with the electrical frequency, which can be obtained multiplying the mechanical frequency by the number of pole pairs.
Focusing on the rotor, as aforementioned, the magnetic flux crossing the surface of the rotor in a slotless machine in open-circuit mode is constant, so the magnetic pressure on the same surface is also constant. Nevertheless, in a slotted stator, the air-gap reluctance is not constant, producing a variation in the magnetic flux that crosses the surface of the rotor. The frequency of this variation is related to the number of slots and pole pairs, as described in (6).
where Q s is the number of slots. In Table 3, the relationship between harmonics of the magnetic flux density and harmonics of the magnetic pressure obtained by (2) is summarized. In this case, it can be deduced that the risky harmonics cannot be removed, but they can be displaced from the problematic frequencies by modifying the Q s /p ratio. The principal component of the magnetic flux created by the winding of the stator propagates along the air gap in synchronism with the main component of the magnetic flux produced by the magnets. That means that it propagates in the same direction and with the same mechanical speed as the rotor, Ω m . Meanwhile, the rest of the spatial components propagate asynchronously with the rotor at Ω n .
Ω n = w s s · n · t p (8) where w s is the electrical speed of the rotor, t p is the periodicity, which is defined as the greatest common divisor of Q s and p, and s = ±1 defines the direction of the propagation. The pulsation of the n-order harmonic of the stator flux density in the rotor reference frame, Ω rn , is the difference between the speed of the rotor, Ω m , and the propagation velocity of this component, Ω n .
Consequently, the frequency of this component referenced to the rotor, f rn , is as follows: where f s is the supplying frequency. Consequently, the order of the harmonics of the magnetic flux density can be computed as follows: From (11), it can be deduced that the harmonic orders of the stator winding flux density referenced to the rotor depend on the p/t p relationship. Therefore, it could be stated that through the modification of this ratio, undesired or critical frequencies could be eliminated. Furthermore, as demonstrated in Section 3.2, the magnitudes of undesired harmonics could also be minimized or even eliminated by adjusting the winding factor, which depends on the chosen winding configuration.
Once the harmonics of the stator winding flux density referenced to the rotor frame are identified (see Tables 4 and 5), the harmonics of the magnetic pressure on the rotor surface can be computed by applying Maxwell's tensor law (2) (see the results in Tables 6 and 7). Table 4. Flux density harmonics created by the armature referenced to the stator and to the rotor in the Qs48p8 machine.

Fixed Reference (Stator)
Moving Reference (Rotor) Order Sign Order    In accordance with the obtained results, it needs to be remarked that in fractional machines (for instance the Qs36p15 machine), fractional pressure harmonics appear on the surface of the rotor. On the other hand, in integer machines (for instance the Qs48p8 machine), only integer and even harmonics appear. Furthermore, in the case of fractional machines, the spectrum of the pressure on the rotor surface presents more harmonics than in the case of integer machines.

Magnetic Forces Due to Eccentricities
It is well known that three types of eccentricities can appear in electrical machines: static, dynamic, or mixed. It is known as static eccentricity when the rotating shaft of the rotor is the center of the rotor, but it is not aligned with the center of the stator (Figure 2a). Dynamic eccentricity occurs when the rotating shaft of the rotor is not aligned with its own center, but it is aligned with the center of the stator (Figure 2b). Finally, mixed eccentricity appears when the previous two cases are mixed. In this case, the rotating shaft of the rotor is aligned with neither the center of the rotor nor the center of the stator (Figure 2c). The effect of the eccentricities can be taken into account by multiplying the magnetic flux density of the air gap in a machine without eccentricity by (12). This distribution function, (13), defines the variation of the air-gap reluctance due to the eccentricity.
where Ω represents the velocity of change of the air-gap length. In the case of dynamic eccentricity, this velocity is the same as the mechanical rotation velocity of the rotor, Ω = Ω , and when there is a static eccentricity this velocity is zero. When there is not any eccentricity the sum term is zero, so the air gap is constant along the whole perimeter, Δ = Δ . Nevertheless, when there is an eccentricity, the complete expression has to be considered. Substituting (13) into (12), (14) is obtained and can be developed to give (15).
From (15), the order of the new magnetic flux density harmonics with both static and dynamic eccentricity are obtained (see Table 8). Table 8. Harmonics of the air-gap flux density in a machine with eccentricity.

Static Eccentricity
Dynamic Eccentricity The effect of the eccentricities can be taken into account by multiplying the magnetic flux density of the air gap in a machine without eccentricity by (12). This distribution function, (13), defines the variation of the air-gap reluctance due to the eccentricity.
where Ω represents the velocity of change of the air-gap length. In the case of dynamic eccentricity, this velocity is the same as the mechanical rotation velocity of the rotor, Ω = Ω m , and when there is a static eccentricity this velocity is zero. When there is not any eccentricity the sum term is zero, so the air gap is constant along the whole perimeter, ∆ = ∆ 0 . Nevertheless, when there is an eccentricity, the complete expression has to be considered. Substituting (13) into (12), (14) is obtained and can be developed to give (15).
From (15), the order of the new magnetic flux density harmonics with both static and dynamic eccentricity are obtained (see Table 8). By applying the simplified equation of Maxwell's stress tensor law (2) to the expression (15) of the flux density, the magnetic pressure including the eccentricities can be computed.
Analyzing (8), it can be deduced that the new harmonics of the magnetic pressure originating from the eccentricities are the ones shown in Table 9. Table 9. Harmonics of the magnetic pressure in a machine with eccentricity.

Static Eccentricity Dynamic Eccentricity
Spatial distribution The magnetic forces on a region can be calculated by (37), integrating magnetic pressures along the region. The resulting magnetic force is a temporal waveform. Therefore, it can be determined that the harmonics of the magnetic forces on a tooth base of the stator present the orders shown in Table 10. As shown, in the case of dynamic eccentricity, the order of the harmonics is related to the number of poles. This is why we can avoid some problematic frequencies through modification of the number of poles. In static eccentricity, this relationship does not exist. Nevertheless, by modifying the magnet pitch, as explained in Section 3.1, we can reduce or eliminate these harmonics. Table 10. Harmonics of the magnetic force in a machine with eccentricity.

Analytical Model
In this section, the analytical model developed for magnetic force computation is described. The proposed analytical model is developed in four steps by means of a Fourier series. In the first two steps, the definition of the flux density created by magnets and the coils of the stator in the air gap of a slotless machine is addressed. Next, the effect of the slots is included. Finally, in the fourth stage, two flux densities (flux densities of the magnets and armature) are superimposed to obtain the overall flux density in the air gap under load conditions. In addition, to ease the mathematical development, a section to explain the complex Fourier series employed is included.

Complex Fourier Series
As is widely known, Fourier series can be employed to represent and realize spectral analysis of periodic functions using (18).
where a 0 , a n , and b n are the Fourier coefficients.
To reduce the mathematical complexity, (18) can be developed into (20) via (19). In this case, it should be stated that negative Fourier coefficients appear. However, these negative coefficients are not harmonics, as can be deduced from (19); they are mathematical phasors.

Magnetic Flux Density Created by Magnets in a Slotless Machine
In PMSMs with a rotor of nonconsecutive poles, the magnitude of the flux density created by the magnets in the air gap can be obtained by (21).
where B r is the remnant value of the magnetic fields, h m is the height, g is the length of the air gap, and µ rm is the relative permeance of the magnets.
After the magnitude is computed, the distribution of this flux density can be represented by means of (22).
where B m gn are the coefficients of the Fourier series, which are obtained by the following equation: where β m is the magnet pitch and a is a parameter known as the fringing coefficient. This is employed to adjust the ideal rectangular distribution with the objective of achieving a more realistic waveform, as shown in Figure 3a. In this case, the fringing coefficient has been parameterized as a function of the different dimensions of the machine in order to provide a physical meaning [17][18][19].
where D ri is the inner diameter of the rotor (Figure 3b).
FOR PEER REVIEW 11 of 24 where is the inner diameter of the rotor (Figure 3b). According to the expression (23), adjusting the magnet pitch can reduce or even eliminate certain harmonics of the flux density created by the magnets. Thus, it is possible to reduce the magnitude of some force harmonics in order to improve the vibroacoustic performance of the machine.

Magnetic Flux Density Created by Coils in a Slotless Machine
The magnetic flux density produced by the armature in the air gap is computed by (25), and that on the rotor surface by (26). These expressions are valid only when magnets are placed on the surface of the rotor. In case of inner magnet configurations, the length of the magnetic path seen by the coils is not constant along the air-gap perimeter, which depends on the rotor position. Furthermore, the problem becomes more nonlinear due to the fact that this length of the magnetic path depends strongly on saturations. Therefore, as the nonlinearities are not considered in this model, rotors with inner magnets placed were not analyzed.
, = ⋅ , where is the magnetomotive force created by the coils in the air gap. To obtain the , the distribution of the conductors must be known. This distribution can be defined using the widely known star of slots method. This method is valid for all types of windings: distributed or concentrated, single layer or double layer, and so on. In Figure 4, the stars of slots corresponding to the machines considered in this work are shown. As shown, the winding of the integer machine is defined by distributed coils in a single-layer configuration, even though a double-layer configuration would also be possible. On the other According to the expression (23), adjusting the magnet pitch can reduce or even eliminate certain harmonics of the flux density created by the magnets. Thus, it is possible to reduce the magnitude of some force harmonics in order to improve the vibroacoustic performance of the machine.

Magnetic Flux Density Created by Coils in a Slotless Machine
The magnetic flux density produced by the armature in the air gap is computed by (25), and that on the rotor surface by (26). These expressions are valid only when magnets are placed on the surface of the rotor. In case of inner magnet configurations, the length of the magnetic path seen by the coils is not constant along the air-gap perimeter, which depends on the rotor position. Furthermore, the problem becomes more nonlinear due to the fact that this length of the magnetic path depends strongly on saturations. Therefore, as the nonlinearities are not considered in this model, rotors with inner magnets placed were not analyzed.
where F MM is the magnetomotive force created by the coils in the air gap. To obtain the F MM , the distribution of the conductors must be known. This distribution can be defined using the widely known star of slots method. This method is valid for all types of windings: distributed or concentrated, single layer or double layer, and so on. In Figure 4, the stars of slots corresponding to the machines considered in this work are shown. As shown, the winding of the integer machine is defined by distributed coils in a single-layer configuration, even though a double-layer configuration would also be possible. On the other hand, in the fractional machine the chosen winding consists of concentrated coils in a double-layer configuration, even though a single-layer configuration would also be possible.  After the physical distribution of the conductors is obtained, the winding factor, , can be computed using (27). Therefore, the magnetomotive force per phase and per current unit can be calculated using (28) and (29).
where ℎ is the number of coils per phase, the position angle of the ℎth coil, the coil pitch angle, the number of turns, and the angle at which the conductors are distributed ( Figure 3b). The magnitude of the winding factor for the n-order harmonic depends on the distribution of the conductors, so adjusting the winding configuration can modulate the winding factor value (in some cases, it is also possible to eliminate the winding factor for certain harmonics). Furthermore, for a given winding factor, adjusting the distribution angle can also minimize or eliminate certain harmonics, which is very useful in optimizing the vibroacoustic performance of electrical machines.
Once is calculated, it is multiplied by the current phase (30) and the magnetomotive forces are subsequently summed (31) in order to compute the overall in the air gap.
Finally, (25) is applied to calculate the magnetic flux density created by the armature in the air gap. It is important to note that in three-phase machines, the spatial shift between magnetomotive forces of different phases is 2 /3 radians, as for the currents. This is the reason why (32) can be deduced from (31).
The expression (32) is referenced to a fixed stator frame. In order to compute the magnetic flux density created by the coils on the rotor surface, it is necessary to change from a fixed stator frame to a moving rotor frame. The new expression for the magnetomotive force in the moving rotor frame is the following: After the physical distribution of the conductors is obtained, the winding factor, Kw n , can be computed using (27). Therefore, the magnetomotive force per phase and per current unit can be calculated using (28) and (29).
where h is the number of coils per phase, ϕ ph the position angle of the hth coil, ϕ the coil pitch angle, N the number of turns, and α the angle at which the conductors are distributed (Figure 3b). The magnitude of the winding factor Kw n for the n-order harmonic depends on the distribution of the conductors, so adjusting the winding configuration can modulate the winding factor value (in some cases, it is also possible to eliminate the winding factor for certain harmonics). Furthermore, for a given winding factor, adjusting the distribution angle can also minimize or eliminate certain harmonics, which is very useful in optimizing the vibroacoustic performance of electrical machines. Once F a is calculated, it is multiplied by the current phase (30) and the magnetomotive forces are subsequently summed (31) in order to compute the overall F MM in the air gap.
Finally, (25) is applied to calculate the magnetic flux density created by the armature in the air gap. It is important to note that in three-phase machines, the spatial shift between magnetomotive forces of different phases is 2π/3 radians, as for the currents. This is the reason why (32) can be deduced from (31). The expression (32) is referenced to a fixed stator frame. In order to compute the magnetic flux density created by the coils on the rotor surface, it is necessary to change from a fixed stator frame to a moving rotor frame. The new expression for the magnetomotive force in the moving rotor frame is the following: F an e −j(nt p θ+nW e t) ie jW e t · 1 + 2 cos 2π 3 (n − 1) As can be seen from (32) and (33), there are some harmonics of the F MM which are annulled when all the magnetomotive forces are summed. It can be easily deduced that in three-phase machines, all harmonics of the total FMM that are multiples of 3 do not exist.

Effect of the Stator Slots
The above equations were developed for a slotless machine. However, in a machine with slots, the air-gap permeance variation due to slots must be included. This effect was studied in [23], where a conformal transformation from the Z plane to the W plane was proposed to obtain the relative permeance of a slot.
When the permeance variation of a slot is obtained (see Figure 5a), the effect for the overall machine can be modeled by (34).
where λ 0 is the mean value of the permeance variation, λ n is the nth-order Fourier coefficient, and Θ is the initial position of the slot. In slotless machines λ 0 = 1 and λ n = 0. In order to represent the flux density in a slotted machine, the expressions of a slotless machine for the flux density due to magnets, (22), and due to armature, (32), must be multiplied by the expression of the reluctance variation, (34), as in the following expression for the armature flux density.
As can be seen from (32) and (33), there are some harmonics of the which are annulled when all the magnetomotive forces are summed. It can be easily deduced that in three-phase machines, all harmonics of the total FMM that are multiples of 3 do not exist.

Effect of the Stator Slots
The above equations were developed for a slotless machine. However, in a machine with slots, the air-gap permeance variation due to slots must be included. This effect was studied in [23], where a conformal transformation from the Z plane to the W plane was proposed to obtain the relative permeance of a slot.
When the permeance variation of a slot is obtained (see Figure 5a), the effect for the overall machine can be modeled by (34).
where is the mean value of the permeance variation, is the nth-order Fourier coefficient, and Θ is the initial position of the slot. In slotless machines = 1 and = 0. In order to represent the flux density in a slotted machine, the expressions of a slotless machine for the flux density due to magnets, (22), and due to armature, (32), must be multiplied by the expression of the reluctance variation, (34), as in the following expression for the armature flux density.

Computation of the Magnetic Forces on the Stator
Magnetic radial pressures are calculated by applying the well-known Maxwell tensor's law, as was introduced in Section 2.1. For this purpose, flux densities computed in the inner perimeter of the stator and in the outer perimeter of the rotor are introduced in Equation (2). Once the magnetic pressures are estimated, pressures are integrated over them to obtain the magnetic radial forces in a region. For instance, the magnetic force in a tooth is obtained by (36), while the magnetic force in a rotor pole is achieved by applying the expression (37).
where is the axial length of the machine and is the interior diameter of the stator.

Computation of the Magnetic Forces on the Stator
Magnetic radial pressures are calculated by applying the well-known Maxwell tensor's law, as was introduced in Section 2.1. For this purpose, flux densities computed in the inner perimeter of the stator and in the outer perimeter of the rotor are introduced in Equation (2). Once the magnetic pressures are estimated, pressures are integrated over them to obtain the magnetic radial forces in a region. For instance, the magnetic force in a tooth is obtained by (36), while the magnetic force in a rotor pole is achieved by applying the expression (37).
where L e is the axial length of the machine and D si is the interior diameter of the stator.

Validation
The analytical model was validated by both FEM simulations and experimental tests performed over two prototypes. The configurations of the prototypes were the ones listed in Table 1.

Analytical Results vs. FEM Results
FEM simulations were carried out using the software Flux from ALTAIR ®® . In this software, three paths are defined in order to compute flux densities, pressures, and forces: one in the middle of the air gap, one on the base of a stator tooth, and the last one on the outer surface of the rotor ( Figure 6). With the first path, the flux density in the air gap is calculated, and with the other two paths, flux densities and magnetic pressures are computed. Once the pressures are estimated, the magnetic forces can be calculated by integrating pressures along the region of interest.

Validation
The analytical model was validated by both FEM simulations and experimental tests performed over two prototypes. The configurations of the prototypes were the ones listed in Table 1.

Analytical Results vs. FEM Results
FEM simulations were carried out using the software Flux from ALTAIR ®® . In this software, three paths are defined in order to compute flux densities, pressures, and forces: one in the middle of the air gap, one on the base of a stator tooth, and the last one on the outer surface of the rotor ( Figure 6). With the first path, the flux density in the air gap is calculated, and with the other two paths, flux densities and magnetic pressures are computed. Once the pressures are estimated, the magnetic forces can be calculated by integrating pressures along the region of interest.  Figure 8, the same flux densities but in the load condition are plotted. It can be observed that the effect of the slots as considered by the analytical model was very similar to the effect computed by FEM.
In Figure 9, the armature flux density calculated in the rotor frame is plotted for the Qs48p8 machine. In Figure 9a, the spatial distribution of this variable is shown, while in Figure 9b, the temporal waveform at a point of the surface of the rotor over time is plotted. In Figure 10, the same results are shown for the Qs36p15 machine.
In the fractional machine, there were fractional harmonics in the spectrum of the flux density in the rotor frame (Figure 11b), while in the spectrum of the flux density in the stator frame, there were only integer harmonics. In the case of the integer machine, both spectra of flux densities (in the rotor frame, Figure 11a, and in the stator frame) consisted of integer harmonics.
According to the results obtained, it can be stated that the analytical model is able to compute flux densities in both a fixed stator frame and a rotor moving frame, with rather high accuracy. Considering that the model does not account for nonlinearities, the achieved accuracy can be considered very good. In Figure 9, the armature flux density calculated in the rotor frame is plotted for the Qs48p8 machine. In Figure 9a, the spatial distribution of this variable is shown, while in Figure 9b, the temporal waveform at a point of the surface of the rotor over time is plotted. In Figure 10, the same results are shown for the Qs36p15 machine.
In the fractional machine, there were fractional harmonics in the spectrum of the flux density in the rotor frame (Figure 11b), while in the spectrum of the flux density in the stator frame, there were only integer harmonics. In the case of the integer machine, both spectra of flux densities (in the rotor frame, Figure 11a, and in the stator frame) consisted of integer harmonics.

Magnetic Pressure
In this section, magnetic pressures calculated in both the stator fixed frame and rotor moving frame are reported. In Figure 12, temporal and spatial waveforms of the magnetic pressure on the rotor are shown for the Qs48p8 machine. In Figure 13, the same results are plotted for the Qs36p15 machine. In the case of the Qs48p8 machine, the results obtained by the model were very similar to the ones obtained by FEM. In the case of the Qs36p15 machine, there were a few differences in the magnitudes, but the results can be considered very precise.
In Figure 14, the spectra of magnetic pressures on the rotor are plotted. As was expected after the analysis performed in Section 2, fractional harmonics appeared in the case of the fractional machine, whereas for the integer machine, all harmonics of the pressure were integers.
Looking at the spectra of the magnetic pressures on the stator plotted in Figure 15, it can be concluded that all harmonics in both machines were integers and even. As was previously pointed out, the spectrum of the pressures on the rotor is modified depending on the winding configuration, whereas the spectrum of pressures on the stator does not depend on the type of winding. Hence, frequencies of harmonics are always the same and they do not depend on whether the machine is fractional or not.
In Figures 16 and 17, the magnetic pressures on the stator are shown for the Qs48p8 and Qs36p15 machines, respectively. Notice that pressures were analyzed at different points of the tooth base, in the edge and in the middle of the slot.
In the case of the Qs48p8 machine, in open-circuit mode the pressures were well computed along all the tooth. However, in load conditions, due to the fact that the magnetic field was higher, the edge of the slot was saturated, leading to an error in the magnitude of the magnetic pressure at the edge. Nevertheless, the saturation did not affect the middle of the slot, so the pressure was accurately estimated in that region.
In the case of the Qs36p15 machine, it seemed to be rather saturated even in opencircuit mode so there was an error in the magnitude of the pressures calculated in the edge of the slot. This error was accentuated when load was applied to the machine. However, the magnetic pressure was estimated with high accuracy in the middle of the slot, where the effect of saturation was smaller.
The edge effect due to saturations can be clearly observed in Figures 18 and 19. As shown, the magnitude of flux density increased considerably in both edges of the slot, leading to the errors mentioned in the estimation of magnetic pressures and forces on the stator tooth.
According to the reported results, it is concluded that the model is able to compute magnetic pressures with rather high accuracy, especially taking into account that satura- According to the results obtained, it can be stated that the analytical model is able to compute flux densities in both a fixed stator frame and a rotor moving frame, with rather high accuracy. Considering that the model does not account for nonlinearities, the achieved accuracy can be considered very good.

Magnetic Pressure
In this section, magnetic pressures calculated in both the stator fixed frame and rotor moving frame are reported. In Figure 12, temporal and spatial waveforms of the magnetic pressure on the rotor are shown for the Qs48p8 machine. In Figure 13, the same results are plotted for the Qs36p15 machine. In the case of the Qs48p8 machine, the results obtained by the model were very similar to the ones obtained by FEM. In the case of the Qs36p15 machine, there were a few differences in the magnitudes, but the results can be considered very precise.  In Figure 14, the spectra of magnetic pressures on the rotor are plotted. As was expected after the analysis performed in Section 2, fractional harmonics appeared in the case of the fractional machine, whereas for the integer machine, all harmonics of the pressure were integers. Looking at the spectra of the magnetic pressures on the stator plotted in Figure 15, it can be concluded that all harmonics in both machines were integers and even. As was previously pointed out, the spectrum of the pressures on the rotor is modified depending on the winding configuration, whereas the spectrum of pressures on the stator does not depend on the type of winding. Hence, frequencies of harmonics are always the same and they do not depend on whether the machine is fractional or not. In Figures 16 and 17, the magnetic pressures on the stator are shown for the Qs48p8 and Qs36p15 machines, respectively. Notice that pressures were analyzed at different points of the tooth base, in the edge and in the middle of the slot.

Magnetic Forces
Finally, in this section, the magnetic forces calculated in a tooth of the stator and on the surface of the rotor are compared. In Figure 20, waveforms, and in Figure 21

Magnetic Forces
Finally, in this section, the magnetic forces calculated in a tooth of the stator and on the surface of the rotor are compared. In Figure 20, waveforms, and in Figure 21 In the case of the Qs48p8 machine, in open-circuit mode the pressures were well computed along all the tooth. However, in load conditions, due to the fact that the magnetic field was higher, the edge of the slot was saturated, leading to an error in the magnitude of the magnetic pressure at the edge. Nevertheless, the saturation did not affect the middle of the slot, so the pressure was accurately estimated in that region.
In the case of the Qs36p15 machine, it seemed to be rather saturated even in opencircuit mode so there was an error in the magnitude of the pressures calculated in the edge of the slot. This error was accentuated when load was applied to the machine. However, the magnetic pressure was estimated with high accuracy in the middle of the slot, where the effect of saturation was smaller.
The edge effect due to saturations can be clearly observed in Figures 18 and 19. As shown, the magnitude of flux density increased considerably in both edges of the slot, leading to the errors mentioned in the estimation of magnetic pressures and forces on the stator tooth.

Magnetic Forces
Finally, in this section, the magnetic forces calculated in a tooth of the stator and on the surface of the rotor are compared. In Figure 20, waveforms, and in Figure 21, spectra  Figure 19. Flux density on the tooth base in load condition: (a) Qs48p8, (b) Qs36p15.
According to the reported results, it is concluded that the model is able to compute magnetic pressures with rather high accuracy, especially taking into account that saturation effects were not considered.

Magnetic Forces
Finally, in this section, the magnetic forces calculated in a tooth of the stator and on the surface of the rotor are compared. In Figure 20, waveforms, and in Figure 21 Regarding the results shown in Figure 20, it can be concluded that the waveforms were properly calculated, but there were some differences in magnitudes, especially in load operation mode. This error in magnitudes was due to nonlinearities. As explained in Section 4.1.3, there was a significant concentration of magnetic flux at the edge of the tooth which led to errors in the computation of the magnitude of the flux density in these areas, and, in consequence, it also led to errors in the pressures computed in the edges. However, in the middle of the tooth, the flux density was estimated with rather high accuracy. Thus, this error in the magnitude of the magnetic forces was due to saturation effects that appeared at the side ends of the tooth of the stator. The spectra in Figure 21 show that the model was able to identify the main force harmonics.
In Figure 22, the magnetic forces calculated on the rotor surface are shown. In this case, magnitudes of oscillation were very small in comparison with the average value of the magnetic field, and also in comparison with the forces in the stator. This makes sense because the main magnetic field created by the magnets is constant in the rotor frame, and harmonics in pressures and forces are produced by the armature field, which is much weaker than the field of the magnets.
To summarize, the magnetic forces computed by the developed analytical model can be considered accurate enough, considering that nonlinearities were not taken into account. Regarding the results shown in Figure 20, it can be concluded that the waveforms were properly calculated, but there were some differences in magnitudes, especially in load operation mode. This error in magnitudes was due to nonlinearities. As explained in Section 4.1.3, there was a significant concentration of magnetic flux at the edge of the tooth which led to errors in the computation of the magnitude of the flux density in these areas, and, in consequence, it also led to errors in the pressures computed in the edges. However, in the middle of the tooth, the flux density was estimated with rather high accuracy. Thus, this error in the magnitude of the magnetic forces was due to saturation effects that appeared at the side ends of the tooth of the stator. The spectra in Figure 21 show that the model was able to identify the main force harmonics.
In Figure 22, the magnetic forces calculated on the rotor surface are shown. In this case, magnitudes of oscillation were very small in comparison with the average value of the magnetic field, and also in comparison with the forces in the stator. This makes sense because the main magnetic field created by the magnets is constant in the rotor frame, and harmonics in pressures and forces are produced by the armature field, which is much weaker than the field of the magnets.
To summarize, the magnetic forces computed by the developed analytical model can be considered accurate enough, considering that nonlinearities were not taken into account.   Regarding the results shown in Figure 20, it can be concluded that the waveforms were properly calculated, but there were some differences in magnitudes, especially in load operation mode. This error in magnitudes was due to nonlinearities. As explained in Section 4.1.3, there was a significant concentration of magnetic flux at the edge of the tooth which led to errors in the computation of the magnitude of the flux density in these areas, and, in consequence, it also led to errors in the pressures computed in the edges. However, in the middle of the tooth, the flux density was estimated with rather high accuracy. Thus, this error in the magnitude of the magnetic forces was due to saturation effects that appeared at the side ends of the tooth of the stator. The spectra in Figure 21 show that the model was able to identify the main force harmonics.
In Figure 22, the magnetic forces calculated on the rotor surface are shown. In this case, magnitudes of oscillation were very small in comparison with the average value of the magnetic field, and also in comparison with the forces in the stator. This makes sense because the main magnetic field created by the magnets is constant in the rotor frame, and harmonics in pressures and forces are produced by the armature field, which is much weaker than the field of the magnets.

Experimental Validation
The analytical model was also validated by experimental measurements. Vibrations were measured in both prototypes shown in Figures 23 and 24 by placing an accelerometer on the outer surface of the stator. The Figures 25 and 26 show that the developed analytical expressions characterized the order and the origin of the various electromagnetic vibrations with accuracy. In the integer machine, all harmonics were integer, and in the fractional machine, fractional harmonics were clearly seen. These harmonics originate on the rotor surface and they are transmitted to the stator surface, for example through the bearings.
The magnitudes of the harmonics could not be validated with the obtained experimental results. This issue is a future prospect; it would be interesting to conduct additional research in order to experimentally validate the magnitudes of these magnetic forces.   To summarize, the magnetic forces computed by the developed analytical model can be considered accurate enough, considering that nonlinearities were not taken into account.

Experimental Validation
The analytical model was also validated by experimental measurements. Vibrations were measured in both prototypes shown in Figures 23 and 24 by placing an accelerometer on the outer surface of the stator. The Figures 25 and 26 show that the developed analytical expressions characterized the order and the origin of the various electromagnetic vibrations with accuracy. In the integer machine, all harmonics were integer, and in the fractional machine, fractional harmonics were clearly seen. These harmonics originate on the rotor surface and they are transmitted to the stator surface, for example through the bearings.

Experimental Validation
The analytical model was also validated by experimental measurements. V were measured in both prototypes shown in Figures 23 and 24 by placing an acc ter on the outer surface of the stator. The Figures 25 and 26 show that the develo lytical expressions characterized the order and the origin of the various electrom vibrations with accuracy. In the integer machine, all harmonics were integer, an fractional machine, fractional harmonics were clearly seen. These harmonics orig the rotor surface and they are transmitted to the stator surface, for example thro bearings.
The magnitudes of the harmonics could not be validated with the obtained mental results. This issue is a future prospect; it would be interesting to conduct ad research in order to experimentally validate the magnitudes of these magnetic fo

Conclusions
In this article, an analytical model based on spatial Fourier series for the calculation of radial magnetic forces and pressures in PMSMs is developed. First, magnetic flux densities created in both open-circuit and load operation conditions are represented by Fourier series. The radial magnetic forces and pressures originating from these magnetic flux densities on both the base of the stator tooth and on the rotor surface are then calculated by applying Maxwell's tensor law. The proposed analytical expressions for the calculus of flux densities, magnetic pressures, and magnetic forces are defined as a function of several design parameters, such as dimensions, magnetic properties of the materials, and so on. Hence, the model developed gives insight into the functioning principle of the electrical machines and the expressions developed enable design rules to be defined. Furthermore, the proposed analytical model makes it possible to face iterative calculus, which is very common in the designing process of electrical machines, within a relatively short time while ensuring accurate enough results.
The proposed model was validated by both FEM simulations and experimental tests. Two types of machine were considered: a fractional machine comprising 36 slots in the stator and 30 poles in the rotor, and an integer machine comprising 48 slots in the stator and 16 poles in the rotor. Results given by the model showed good agreement with the FEM results. The results for flux densities, pressures, and forces were compared with The magnitudes of the harmonics could not be validated with the obtained experimental results. This issue is a future prospect; it would be interesting to conduct additional research in order to experimentally validate the magnitudes of these magnetic forces.

Conclusions
In this article, an analytical model based on spatial Fourier series for the calculation of radial magnetic forces and pressures in PMSMs is developed. First, magnetic flux densities created in both open-circuit and load operation conditions are represented by Fourier series. The radial magnetic forces and pressures originating from these magnetic flux densities on both the base of the stator tooth and on the rotor surface are then calculated by applying Maxwell's tensor law. The proposed analytical expressions for the calculus of flux densities, magnetic pressures, and magnetic forces are defined as a function of several design parameters, such as dimensions, magnetic properties of the materials, and so on. Hence, the model developed gives insight into the functioning principle of the electrical machines and the expressions developed enable design rules to be defined. Furthermore, the proposed analytical model makes it possible to face iterative calculus, which is very common in the designing process of electrical machines, within a relatively short time while ensuring accurate enough results.
The proposed model was validated by both FEM simulations and experimental tests. Two types of machine were considered: a fractional machine comprising 36 slots in the stator and 30 poles in the rotor, and an integer machine comprising 48 slots in the stator and 16 poles in the rotor. Results given by the model showed good agreement with the FEM results. The results for flux densities, pressures, and forces were compared with those obtained by FEM and it was demonstrated that the magnitudes and frequencies of the main harmonics were computed with rather good accuracy. However, some differences were observed in the calculated forces and pressures in the stator. This was probably due to the nonlinearity of the problem, which was not considered in the proposed model. This nonlinearity is greater in slotted stators because the flux is concentrated at the edges of the slots saturating these regions. For this reason, significant differences in the pressures calculated at the edge of the slots were obtained, while the pressures computed in the middle of the slots were very accurate.
In addition, two prototypes were built and vibrations were measured under load conditions. According to the results, it can be stated that the proposed model is able to identify the frequencies of the main force harmonics. Nevertheless, with these experimental results, the magnitudes of the harmonics could not be validated. This issue is challenging and as a future research prospect, it would be interesting to conduct additional research in order to experimentally validate the magnitudes of these magnetic forces.
Furthermore, it is important to notice that the developed model estimates magnetic forces not only on the stator surface but also on the rotor surface. Interesting results were obtained when comparing the performances of both machines. In the integer machine, all estimated harmonics were even. In the case of the fractional machine, the harmonics of the force on the stator were even, but those on the rotor were fractional. These results were corroborated by the experimental tests. Vibrations were measured in both machines by placing an accelerometer on the surface of the stator. In the case of the fractional machine, fractional harmonics were clearly identified.
The main advantages of the model compared to the FEM method are the rapidity of the calculi and the ease of relating components of the flux density and magnetic forces to design variables. Therefore, by employing the presented model, an iterative optimization of machine design may be accomplished in a relatively easy and quick way.
Finally, as future prospect, it would be interesting to include in the model the amplitude of the forces created by the eccentricities. In this way, the model would be able to compute all magnetic forces that appear in PMSMs.