Analytical Multiphysics Methodology to Predict Vibroacoustics in PMSMs Combining Tangential Electromagnetic Excitation and Tooth Modulation Effects

The vastly increased range of operating speeds and loads of electric vehicle (EV) traction motors implies that the use of traditional finite element (FE)-based noise vibration and harshness (NVH) optimization methodologies becomes challenging due to extensive computational loads. Thus, reduced-order analytical methodologies have become of crucial significance for fast decision-making at predesign stages. In this work, a novel combination of analytical techniques is utilized, to formulate a multiphysics methodology for electromagnetic (EM) NVH prediction of typical high-speed surface-mounted permanent magnet synchronous machines (S-PMSMs). Calculation of the EM stresses in open-circuit conditions was performed using the 2-D complex permeance (CP) methodology. Vibroacoustic predictions were made using analytical expressions for an equivalent 2-D representing the stator. Unlike traditional analytical methodologies, further refinements were implemented to improve the accuracy of the vibroacoustic calculations without sacrificing computational efficiency, through a set of force transformation techniques. These allowed for consideration of the tooth modulation and tangential excitation effects, which are typically neglected in similar studies. The methodology was applied on a 48-slot 8-pole S-PMSM with EM and vibroacoustic results validated numerically. A fast parametric study was performed on the design parameters for the optimization of the generated force harmonics, achieving significant reductions in the sound power levels (SWL) at specific frequencies.


Analytical Multiphysics Methodology to Predict Vibroacoustics in PMSMs Combining Tangential Electromagnetic Excitation and Tooth Modulation Effects
Panagiotis Andreou , Amal Z. Hajjaj, Mahdi Mohammadpour, and Stephanos Theodossiades Abstract-The vastly increased range of operating speeds and loads of electric vehicle (EV) traction motors implies that the use of traditional finite element (FE)-based noise vibration and harshness (NVH) optimization methodologies becomes challenging due to extensive computational loads.Thus, reduced-order analytical methodologies have become of crucial significance for fast decision-making at predesign stages.In this work, a novel combination of analytical techniques is utilized, to formulate a multiphysics methodology for electromagnetic (EM) NVH prediction of typical high-speed surface-mounted permanent magnet synchronous machines (S-PMSMs).Calculation of the EM stresses in open-circuit conditions was performed using the 2-D complex permeance (CP) methodology.Vibroacoustic predictions were made using analytical expressions for an equivalent 2-D representing the stator.Unlike traditional analytical methodologies, further refinements were implemented to improve the accuracy of the vibroacoustic calculations without sacrificing computational efficiency, through a set of force transformation techniques.These allowed for consideration of the tooth modulation and tangential excitation effects, which are typically neglected in similar studies.The methodology was applied on a 48-slot 8-pole S-PMSM with EM and vibroacoustic results validated numerically.A fast parametric study was performed on the design parameters for the optimization of the generated force harmonics, achieving significant reductions in the sound power levels (SWL) at specific frequencies.

I. INTRODUCTION
W ITH demands for electric vehicles (EVs) being at an all-time high [1], the efficiency and performance of electric traction motors are at the center of attention for EV manufacturers.Among others, power density and system weight are two optimization factors with significant importance [2].However, increased power density and weight minimization typically result in adverse effects on the vibroacoustic behavior of electric motors, due to amplifications of the excitation forces in conjunction with reduced system stiffness.Moreover, operation at higher e-motor speeds can introduce higher frequency excitations and induce aggressive vibration and noise.Thus, noise, vibration, and harshness (NVH) behavior and characterization become increasingly important criteria, which must be considered from early design stages.
While aerodynamic and mechanical excitations can contribute to noise generation, the main components of noise generated by EV traction motors are of electromagnetic (EM) origin [3], primarily due to the Maxwell stresses induced at the interface between the air gap and the ferromagnetic materials [4].These stresses are acting on both the stator and the rotor surfaces and are a result of the varying distribution of the magnetic flux density within the e-motor, necessary for electromechanical conversion.
Therefore, when performing a vibroacoustic study on an electric motor, it is necessary for a multiphysics modeling methodology to be adopted, using the topology and operating conditions of the machine as inputs, hence predicting the excitations and the NVH response of the system.It is common for the workflow to be divided into several modules, with data transferred between modules in a sequential, weakly coupled approach.A typical workflow for this purpose is shown in Fig. 1, consisting of four modules that follow the causal path from EM excitation to noise generation [5].
There have been several reported studies in the literature investigating the NVH of electric motors, with varying degrees of complexity concerning each of the four modules.Recently, the use of numerical techniques both for the EM resolution through finite element (FE) modeling as well as acoustic FE or boundary element (BE) modeling has become widespread [6], [7], [8], [9], [10], [11], [12], [13], [14] for all motor topologies [6], [7], [8], [9], [10], [11], [12], [13], [14].This is due to their higher accuracy and ability to account for more complex effects such as local saturation, hysteresis, material nonhomogeneities, and anisotropies [4].Studies have also been performed using numerical techniques to observe the effects of dysfunctions including eccentricities [15], uneven magnetization of the permanent magnets [16], as well as short circuits between coil turns [17], illustrating their impacts on the spatial and temporal distribution of the excitation forces on the structure, as well as the variations in acoustic behavior.
Nevertheless, at a predesign stage, where fast decisions regarding design parameters can give a competitive edge to manufacturers, numerical methods can be cumbersome due to their computational demand for three reasons [4], [18].
1) The time discretization should be small enough to capture higher time harmonics, as well as satisfy the Nyquist criterion, which requires that the sampling frequency is at least twice the maximum frequency that is to be observed, to avoid aliasing [4], [18].
2) The space discretization should follow the same reasoning, to avoid aliasing in the higher wavenumbers in the spatial domain [4], [18].3) In speed sweeps, the speed curve should be discretized in such a way that is fine enough such that any potential resonance peaks can be captured correctly [4], [18].
For these reasons, efforts have been made to replace some or all the numerical steps in the workflow with analytical or semianalytical techniques that can generate accurate solutions and give meaningful insights into the origins of the main noise components while remaining computationally efficient.Le Besnerais [19] developed a fully analytical methodology for the prediction of the vibroacoustic behavior of an induction machine (IM).The permeance-magnetomotive force (P-MMF) methodology was utilized for the resolution of the EM problem, allowing for the computation of the radial EM stresses acting on the stator.The stator's structure was modeled as a 2-D ring and the sound radiation as a result of the radial EM stresses was estimated.The methodology was later extended to consider the effects of the pulsewidth modulation (PWM) harmonics on airborne noise generation [20].A similar methodology was used by Zhu [21] for the prediction of noise and vibration of brushless permanent magnet dc motors, with a particular interest in how the masses of the teeth, windings, and endplates affect the prediction of the natural frequencies through the 2-D ring model.More recently, the P-MMF approach has been adapted to consider saturation effects [22], skewing [23], and eccentricities [24].In these earlier models, the main limitation was the inability to compute the tangential component of the magnetic flux density and, hence, the vibrations of the stator induced by the tangential stresses.However, other studies have identified that the effects of the tangential stresses are nonnegligible [5], [25].This led to the development of more advanced methods, which allowed for the calculation of both the radial and tangential components of the magnetic fields.Devillers [4] developed a methodology utilizing a subdomain model to improve the efficiency of the EM resolution of both a surface permanent magnet synchronous machine (S-PMSM) and a squirrel-cage induction machine (SCIM) while maintaining high accuracy compared to EM FE modeling.Other subdomain models allowed for the effects of skewing [26] and finite iron permeability [27] on the magnetic fields to be investigated.Magnetic equivalent circuit (MEC) models are another powerful tool used for the derivation of the EM fields in electric motors and have been applied to several machine topologies reported in the literature, including S-PMSMs [28], [29], [30], Interior-PMSMs [30], [31], switched reluctance motors (SRMs) [32], and IMs [33].Despite their advantages, the computational efficiency of subdomain and MEC models is significantly reduced when higher spatial and temporal harmonics are to be considered, as the number of unknowns that need to be numerically determined becomes rather large.
Zarko et al. [34] developed another methodology for the calculation of EM fields in brushless dc (BLDCs) motors and S-PMSMs, referred to as the complex permeance (CP) method, where the resultant magnetic flux density in the e-motor is computed as the product of the analytically calculated magnetic flux density in an equivalent slot-less topology and a complex modulating permeance function [35].The CP function was obtained through the use of Schwartz-Christoffel conformal transformations and generated both radial and tangential components of the magnetic flux density.Further developments in CP have allowed for the application of the method to I-PMSMs [36], [37], [38], [39], [40] showing good accuracy when compared to finite element models (FEMs) and are also within acceptable deviation from experimental data.In more recent work by Tessarolo and Olivo [41] and Wang et al. [42], analytical expressions were developed to allow for the computation of the CP radial and tangential components for PMSMs, based on exact motor geometrical features, in Cartesian and cylindrical coordinates, respectively.This eliminated the need for numerical conformal transformation procedures, making the method substantially more computationally efficient.Although the CP method proved its potential in accurately predicting the magnetic flux density components in PMSMs, no studies have been identified to adopt the vibroacoustic predictions method, other than Fakham et al. [43] who utilized a hybridized version of the approach, combining FEM for the determination of the permeance functions.
Another important assumption typically found in analytical methodologies is that each air-gap force density harmonic could only excite a stator's vibration mode with the equivalent spatial order.This suggests that higher force harmonics could be neglected as vibration amplitudes would be significantly smaller [44].However, recent attention to this matter has shown that similar to the Nyquist theorem, when higher force harmonics are projected onto the stator structure, an aliasing effect occurs, referred to as tooth modulation effect.This results in unexpected vibroacoustic contribution by higher force harmonics, as they result in the excitation of lower vibration modes [45].
To improve the accuracy of vibroacoustic predictions of high-speed S-PMSMs performed using computationally efficient approaches, this study proposes a purely analytical multiphysics methodology considering both the additional contribution of the tangential forces, as well as tooth modulation effects.In Section II, the CP method employed for the EM resolution of an S-PMSM is presented, along with the main modeling assumptions.In addition, the techniques implemented to consider the tangential and tooth modulation effects on the system's dynamic response are outlined.The expressions used to analytically predict the natural frequencies and the vibroacoustic response of the stator under harmonic excitation are presented in Section III.The methodology is validated in Section IV, using experimental and numerical results from literature, for a 12-slot 10-pole S-PMSM.In Section V a 48-slot, 8-pole S-PMSM is utilized for the demonstration of the aforementioned effects on the force content and noise generation, and the results are assessed against numerical EM and acoustic predictions.Moreover, a parametric study is performed based on a literature case study to demonstrate the capabilities of the proposed methodology.This combination of analytical techniques for the formation of a computationally efficient multiphysics methodology with considerations for tangential excitation and tooth modulation effects on the vibroacoustic response of high-speed S-PMSMs has hitherto not been identified in the literature, to the best knowledge of the authors.Finally, the main conclusions of the research are summarized in Section VI.

A. Complex Permeance (CP)
To formulate the CP methodology, a set of necessary assumptions are introduced, which allows for the Schwartz-Christoffel conformal transformations between the complex planes to be implemented, such that the magnetic flux density in the middle of the slotted air gap can be calculated [34].This procedure involves transformations from an infinite space into a finite space representing the slotted stator, which inevitably deforms the local space around the slot openings.Nevertheless, the introduced errors because of this appear to be minimal as suggested by Zarko et al. [46].The assumptions are given as follows [34]: 1) The flux variation in the axial dimension and end effects is negligible, and thus, the geometry is considered to be 2-D and is described in polar coordinates.
2) The physical properties of the materials are constant and isotropic.
3) The iron core is infinitely permeable and magnetic saturation effects are neglected.4) The stator slots are infinitely deep.5) The mutual influence between slots is ignored.Based on the above, the radial and tangential components of the flux density are first evaluated in an equivalent slot-less structure, and the slotting harmonics are introduced through the complex modulating function, λ * , evaluated at the middle of the air gap.More precisely [42] where B stands for magnetic flux density in units of tesla (T), λ a and λ b are the radial and tangential components of the CP modulating function, respectively, and j is the imaginary unit.The slot-less magnetic flux density can be expressed as [47] where for np ̸ = 1 and while the coefficient M n can be expressed as [21], [48] M Regarding ( 4)-( 9), p is the number of pole pairs, θ is the mechanical angle in the air gap, ω r is the angular velocity of the rotor in rad s , t is the time, µ 0 is the vacuum magnetic permeability, µ r is the relative recoil permeability of the magnets, R r is the outer radius of the rotor, R m is the radius at the magnets' outer surface (R m = R r + h m , where h m is the magnet thickness), R s is the stator's inner radius, r is the radius at the middle of the air gap, B 0 is the magnet's remanent flux density, and α p is the pole-arc-to-pole-pitch ratio.
For the calculation of the radial and tangential components of the CP, the following analytical expressions were used, derived from the work of Wang et al. [42] where In ( 10)-( 15), coefficient α s0 stands for the slot opening angle, Z s is the number of teeth/slots, g 1 is the perpendicular distance between the rotor's outer surface and the stator's inner surface (g 1 = R s − R r ), and k c is the Carter coefficient with t s being the slot pitch distance and b s the slot opening width.The radial and tangential components of the Maxwell force densities (i.e., with units N/m 2 ) in the middle of the air gap, P r and P t , can then be calculated according to the Maxwell stress tensor (MST) as [49], [50]

B. Tooth Modulation
In conventional analytical methodologies, the Maxwell force densities are directly decomposed into spatial and temporal harmonics using the 2-D fast Fourier transform (FFT) technique.Then, the harmonics are applied onto the cylindrical structure under the assumption that a vibration mode of certain spatial wavenumber can only be excited by a force of the same spatial harmonic.As explained in Section I, this may lead to significant inaccuracies in vibroacoustic predictions, as forces of higher spatial harmonics can indeed excite vibration modes of lower orders.
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.To alleviate this issue and allow for consideration of the tooth modulation effects, the force densities are first transformed into forces acting on the teeth surface, as shown in Fig. 2, similar to the method proposed in [45] where F r,Z and F t,Z represent the radial and tangential forces at the zth tooth, respectively, L is the effective length of the motor, θ is the slot pitch angle, and θ z is the location of the zth tooth.Unlike in [45], the derived expression regarding the twisting moment acting on the tooth surface is not utilized in this study, and as such, a moment could not be applied on the simplified structural model representing the stator.However, the application of the force couple through the transformation explained in Section II-C should theoretically result in an equivalent effect and is thus assumed to be a reasonable approximation.

C. Application of Tangential Forces
The application of tangential forces onto a numerical structural model is generally a straightforward task.However, for the scope of this work, such an approach was ruled out for computational efficiency.Thus, to allow for consideration of the vibrations resulting from the tangential forces, a transformation technique based on the work of Roivanen [25] was implemented.Assuming that the stator yoke curvature and the inertia of the teeth are negligible, the concentrated tangential teeth forces evaluated by (20) can be transformed into a radial force couple acting on the middle surface of the stator yoke where h t and h y are the tooth and yoke heights, respectively, and R yoke is the radius at the center of the stator yoke.The transformation is shown in Fig. 3. Subsequently, the derived forces can be applied in a similar manner on the stator.

III. VIBROACOUSTIC MODELING A. Calculation of Natural Frequencies
Literature suggests that if the stator length-to-mean diameter (d mean ) ratio is L d mean ≤ 1, the stator can be considered "short," and can be represented as a 2-D freely supported ring [51].Hence, classical formulas can be utilized for the calculation of its natural frequencies, considering only circumferential vibration modes.Although the free boundary conditions are not representative of a real setup, the 2-D nature of the model prevents any other conditions from being applied as any displacements and oscillatory motions would be prohibited.Hence, other studies utilizing similar models have also used the same assumption regarding the boundary conditions [4], [21], [51].Considering the frame attached around the stator, the breathing mode frequency f 0 was derived as [19] where E s , K f s , ρ s , and d so denote, respectively, the stator's Young's modulus, stacking factor, material density, and outer diameter, while h f is the frame thickness.The coefficients W st , W sw , W f , and W sy represent the masses of the stator teeth, windings, frame, and stator yoke, respectively.Due to the nature of the system, no excitation of the first bending mode is expected.For higher circumferential modes, m > 1, the natural frequencies were computed by [51] where and K m is a correction factor described in [52] and [53].Note that the modal damping, ξ m , was calculated through the empirical expression [51] Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

B. Stator Vibrations
Following the 2-D decomposition of the concentrated forces, a set of force harmonics with distinct spatial (m) and temporal (ω) orders was obtained, F m,ω .However, in order to compute the circumferential displacements resulting from each harmonic, it was necessary to switch back to force densities.This was performed by assuming that Hence, the static displacements resulting from each force harmonic, S static m,ω were computed as [19] To account for resonance effects between the force harmonics and the stator/frame structure, dynamic deflections were computed by The dynamic displacements were then used to obtain the spatially averaged mean square velocity over the radiating surface [54], where f ω is the temporal frequency of the force harmonic in Hz.

C. Acoustic Modeling
One of the main assumptions imposed at the beginning of the study was that the axial variations in the flux density were negligible, allowing for a 2-D resolution to be carried out.Consequently, the stator displacements, and hence acoustic power calculations, were also performed considering only the effect of circumferential vibrations.Thus, the sound power levels (SWL) generated from the vibrating cylindrical structure due to each vibration mode m were computed as where ρ 0 is the density of air, c 0 is the speed of sound, and S f is the outer surface of the vibrating structure (S f = 2π R frame L).The function σ m defines the effectiveness of the vibrating structure to generate airborne noise and is referred to as modal radiation efficiency.In this case, the modal radiation efficiency of a simply supported and infinitely long cylinder will be used (see Appendix), although similar studies have also used expressions for a spherical structure [19].

IV. METHODOLOGY VALIDATION
For validation purposes, the methodology was used to model S-PMSMs from case studies found in the literature and hence assess the accuracy of the EM and vibroacoustic modules outlined above.In Devillers's work [4], a hybrid subdomain model was developed and used to predict the EM distribution of a 12-slot, 10-pole prototype S-PMSM.The results were assessed against nonlinear EM FEM.In addition, a test rig was built to examine the vibroacoustic response of the system experimentally under open-circuit conditions.
A comparison of the magnetic flux density components in the middle of the air gap obtained by the proposed CP modeling and the nonlinear FEM performed in [4] is presented in Fig. 4. The analytical predictions of the radial and tangential show good agreement with the numerically predicted data, with a mean deviation of 3.3% and 7%, respectively.It must be noted that the maximum deviations occur near the peaks in the radial and tangential flux density components and can reach up to 11% and 45%, respectively.These are a result of the spatial deformations arising from the conformal mapping techniques used in the CP methodology, in the sensitive areas near the interaction of the magnet edges and the stator slot openings, as suggested by Žarko et al. [46].Nonetheless, while the maximum errors are significant, the numerically and analytically predicted components of the flux density show for the biggest part a very high conformity and, hence, were considered acceptable for the purposes of this study, especially when taking into account that the radial component is considered the primary contribution to airborne noise generation.
Furthermore, the analytical methodology was used to estimate the acoustic response of the prototype motor at 1200 RPM such that comparisons could be made with experimental and numerical data provided in [4].Initially, the 2-D  ring model was used to predict the natural frequencies of the system.These are presented in Table I and are assessed against structural FEA predictions from [4].The results suggest that the methodology can predict the natural frequencies of the system with an acceptable accuracy (<10%).
The method can correctly predict the main resonance at 720 Hz corresponding to the second circumferential vibration mode of the stator, as well as the two noise peaks at 200 Hz (2fs) and 1000 (10fs) that appear in the experimental sound pressure level (SPL) measurements presented in [4].However, a direct comparison between the amplitudes of the peaks cannot be performed, as the analytical model predicts SWL, whereas the experimental data show SPL.Therefore, this can be considered a qualitative validation.One should note that the experimental measurements presented in [4] (represented as red stars in Fig. 5) suggest some contribution at 850 and 1150 Hz due to eccentricity and axial vibrations, respectively.These cannot be captured by the adopted 2-D symmetrical model.

V. APPLICATION OF METHODOLOGY ON 48-SLOT 8-POLE S-PMSM
The methodology presented is applicable to any S-PMSM in open-circuit conditions and is not limited to the 12-slot 10-pole motor used for validation purposes.However, the slot/pole combination of the aforementioned machine is not optimal, as tooth modulation can introduce low-order excitations and result in aggressive NVH behavior.For EVs, it is typical for the traction motor to have a high number of stator slots to minimize torque ripple, as well as for the number of stator slots to be a multiple of the number of poles for NVH suppression [55].For this reason, the methodology was applied to a 48-slot, 8-pole S-PMSM to further investigate the tooth modulation and tangential excitation effects on the vibroacoustic response of the typical traction electric motor, which is shown in Fig. 6.
In parallel, a numerical multiphysics modeling workflow was developed, aiming to first assess the quality of the results predicted using the proposed methodology and second provide a direct comparison of the computational demands of the two workflows.Ansys Maxwell 1 [56] was utilized for the prediction of the induced EM stresses in 2-D, under open-circuit conditions and a rotor speed of 5000 RPM.The prediction of the system's vibration modes and frequencies was performed using Abaqus CAE [57].The parameters used for the development of these models as well as the results for the first five natural frequencies predicted both analytically and numerically are shown in Tables II and III, and shown in Fig. 7. To account for the different materials of the stator/frame structure in the analytical prediction of the system's natural frequencies, equivalent material properties were utilized based on the volume fractions of the two components, as suggested in [58] where V s , V f , and V total represent the stator, frame, and total volume, respectively.The dynamic response of the structure under the computed EM stresses was evaluated using AVL Excite 1 Power Unit [59], and by utilizing the predicted surface velocities, the generated SWL was estimated in AVL Excite 1 Acoustics [60].The 3-D numerical model generated for the structural and vibroacoustic simulations consisted of 378 000 nodes and 71 500 elements.Furthermore, the distribution of the radial and tangential components of the magnetic flux density in the middle of the air gap is shown in Fig. 8(a) and (b), respectively.The CP predictions show excellent conformity with the EM FEA data and highlight the methods' capabilities regarding the scope of this investigation.
To illustrate the effects of tooth modulation, the 2-D decomposition of the predicted radial Maxwell force densities is presented in Fig. 9(a), while Fig. 9(b) presents the 2-D decomposition of the forces having used (19), (20), and (30).The most prominent difference is that the spatial harmonic content becomes much lower as a result of the integral.More specifically, the force harmonics are limited to half the 1 Trademarked  number of teeth, in this case ±24.In addition, the 48th spatial harmonic at 12 f s , which is directly linked with cogging torque, as well as its multiples, is modulated to a zero-order spatial harmonic, as they act equally on each tooth in the space and time domains.Consequently, these result in the excitation of the stator's breathing mode.Note that in Fig. 9(a), there Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.exist force contributions at multiples of 2 f s throughout all spatial harmonics.This is not a characteristic of such motor topologies and therefore can only be considered numerical error.Similar observations can be made for the tangential stresses.
Although small in magnitude, it appears that these numerical errors in the force spectrum result in significant overestimations in the SWL predictions using the conventional analytical method, especially at lower frequencies where smaller excitations tend to have a higher impact on noise generation.As shown in Fig. 10, when the force densities are used for the prediction of noise generation, prominent peaks appear at approximately 110 and 530 Hz, which is a result of excitation of the first two vibration modes of the stator.However, in this ideal representation of the system where eccentricities and asymmetries are not considered, such excitations, and hence noise components, should not exist or be as pronounced.In addition, strong contributions are present at all multiples of 2 f s (= 666.67Hz at 5000 RPM).In particular, for the fundamental harmonic at 666.67Hz, caused by an 8th spatial order excitation, it is expected that the sound power generation should be much smaller in magnitude, due to the low sound radiation efficiency of the specific vibration mode in that frequency region.When tooth modulation effects and subsequently the correct spatial harmonics corresponding to those excitations are considered, these noise components are reduced/eliminated.Finally, when tangential forces are also included in the acoustic predictions through the transformation technique outlined in Section II, no additional noise components are introduced as expected.However, the noise amplitudes of the harmonics are altered depending on the phase difference between the tangential and radial force waves, which may result in constructive or destructive interaction.Based on the above, a speed sweep was performed, and acoustic predictions were produced from 200 to 15 000 RPM, at intervals of 25 RPM.Simulations were performed on a 40core computer for a run time of 1 s per speed increment, with a sampling rate of 100 kHz.The total computational time was 6:25 h, averaging only 39 s per speed interval, for the resolution of the EM and vibroacoustic problems.On the contrary, the equivalent numerical workflow consisting of the EM, structural, and vibroacoustic modules necessitated approximately 12 h of total simulation time-excluding model development-for a single speed point.This accentuates the capabilities and effectiveness of the reduced-order methodology proposed in this study for fast predictions in early development stages.The predicted SWL is shown in Fig. 11, highlighting the presence of the expected harmonics in the airborne noise spectrum [3].In addition, the natural frequencies of the system appear as expected, as vertical lines on the spectrum.

A. Parametric Study
Finally, a parametric study was performed along the lines of the work of Wang et al. [61] as an attempt to optimize the airborne noise generated by the motor.The two parameters under investigation were the pole-arc-to-pole-pitch ratio and the slot opening width.The parameters were examined individually to replicate the results of [61], with the pole-arcto-pole-pitch ratio varied from 0.5 to 0.94 and the slot opening width varied from 1 to 6 mm.
Examination of the harmonic content of the radial force density shown in Fig. 12(a) suggests that an increase in the ratio results in an almost linear decrease in the amplitude of the 8th harmonic.This harmonic dominates the vibration and noise, implying that a greater ratio will significantly reduce noise generation.Nevertheless, almost all other harmonics follow an opposing trend, and thus, a compromise has to be made in order to avoid higher noise generation at higher frequencies.Noise generation is not affected by the increasing trend of the zeroth harmonic, as it only produces static deformations of the stator.As the simulations were performed in opencircuit conditions, the constant component of the tangential  force density shown in Fig. 12(b) is always zero.Therefore, the effects of varying the ratio on torque production cannot be observed.However, it appears that by increasing the ratio from 0.8 to 0.94, the 48th harmonic, which is directly linked to cogging torque, increases dramatically.
The effect of the slot opening width on the radial and tangential force density harmonics is shown in Fig. 13(a) and (b), respectively.The most significant effect that can be Fig.14.
Effect of varying the slot opening width on cogging torque peak-to-peak amplitude.Fig. 15.
SWL before and after parametric optimization of the pole-arc-to-pole-pitch ratio and slot opening width, predicted using (a) proposed reduced-order methodology and (b) numerically.
observed is that reducing the slot opening width results in a drastic decrease of the 48th tangential force harmonic.This results in a significant decrease in the cogging torque, which is better visualized in Fig. 14.The effects on the other harmonics are minimal.However, a very small slot opening width may result in adverse effects on the average torque as well as impose manufacturing difficulties and increased costs.
Based on the above results, it was decided that a pole-arcto-pole-pitch ratio of 0.88 and a slot opening width of 2 mm were the optimal design points.The SWL predictions for the original and optimized designs at 5000 RPM are presented in Fig. 15(a).The employed design modifications contributed to a considerable decrease in noise generation at the higher temporal harmonics, with a suggested 4-dB reduction at 10 f s (3.33 kHz), 19-dB reduction at 12 f s (4 kHz), and approximately a 2-dB reduction to the SWL floor above 1500 Hz.
To further assess the validity of the above, the optimization procedure has been replicated using the numerical workflow developed, and the respective results are shown in Fig. 15(b).Upon inspection, one can observe similar trends between the two sets of data, as once again the numerical results Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.
suggest reductions in the SWL at 10 f s and 12 f s , both of 12-dB magnitude, compared to 4 and 19 dB from Fig. 15(a).In addition, the numerical model predicts a 38-dB increase in SWL generation at 4 f s , coinciding with the suggested 22-dB increase from Fig. 15(a).This discrepancy between the magnitudes of the two sets of data mainly arises from the assumptions introduced in developing the structural model.The geometry simplifications, material property approximations, as well as the 2-D nature of the model may not be able to fully capture the complexity of the dynamic model and, however, still provides meaningful insights in the behavior of the system with regard to the main excitation and acoustic harmonics, and the effects of modifying key parameters on their amplitude.
Some resemblance can also be identified with the predictions of Wang et al. [61], and however, an opposing trend can be observed for lower harmonics, with significant amplifications at 4 f s and 8 f s relative to the original design.This phenomenon was not present in [61] and could be a result of the inclusion of the tangential excitations in combination with the tooth modulation effects to the methodology.Nevertheless, the above results highlight the capabilities of the reduced-order methodology proposed, which allows for qualitative predictions to be made at early design stages.This is of crucial importance to the NVH engineer, as it generates an opportunity for preliminary fine-tuning of the electrical machine, such that any resonances between the excitation harmonics and vibration modes of the machine are avoided/minimized.Finally, it must be noted that under open-circuit conditions, it is assumed that the effects of modifying these parameters on the amplitudes of force harmonics are independent of the motor speed, and thus, similar observations are expected throughout the whole speed range of the electric motor.

VI. CONCLUSION
In this work, a reduced-order multiphysics methodology for the EM and vibroacoustic resolution of typical traction S-PMSMs in open-circuit conditions has been presented.The novelty of this work lies in the combination of the analytical EM technique with the refined analytical vibroacoustic model, allowing for consideration of tooth modulation and tangential EM excitations through the implementation of the force transformation methods presented.These advancements resulted in an extremely computationally efficient algorithm with higher capabilities and accuracy than traditional analytical techniques.Once again, it was confirmed that both of these phenomena can have a significant impact on e-motor NVH predictions and should therefore not be neglected in such studies.The presented results were verified through numerical, as well as experimental results from the literature and showed very good conformity.Furthermore, the analytical methodology gives direct insight on the effect of the parameters of interest and is orders of magnitude faster than 3-D numerical methods.
The proposed methodology is orders of magnitude faster when compared to FE-based methodologies.This allowed for a speed sweep up to 15 000 rev/min to be performed on the examined S-PMSM, identifying the main resonances and noise harmonics with a frequency range in excess of 10 kHz.In addition, a fast parametric study was performed, where through optimization of the pole-arc-to-pole-pitch ratio and slot opening widths of the original design, significant SWL reductions were achieved at specific dominant harmonics.The results of the parametric study were also assessed against FEA predictions and showed good qualitative agreement at those harmonics.However, the modifications resulted in adverse effects in lower frequencies.Thus, for fine-tuning of a machine's parameters, the optimization procedure must consider all excitation harmonics as well as any potential resonant excitations over the complete speed range of the e-motor.
To further improve the accuracy and capabilities of the methodology, efforts will be made to include the effects of the armature magnetic field generated under loaded conditions.Additional work will be performed for the prediction of the system's vibration in the time domain, with consideration for the system's boundary conditions and axial vibrations, as they appear to have a significant impact on noise generation.Finally, additional experimental validations would be of crucial importance, especially considering PMSMs with different sizes and pole-slot combinations to further validate the methodology.APPENDIX Modal radiation efficiency for simply supported infinite-length cylinder [3], as shown in the equation at the bottom of the page, where k 0 is the acoustic wavenumber, J m (k 0 R frame ) and J m+1 (k 0 R frame ) are the Bessel functions of the first kind of the mth and (m+ 1)th order, respectively, and Y m (k 0 R frame ) and Y m+1 (k 0 R frame ) are the Bessel functions of the second kind (Neumann functions) of the mth and (m + 1)th order, respectively.
Authorized licensed use limited to the terms of the applicable license agreement with IEEE.Restrictions apply.

Fig. 2 .
Fig. 2. Transformation of air gap force density to concentrated radial and tangential tooth forces.

Fig. 4 .
Fig. 4. Assessment of (a) radial and (b) tangential components of flux density produced by the CP model and compared to the nonlinear FEM results presented in [4], for 12-slot, 10-pole S-PMSM.

Fig. 7 .
Fig. 7. Illustrations of the first five circumferential vibration modes of the stator/frame structure predicted using FEA.

Fig. 8 .
Fig. 8.Comparison of analytically (CP) and numerically predicted (a) radial and (b) tangential components of the air-gap magnetic flux density distribution.

Fig. 10 .
Fig. 10.Comparison of predicted SWL when tooth modulation and tangential effects are considered from (a) 0 to 8000 Hz and (b) 0 to 3000 Hz.

Fig. 13 .
Fig. 13.Effect of varying the slot opening width on (a) radial and (b) tangential force harmonics.

ACKNOWLEDGMENT
The authors would like to express their gratitude to the Engineering and Physical Sciences Research Council (EP/T518098/1-DTP 2020-2021 Loughborough University), the U.K. Research and Innovation (2585407-Electric motor and transmission coupled vibro-acoustics of electric powertrains), and the Loughborough University Doctoral College for their support.Additionally, they would like to thank AVL List GmbH for generously providing access to their software, AVL Excite 1 during the course of this research.For the purpose of open access, the authors have applied a Creative Commons Attribution (CC BY) license to any Author Accepted Manuscript version arising.

TABLE I COMPARISON
OF NUMERICALLY AND ANALYTICALLY PREDICTED NATU- RAL FREQUENCIES FOR 12S-10P MOTOR

TABLE II MOTOR PARAMETERS TABLE III COMPARISON
OF NUMERICALLY AND ANALYTICALLY PREDICTED NATU-RAL FREQUENCIES FOR 48S-8P MOTOR