A novel sparse reduced order formulation for modeling electromagnetic forces in electric motors

Abstract

A novel model order reduction (MOR) technique is presented to achieve fast and real-time predictions as well as high-dimensional parametric solutions for the electromagnetic force which will help the design, analysis of performance and implementation of electric machines concerning industrial applications such as the noise, vibration, and harshness in electric motors. The approach allows to avoid the long-time simulations needed to analyze the electric machine at different operation points. In addition, it facilitates the computation and coupling of the motor model in other physical subsystems. Specifically, we propose a novel formulation of the sparse proper generalized decomposition procedure, combining it with a reduced basis approach, which is used to fit correctly the reduced order model with the numerical simulations as well as to obtain a further data compression. This technique can be applied to construct a regression model from high-dimensional data. These data can come, for example, from finite element simulations. As will be shown, an excellent agreement between the results of the proposed approach and the finite element method models are observed.

Introduction

The electric powertrain is drastically growing its importance in industry, specifically in the automotive one, because of different reasons. One of them is the environmental and energy regulations laid out in United Nations Climate Change Conference held in Paris in 2015 [1, 14]. In addition, the development of technologies such as electric vehicles (EVs), hybrid electric vehicles (HEVs) or self-driving cars encourages even more the research and development in this area.

Furthermore, an important effort of research in this direction of both companies and national governments (such as USA or China) is carried out in the last years. For example, in Europe, the creation of The European Technology & Innovation Platform (ETIP) on batteries, named BatteRIes Europe [2], clarifies and envisages even more the importance of the electric powertrain and the EVs in our next years.

Due to this interest, it is mandatory to analyze the vehicle noise, vibration and harshness (NVH) because contrary to popular costumers’ belief, EVs are not silent at all. It is true that they can present lower overall levels comparing to internal combustion engines (ICEs) but unfortunately, they have high-frequency and tonal content that makes the electric motor noise annoying if this issue is not addressed correctly [8, 15].

The noise in the EV can be divided in 4 main sources: powertrain noise, wind noise, tire/road noise, and ancillary noise (where the first is one of the most significant). For this reason, it is so important to analyze the electric powertrain, specially, the electric motor in the NVH studies.

As it is shown in different sources [3, 9] , vibration and noise produced by electrical machines can be divided into three categories [10] : electromagnetic vibration and noise, mechanical noise (related to the mechanical assembly, in particular bearings) and aerodynamic noise (they are mainly caused by aerodynamic forces in ventilation components of the motor).

Adding complexity to the analysis, the electromagnetic vibration and noise of electric motors is a multi-physics problem, involving multiple fields including electromagnetism, structural dynamics and acoustics.

Electromagnetic-induced vibration and noise of an electric motor is mainly caused by radial force waves [3] on the stator surface as said in [10] and [7, 11, 24] also agree on this. In addition, this calculation is one of the most highly-time consuming and challenging. Consequently, reduced order models (ROMs) are an appealing alternative to compute these force waves which will be the input for the other stages of the NVH analysis, as illustrated in Fig. 1.

Fig. 1 [Images not available. See PDF.]

Conversion of electric energy into acoustic energy

Contemporary electric motors are designed with higher magnetic flux density in the air gap which produces higher radial magnetic forces acting on the stator. This can lead to a rise of vibration and acoustic problems [3, 10, 24].

Even worse, the relatively small-size and lightweight design for motors, the large electromagnetic forces as well as the poor rigidity of motor structures increase this serious problem of electromagnetic vibration and noise which will affect the riding comfort. Therefore, it is clearly justified to incorporate in the design the requirements of noise and vibration to avoid large retrofit expenses when the overall performance is being optimized/balanced as said in [10] and [7, 9, 14] also agree on this.

A possible gateway for enabling more efficient designs could be the simulation of complex models of the electric machine. However, although it is a valid option in some cases, it can be problematic when a detailed analysis of the machine is being carried out or the machine is beeing simulated with the other physical systems which interact with it.

The reasons are the following ones. On the one hand, a lot of simulations are needed to analyze the NVH behavior of a motor, thus being required faster simulations without losing excessively accuracy. On the other hand, these complex models cannot be used to analyze the whole EV system because of the difficulties to couple the machine model with the other subsystems. Nevertheless, this system simulation is important in the analysis because the real inputs for the motor are computed and the entire vehicle system is modeled (within the so-called system engineering) obtaining the predictive responses of the different subsystems when interacting with each other.

Therefore, the main aim of this paper is to pave the way towards the use of simple, accurate and fast ROMs to predict the electromagnetic forces on the stator surface which is one of the most challenging and highly time consuming steps of the NVH analysis. For this reason, the proposed methodology can be used to compute these force waves in almost real-time (because of the simple algebraic expression to be manipulated) and with accuracy respect to the finite element model (FEM) where it is based.

Concretely, the methodology proposed in this work is a novel formulation of the sparse proper generalized decomposition (sPGD) [12, 19], combining it with a reduced basis (RB) approach, which is used to fit correctly the ROM with the numerical simulations and to perform a further data compression. This procedure is described in Sect. 3.3 and can be extended and used in other industrial problems without issues when appropriate.

The proposed technique can be applied to construct a regression model constructed from data. This data can come for example from a finite element method (FEM) software. Moreover, it is perfectly suitable to create a high-dimensional function to give us the electromagnetic force or pressure considering many parameters, including the geometric ones.

One of the main advantages of this technique is that it deals with the well-known course of dimensionality issue allowing to have good results with only few snapshots, that is, high fidelity solutions. This way, given some snapshots, the regression model can be constructed using mainly polynomials, therefore once the model is constructed the responses are computed immediately because of the simple expression obtained.

Moreover, the accuracy of the regression model will be strongly related to the data from which the regression is carried out. For example, if the data comes from a 3D model, the regression results will be better than the results obtained from a 2D model. Furthermore, the methodology proposed in this paper can easily be extended to other type of motors and modeling frameworks.

To illustrate the process, the sPGD-RB procedure is applied to construct 3 regression models from the 2D FEM models of an induction and synchronous motor presented in Sect. 2. The open source software finite element method magnetics (FEMM) [16] is used to obtain the pseudo experimental data. In Sect. 3, the proposed approach is presented and discussed.

The results obtained as well as a comparison study between the sPGD and FEMM predictions are presented in Sect. 4 to evaluate the accuracy of the ROM.

Finally, in Sect. 5 the general conclusions of the present work are discussed as well as works in progress.

Analyzed motor

Induction motor

We focus on a 2D current-based formulation for a squirrel-cage induction motor. According to [4, 18], the 2D problem is given by the following PDE:

1

$$\begin{array}{cc}& \frac{\partial }{\partial x}(\frac{1}{\mu }\frac{\partial A_z}{\partial x})+\frac{\partial }{\partial y}(\frac{1}{\mu }\frac{\partial A_z}{\partial y})\hfill \\ \hfill & =-J_0+\sigma \frac{\partial A_z}{\partial t}-\sigma \overrightarrow{v}\times (\mathrm{\nabla }\times \overrightarrow{A}),\hfill \end{array}$$where $\overrightarrow{A}$ is the magnetic vector potential, $\overrightarrow{v}$ is the velocity, $A_z$ is the z-component of the magnetic vector potential, $J_0$ is the applied density current source, $\sigma$ is the electric conductivity and $\mu$ represents the permeability.

Considering $J_0(t)=\text{Re}(\widehat{J_0}{e}^{j\omega t})$ , $A_z(t)=\text{Re}(\widehat{A_z}{e}^{j\phi t}{e}^{j\omega t})$ and $\widetilde{A_z}=\widehat{A_z}{e}^{j\phi t}$ (where $\phi$ is the phase angle between $A_z(t)$ and $J_0(t)$) as well as the assumptions and mathematical procedure shown in [4, 18], the above problem can be simplified to the following expression for an harmonic analysis:

2

$$\begin{array}{c}\hfill \frac{\partial }{\partial x}(\frac{1}{\mu }\frac{\partial \widetilde{A_z}}{\partial x})+\frac{\partial }{\partial y}(\frac{1}{\mu }\frac{\partial \widetilde{A_z}}{\partial y})=-\widehat{J_0}+j\omega {\sigma }_{\mathit{eq}}\widetilde{A_z},\end{array}$$where $\widetilde{A_z}$ is a complex number, j is the unit imaginary number, $\omega =2\pi f$, f is the supply frequency, ${\sigma }_{\mathit{eq}}$ is an equivalent conductivity computed as ${\sigma }_{\mathit{eq}}=\sigma s$ and s is the slip. The slip in induction motors is defined as $s=\frac{n_s-n_r}{n_s}$, where $n_s$ is the synchronous speed and $n_r$ is the rotor speed.

This formulation transforms the magneto-dynamic field problem expressed by Eq. (1) to a magnetostatic complex field problem with induced currents. In the chosen approach the rotor is fixed in the stator reference frame and an equivalent conductivity is assigned to the rotor bars to take into consideration the motional term of the current density, that is, the induced current density due to the movement.

Therefore, we can represent the motor at any operation point by multiplying the rotor conductivities by the slip. This is similar to the procedure used in the standard motor equivalent circuit where the rotor resistance is divided by the slip.

It is true that solving Eq. (1) is more accurate than solving Eq. (2) but we did this choice to adapt us to the free FEMM software and its capabilities

To take into consideration the nonlinear relationship B-H in Eq. (2), FEMM includes a nonlinear time harmonic solver that it is used in this work.

This nonlinear time harmonic analysis seeks to include the effects of nonlinearities like saturation and hysteresis on the fundamental of the response, while ignoring higher harmonic content.

There are several subtly different variations of the formulation that can yield slightly different results, so documentation of what has actually been implement is important to the correct interpretation of the results from this solver. An excellent description of this formulation is contained in [13, 16].

Synchronous machine: brushless motor

In this machine, a two-dimensional steady state analysis is carried out to adapt us to the capabilities of the open source FEMM. In three-phase motors, as in the other polyphase configurations of the synchronous machines, the stator-produced magnetomotive force (MMF) rotates at synchronous speed. Since the rotor is also rotating at synchronous speed in the steady state, an observer on the rotor experiences a constant field ($\frac{\partial B}{\partial t}=0$), and therefore, there are no eddy currents on the rotor.

On the other hand, an observer on the stator experiences a time varying field whose fundamental is at the system frequency. Since the stator is laminated and the stator windings are stranded and transposed, the eddy currents are resistance limited and can be neglected in the field computation. Hence the term ($\sigma \frac{\partial A_z}{\partial t}$) in the diffusion equation is neglected also in this frame since $\sigma$ can be considered zero.

If we take into consideration the above assumptions in Eq. (1) as well as adding the modeling term for the permanent magnets, it will lead us to the Poisson’s equation for a magnetostatic analysis:

3

$$\begin{array}{cc}& \frac{\partial }{\partial x}(\frac{1}{\mu }\frac{\partial A_z}{\partial x})+\frac{\partial }{\partial y}(\frac{1}{\mu }\frac{\partial A_z}{\partial y})\hfill \\ \hfill & =-J_0-\frac{\partial }{\partial x}(\frac{B_{r,y}}{\mu })+\frac{\partial }{\partial y}(\frac{B_{r,x}}{\mu }),\hfill \end{array}$$where $B_{r,x}$ and $B_{r,y}$ are respectively the x and y components of the remanent flux density.

A fixed reference frame is used in the above equation where the PDE is solved for each rotor position. To further details in electric machine modeling, we kindly suggest the reading of [4, 18]. Note also that Eqs. (2) and (3) are solved in Cartesian coordinates.

Post-processing step: computation of the radial force waves

Once the PDE is solved, the B field must be obtained from the expression:

4

$$\begin{array}{c}\hfill \overrightarrow{B}=\mathrm{\nabla }\times \overrightarrow{A},\end{array}$$and then, the Maxwell stress tensor is used to compute the forces.

In this work we use, for the sake of simplicity, the approach applied to the path at constant radius in the middle of the air gap [21]. Other possibilities can be envisaged. For example, in [17], the Maxwell stress tensor is used and compared under different paths. However, it is important to highlight that using a more complicated post-processing will not affect the fast computational features of the final ROM obtained in this work.

The normal and tangential components of magnetic pressure are:

5

$$\begin{array}{cc}\hfill P_n(\alpha ,t)=& \frac{1}{2\mu }(B_n^2(\alpha ,t)-B_{\tau }^2(\alpha ,t)),\hfill \end{array}$$

6

$$\begin{array}{cc}\hfill P_{\tau }(\alpha ,t)=& \frac{1}{\mu }{B}_n(\alpha ,t)B_{\tau }(\alpha ,t),\hfill \end{array}$$where $\alpha$ is the angle of a polar coordinate system pointing to the selected air-gap point, subscript n refers to the radial component in the air-gap midline and the subscript $\tau$ refers to the tangential component in the air-gap midline.

In the vibro-acoustic context, simplifying assumptions are often added neglecting the tangential terms. The reason is that the tangential component of the flux density is much smaller than the normal component [11, 17, 24]. This leads to:

7

$$\begin{array}{cc}& P_n(\alpha ,t)=\frac{1}{2\mu }{B}_n^2(\alpha ,t),\hfill \end{array}$$

8

$$\begin{array}{cc}& P_{\tau }(\alpha ,t)\approx 0.\hfill \end{array}$$

Motor parameters

Units used

We introduce the units used in the present work in the following table (Table 1).

Table 1. Units used

Unit	Symbol	Unit	Symbol
Millimeters	mm	Meters	m
Degrees	deg.	Radians	rad
Horsepower	HP	Ampere	A
Hertz	Hz	Revolutions per minute	rpm
Tesla	T	Megasiemens	MS
Pascal	N/m${}^2$ = Pa	–	–

Induction motor

We employed an example from FEMM website [16], where geometry and further details can be found. Main parameters/features of the motor:

• 2 HP motor, 50 Hz, 3-phase supply. It is a 4-pole machine (i.e., p = 2).

• The winding configuration for one pole of the machine is: A+, A+, A+, C−, C−, C−, B+, B+, B+ (the nine slots from 0 to 90 geometrical degrees).

• There are a total of 36 slots on the stator and 28 slots on the rotor. A total of 44 turns sit inside each stator slot.

• The rotor’s diameter is 80 mm, and the air gap between the rotor and stator is 0.375 mm. The length of the machine in the into-the-page direction is 100 mm.

Materials used:

• Aluminum for rotor bars ($\sigma$ = 34.45 MS/m)

• Air (${\mu }_r$ = 1)

• Stator Winding (${\mu }_r$ = 1)

• For the Case 1. Linear B–H relationship.

  [NEST]

  • Silicon Core Iron for the ferromagnetic materials (${\mu }_r$ = 7000)

• For the Case 2. Nonlinear B–H relationship.

  [NEST]

  • Carpenter Silicon Core Iron “A”, 1066C Anneal (B–H curve taken from FEMM library)

Brushless motor

We employed an example from FEMM website [16] where geometry and further details can be found. Main parameters/features of the motor:

• 3-phase supply. It is a 8-pole machine (i.e. p = 4)

• The winding configuration for one pole of the machine is: A+, B−, C+ (the three slots from 0 to 45 geometrical degrees).

Geometry:

• Rotor inner diameter: 22.8 mm.

• Rotor iron outer diameter: 50.5 mm.

• Rotor outer diameter: 55.1 mm.

• Air gap length: 0.7 mm.

• Stator outer diameter: 100 mm.

• Angle spanned by tooth: 11.9 deg.

• Turns/slot: 46.

• The complete geometry of the machine is reconstructed in Fig. 2.

Materials:

• Winding wire: 4X20AWG copper wire

• Magnet material: Sm2Co17 24MGOe

• Stator material: 24 Gauge M19 NGO Steel @ 98% fill

• Rotor material: 1018 steel

Fig. 2 [Images not available. See PDF.]

Brushless machine: complete geometry drawn in FEMM software

Reduced order model

Introduction

As already mentioned in Sect. 1, this paper aims at proposing a new methodology to obtain accurately the electromagnetic forces on the stator surface (specifically the magnetic pressure) in almost real-time for any choice of a given set of parameters. The reason is that this is one of the most challenging steps of the vibro-acoustic analysis as it was discussed in Sect. 1.

In fact, the proposed approach allow us to obtain the force or magnetic pressure immediately when changing different parameters of the problem such as conductivities or the operation point of the motor, for instance. This methodology opens the door to a more efficient vibro-acoustic analysis during the design and optimization process of the electric machine as well as to improve the prediction capacities when the whole vehicle system is considered.

The ROM’s accuracy will depend on the data used for the regression technique. In the cases analyzed here, the models described in Sect. 2 are used but without loss of generality of the foregoing, the proposed approach in this paper can easily be extended to more complex models increasing more the ratio between accuracy and computational efficiency. But not only that, the proposed approach also can be extended and used to other type of problems and modeling frameworks.

As it was introduced in Sect. 1, the ROM is made of a regression combining both the sPGD and the RB techniques. In Sect. 3.2, the sPGD technique is exposed. This methodology will allow us to achieve excellent results when dealing with high-dimensional spaces and sparse data. This technique is specially convenient because only sparse data is available when dealing with high-dimensional problems [12, 19]. This allows us to cope with the curse of dimensionality.

Then, as it will be discussed in Sect. 3.3, because of the presence of localized behaviors and discontinuities in the computed solutions, more than interpolating the solution itself, we consider the construction of a RB which will be used for the regression procedure. The RB is inserted in the sPGD formulation, creating the ROM described in Sect. 3.3.

Sparse PGD

If needed, standard references concerning the Proper Generalized Decomposition (PGD) method for solving PDEs can be reviewed in [5, 6, 20]. In this section a brief exposition of the sPGD is presented based on a more detailed work [19]. After discussing this basis, we will present the novel approach in Sect. 3.3. First, we are going to define the following function$$\begin{array}{c}\hfill f(s^1,\dots ,s^{n_d}):\Omega \subset {\mathbb{R}}^{n_d}\to \mathbb{R},\end{array}$$which depends on $n_d$ different variables, considered as dimensions of the state space $s^k$, $k=1,\dots ,n_d$.

The sparse PGD (sPGD) technique is based on approximating this unknown function f employing a separated (tensor) representation. As in standard PGD procedures, the function f is decomposed using a sum of products of one-dimensional functions each one involving one independent variable. Each sum is called a mode.

In the context of non-intrusive ROMs, the main objective is to find a function $\tilde{f}$ which minimizes the distance to the sought function

9

$$\begin{array}{c}\hfill {\displaystyle \tilde{f}={\arg \min }_{f^{\ast }}\sum _{\begin{array}{c}i=1\end{array}}^{n_t}∥f,({\overrightarrow{s}}_i),-,f^{\ast },({\overrightarrow{s}}_i)∥,}\end{array}$$and that takes the separated form

10

$$\begin{array}{c}\hfill {\displaystyle \tilde{f}(s^1,\dots ,s^{n_d})=\sum _{\begin{array}{c}m=1\end{array}}^M\prod _{k=1}^{n_d}{\psi }_m^k(s^k),}\end{array}$$where M is the number of modes and ${\psi }_m^k$ is the one-dimensional function for the dimension k and mode m. $n_t$ is the cardinality of the training set used to construct the model and ${\overrightarrow{s}}_i$ are the different vectors which contain the data points to perform the regression. $∥·∥$ is the chosen norm to measure the distance between two points.

The other goal is that the model $\tilde{f}$ has to perform as well in the training set as in other unseen scenarios. This second target is more difficult to reach, yet is more crucial because this indicates the predictive ability of the ROM $\tilde{f}$, that is, the prediction accuracy when it is fed with unseen data. Tis is specially challenging to achieve when facing a high-dimensional problem, which provides sparse data.

As previously introduced, the sPGD technique expresses the function $\tilde{f}$ with the separated form expressed by (10). Then, the functions ${\{{\psi }_m^k\}}_{m=1}^M$ for each k are formed by a linear combination of a set of basis functions:

11

$$\begin{array}{c}\hfill {\psi }_m^k(s^k)=\sum _{j=1}^D{N}_{j,m}^k(s^k)a_{j,m}^k={({\overrightarrow{N}}_m^k)}^{\top }{\overrightarrow{a}}_m^k,\end{array}$$where D represents the degrees of freedom of the selected approximation. Moreover, ${\overrightarrow{N}}_m^k$ is a column vector with the set of basis functions for the k dimension and the m-th mode and ${\overrightarrow{a}}_m^k$ is a column vector with the coefficients for the k dimension and the m-th mode. The important question here is which type of basis functions are the most convenient at hand. For instance, a polynomial basis or a Kriging basis can be chosen.

The computation of the coefficients in each one-dimensional function for each mode $m=1,\dots ,M$ is done by employing a greedy algorithm such that, once the approximation up to order $M-1$ is known, the new M-th order term is found using a non-linear solver (Picard, Newton, for instance):

12

$$\begin{array}{c}\hfill {\displaystyle {\tilde{f}}^M=\sum _{\begin{array}{c}m=1\end{array}}^{M-1}\prod _{k=1}^{n_d}{\psi }_m^k(s^k)+\prod _{k=1}^{n_d}{\psi }_M^k(s^k).}\end{array}$$A standard choice is to select the same basis for each one of the modes:

13

$$\begin{array}{c}\hfill {\overrightarrow{N}}_1^k={\overrightarrow{N}}_2^k=\dots ={\overrightarrow{N}}_M^k\text{, for}k=1,\dots ,n_d.\end{array}$$This choice may seem reasonable, however it may not be appropriate when dealing with non-structured sparse data.

It is known that the cardinality of the interpolation basis must not exceed the maximum rank provided by the training set. Indeed, this constraint, which provides an upper bound to build the interpolation basis, only guarantees that the minimization is satisfied by the training set, without saying anything of the other points. Hence, if there is not an abundance of sampling points in the training set, in the low-data limit, high oscillations may appear out of these measured points because of the increased risk of overfitting. Usually, this is an undesirable effect because it affects the predictive ability of the constructed regression model.

To deal with this problem, the sPGD uses the modal adaptivity strategy (MAS) to take advantage of the greedy PGD algorithm. The idea is to minimize spurious oscillations out of the training set by starting the PGD algorithm looking for modes with low order degree. When it is observed that the residual decreases slowly or stagnates, higher order approximation functions are introduced. By following this procedure, undesired oscillations are decreased, since a higher-order basis will try to capture only the remaining residual.

It is strongly recommended to define an indicator and a stopping criterion to design the above strategy. Many different procedures can be thought. Here, we employ that defined in references [12, 19], where the methodology of the sPGD is deeply described. Following that reference, the following norm is used for the PGD residual in the present work:

14

$$\begin{array}{c}\hfill {\displaystyle {\mathcal{R}}_T^M=\frac{1}{\sqrt{n_t}}\sqrt{\sum _{\begin{array}{c}i=1\end{array}}^{n_t}(f({\overrightarrow{s}}_i)-{\tilde{f}}^M({\overrightarrow{s}}_i){)}^2},}\end{array}$$where ${\mathcal{R}}_T^M$ is the residual of the PGD solution of M modes in the training set T and ${\tilde{f}}^M$ is the PGD solution composed of M modes.

Therefore, we define for each enrichment step, ${\tilde{f}}^M$ as:

15

$$\begin{array}{c}\hfill \Delta {\mathcal{R}}_T^M={\mathcal{R}}_T^M-{\mathcal{R}}_T^{M-1}<{ϵ}_r,\end{array}$$where ${ϵ}_r$ is a parameter defining the resilience of the sPGD to increase the number of elements of the interpolation basis.

Constructing a novel ROM by combining sPGD and RB

Regarding the compactness, robustness and simplicity of the PGD models obtained by the sPGD technique, global polynomial basis are usually selected to use the sPGD explained in Sect. 3.2.

As it is well known, polynomials are differentiable functions for all arguments. Therefore, trying to capture a differentiable function using polynomials is, above all, a consistent idea because the function which wants to be captured has the same properties that the basis where is projected.

In fact, according to the Weierstrass approximation theorem, if f is a continuous real-valued function on [a, b], for a given $ϵ$, then there exists a polynomial p on [a, b] such that:

16

$$\begin{array}{c}\hfill \mid f(x)-p(x)\mid <ϵ,\end{array}$$for all $x\in [a,b]$. In words, any continuous function on a closed and bounded interval can be uniformly approximated on that interval by polynomials to any degree of accuracy. However, sometimes the function which is trying to be captured is not a differentiable function in some points or even presents discontinuities.

In this case, global polynomials are far away to be the best choice to approach this type of functions because they are, in fact, very poor at interpolating discontinuities. To demonstrate this, Fig. 3 shows the interpolation of the unit step function with 16 points using a 15th-degree polynomial.

Fig. 3 [Images not available. See PDF.]

Approximating a unitary step function using global polynomials. In the first case, Chebyshev nodes are used to fit the polynomial. In the second case, equidistant nodes are used

Two different solutions to the problem can be envisaged. The first one is to use piecewise polynomial interpolation. This way, we can use different polynomials to approach the function at the right and at the left of the discontinuity. For example, the step function of the previous example will be composed of:

17

$$\begin{array}{c}\hfill f(x)=\left\{\begin{array}{cc}0& x\in [-1,0],\\ 1& x\in (0,1],\end{array}\right)\end{array}$$which are zero degree polynomials. However, an issue of this type of approaches is its rank deficiency when combined with the sparsity used in the sPGD [12].

Other possible solution is to use a basis (which can contain discontinuous functions) whose linear combination produces the class of discontinuities of our problem in the right places.

For instance, a basis composed of different step functions can be used to approximate the unitary step function discussed during this example.

In the industrial problem that this work is dealing with, discontinuities change their place in space when changing some values in the parameter space. Therefore, the previous discussed basis approach is prefered to deal with the discontinuity problem.

Other issue is that if we try to capture a non-regular function of this type without having preliminary knowledge of the system, a lot of nodes are needed to detect where and how these singular points are present in the dimension where this behaviour happens.

Therefore, to deal with the previous discussed issues, we propose the following approach:

1. Find the spatial dimension(s) where singular points are placed.

2. Detect parameter(s) which can change location of the singular points along spatial dimension/s.

3. Construct a RB considering the non-regular dimensions found in steps one and two. Not sparse sampling will be used along the dimension/s contained in the RB.

To insert the RB in the PGD procedure, we propose to reformulate the regression problem in the following way (where without loss of generality only one dimension—$s^1$—is assumed causing troubles):

18

$$\begin{array}{cc}\hfill \tilde{f}(s^1,\dots ,s^{n_d})=& {\beta }^1(s^2,\cdots ,s^{n_d})·N_1(s^1)\hfill \\ \hfill & \cdots +{\beta }^T(s^2,\cdots ,s^{n_d})·N_T(s^1),\hfill \end{array}$$where:

19

$$\begin{array}{c}\hfill {\displaystyle {\beta }^p(s^2,\cdots ,s^{n_d})=\sum _{\begin{array}{c}i=1\end{array}}^{I_p}\prod _{k=2}^{n_d}{\psi }_{i,p}^k(s^k),p=1,2,\cdots ,T.}\end{array}$$$N_1(s^1),\cdots ,N_T(s^1)$ form the RB obtained with the SVD (see “Appendix”) along the $s^1$ dimension, the ${\beta }^p(s^2,\dots ,s^{n_d})$ terms represent the unknown functions for the sPGD problem for a given $N_p$, $I_p$ is the number of modes used to decompose ${\beta }^p$ and ${\displaystyle {\sum }_{p=1}^T{I}_p}$ is the total number of modes of the ROM.

The training set is then used to obtain the reduced basis as well as the value of the ${\beta }^p$ coefficients in the training points. In addition, once these points are obtained, the sPGD procedure is used to obtain the separated representation of these functions using polynomial basis according to Eq. (19). To obtain the RB, we consider the SVD, revisited in “Appendix”, according to the discussion that follows.

Defining ${\overrightarrow{z}}_i$ as a point in the dimensions $(s^2,\dots ,s^{n_d})$ then, the set $\mathcal{Y}$ of one-dimensional functions, created by the points belonging to the training set $\mathcal{T}$, can be defined as:

20

$$\begin{array}{c}\hfill \mathcal{Y}=\{f(s^1,{\overrightarrow{z}}_i):{\overrightarrow{z}}_i\in \mathcal{T}\}\end{array}$$Therefore, the set $\mathcal{Y}$ is created collecting the one-dimensional functions for the different points which are selected for the training set $\mathcal{T}$ to do the regression via sPGD.

Consequently, the snapshots in $\mathcal{Y}$ are the ones used to construct the matrix $\mathbf{Y}$ (see “Appendix”). Then, the SVD can be used to extract a reduced basis ${\{{\overrightarrow{N}}_i\}}_{i=1}^T$, which best approximates the set $\mathcal{Y}$.

To end up, an important issue is the choice of arguments of the parametric model. They must be independent or poorly correlated to avoid increasing the redundancy and the complexity of the model without necessity.

Industrial and technical knowledge can be used to determine the appropriate choice of variables as it is the case in this work.

If not, a manifold learning (such as kernel-PCA, [22]) or a dimensionality reduction technique (viz. topological data analysis, [23]) can be applied to reduce unnecessary variables. In addition, the ANOVA analysis can also be carried out to determine the importance of each input to keep the most relevant ones.

Results

In this section, the results of the proposed ROMs for each motor are presented.

In addition, to illustrate some of the advantages concerning the separated representation of the PGD, a sensitivity analysis is carried in Sect. 4.3 to measure the impact of each variable.

Furthermore, a study comparing the error between the ROM and the FEM software FEMM is carried out to check the accuracy of the proposed approach.

Induction motor

Linear B–H

The approach shown is based on a current-based model where a balanced three-phase system is supposed for both the fundamental and the harmonic component of the current. The searched PGD function is:

21

$$\begin{array}{c}\hfill B_n(\alpha ,f,s,I_p,f_h,s_h,I_{\mathit{ph}},\gamma );B_n\in \mathbb{C},\end{array}$$where $\alpha$ refers to an angle pointing to a node in the air-gap midline, f is the supply frequency, s is the slip, $I_p$ is the current peak value of the fundamental frequency of the source, $f_h$ is the harmonic frequency of the source, $s_h=\frac{f_h-(1-s)f}{f_h}$ is a redefined slip concerning the harmonic component, $I_{\mathit{ph}}$ is the current peak value of the harmonic component, and $\gamma$ is the relative rotor position in relation to stator.

In this section, the problematic dimensions discussed in Sect. 3.3 are $\gamma$ and $\alpha$. Therefore, the RB procedure will be applied in these dimensions as explained in the reformulation of Sect. 3.3.

To obtain the above function, a multi-PGD procedure is used to decompose the function in more than one PGD solution. Consequently, the searched function is now:

22

$$\begin{array}{cc}\hfill B_n(\alpha ,f,s,I_p,f_h,s_h,I_{\mathit{ph}},\gamma )=& \widetilde{f_1}(\alpha ,f,s,I_p,\gamma )\hfill \\ \hfill & +\widetilde{f_2}(\alpha ,f_h,s_h,I_{\mathit{ph}},\gamma ),\hfill \end{array}$$where the sPGD technique will be used first for $\widetilde{f_1}$ and then for $\widetilde{f_2}$.

The first PGD search will focus on the range of parameters of the fundamental component of the source and the second one regarding the harmonic one.

Considering the extraction of snapshots the following remark must be considered to approach the problem.

As the system is considered linear (Linear B–H relationship), the total response is the sum of the responses obtained from each source considered separately (superposition theorem). Therefore, the chosen approach is to analyse harmonic content separately, considering each source component independently and then adding each time response to $B_n$. Finally, when the total $B_n(t)$ is obtained, the post-processing of Eq. (7) must be carried out.

To compare the results for different ${\overrightarrow{z}}_i$ along the $\alpha$ dimension, the following expressions are used:

23

$$\begin{array}{c}\hfill \begin{array}{cc}\hfill {\text{err}}_i^{\mathrm{real}}& =\frac{{∥\text{real},(B_n^{\mathrm{PGD}}(\alpha ,{\overrightarrow{z}}_i)),-,\text{real},(B_n^{\exp }(\alpha ,{\overrightarrow{z}}_i))∥}_2}{{∥\text{real},(,B_n^{\exp },(\alpha ,{\overrightarrow{z}}_i),)∥}_2},\hfill \\ \hfill {\text{err}}_i^{\mathrm{imag}}& =\frac{{∥\text{imag},(B_n^{\mathrm{PGD}}(\alpha ,{\overrightarrow{z}}_i)),-,\text{imag},(B_n^{\exp }(\alpha ,{\overrightarrow{z}}_i))∥}_2}{{∥\text{imag},(,B_n^{\exp },(\alpha ,{\overrightarrow{z}}_i),)∥}_2}\hfill \end{array}\end{array}$$where the superscript $\text{PGD}$ denotes the results obtained by the sPGD and the superscript “$\text{exp}$” denotes the experimental measurements (pseudo-experimental results obtained in this case by the FEMM software).

Furthermore, to sample the training set, different latin hypercubes (LHs) are taken using a grid composed of Chebyshev nodes along dimensions $f,s,I_p,f_h,s_h,I_{\mathit{ph}}$. In addition, eight hours are used to obtain the training set in the offline stage.

In Fig. 4, a comparison between the expperimental and PGD results for a ${\overrightarrow{z}}_i\notin \mathcal{T}$. For this plot, the error measured as Eq. (23) is ${\text{err}}_i^{\mathrm{real}}=0.000355$, ${\text{err}}_i^{\mathrm{imag}}=0.000341$; and ${\text{err}}_i^{\mathrm{real}}=0.000037$, ${\text{err}}_i^{\mathrm{imag}}=0.000029$ for the induction component caused by the harmonic. On the other hand, the error in the training set measured as Eq. (23) is always lower than ${10}^{-5}$.

In addition, the alternating force (specifically, the alternating magnetic pressure) for the selection of parameters used in Fig. 4 is shown in Fig. 5. This force is obtained combining the different sine/cosine waves of the B-field in Eq. (7).

Fig. 4 [Images not available. See PDF.]

Comparison between the sPGD and FEMM model (linear B–H relationship) for the parameters $f=40$ Hz, $s=2.5$ %, $I_p=3$ A, $f_h=5\times 40$ Hz, $s_h=80.50$ %, $I_{\mathit{ph}}=3/5$ A, $\gamma =1$ degree

Fig. 5 [Images not available. See PDF.]

Alternating force components obtained with the parameters $f=40$ Hz, $s=2.5$ %, $I_p=3$ A, $f_h=5\times 40$ Hz, $s_h=80.50$ %, $I_{\mathit{ph}}=3/5$ A, $\gamma =1$ degre (linear B–H model)

Fig. 6 [Images not available. See PDF.]

Relative error of the sPGD model (linear B–H relationship) to determine the real part of $B_n$ for different untrained ${\overrightarrow{z}}_i$. Error criteria of Eq. (23) is used

Fig. 7 [Images not available. See PDF.]

Relative error of the sPGD model (linear B–H relationship) to determine the imaginary part of $B_n$ for different untrained ${\overrightarrow{z}}_i$. Error criteria of Eq. (23) is used

In Figs. 6 and 7, the error for Eq. (21) for some ${\overrightarrow{z}}_i\notin \mathcal{T}$ can be seen. As it can be noticed, an excellent agreement between FEMM and PGD results is achieved even outside the training set. The main advantage is that the PGD model computes induction and force for a given ${\overrightarrow{z}}_i$ in less than 0.2 s independently of the computational cost of the finite element solutions used for the snapshots.

Nonlinear B–H

The approach shown here is also based on a current based model. In addition, perfect sine wave current functions are supposed for the balanced three-phase system. In this case, the searched PGD function is:

24

$$\begin{array}{c}\hfill B_n(\alpha ,f,s,I_p,\sigma ,\delta ,\gamma );B_n\in \mathbb{C},\end{array}$$where $\sigma$ is the bar conductivity and $\delta$ is the air gap of the electric machine.

In this section, the problematic dimensions described in Sect. 3.3 are $\alpha$ and $\gamma$. Therefore, the RB procedure will be applied in this dimension as explained in the reformulation of Sect. 3.3.

Furthermore, to sample the training set, different Latin Hypercubes (LHs) are taken using a grid composed of Chebyshev nodes along dimensions $f,s,I_p,\sigma ,\delta$. In addition, nine hours are used to obtain the training set in the offline stage.

In Fig. 8, a comparison for Eq. (24) is shown for a ${\overrightarrow{z}}_i$ belonging to the training set. For this plot, the error measured as Eq. (23) is ${\text{err}}_i^{\mathrm{real}}=0.0008$ and ${\text{err}}_i^{\mathrm{imag}}=0.001$.

In Fig. 9, a comparison between the FEMM and PGD results for a ${\overrightarrow{z}}_i\notin \mathcal{T}$. For this plot, the error measured as Eq. (23) is ${\text{err}}_i^{\mathrm{real}}=0.009$ and ${\text{err}}_i^{\mathrm{imag}}=0.012$. In addition, the alternating magnetic pressure related to this B-field can be seen.

In Figs. 10 and 11, the error for Eq. (24) for some ${\overrightarrow{z}}_i\notin \mathcal{T}$ can be observed.

As it can seen in these figures, excellent agreement between FEMM and PGD results is achieved outside the training set.

Fig. 8 [Images not available. See PDF.]

Comparison between the sPGD and FEMM model (nonlinear B–H relationship) for a ${\overrightarrow{z}}_i$ belonging to the training set. Parameters $f=10.4772$ Hz, $s=5$ %, $I_p=3.3858$ A, $\sigma =47.6537$ MS/m, $\delta =0.4876$ mm, $\gamma =0$ degrees. The alternating force obtained during the post-processing step is also shown

Fig. 9 [Images not available. See PDF.]

Comparison between the sPGD and FEMM model (nonlinear B–H relationship) for a ${\overrightarrow{z}}_i\notin \mathcal{T}$. Parameters $f=50$ Hz, $s=3$ %, $I_p=3.5$ A, $\sigma =40$ MS/m, $\delta =0.5$ mm, $\gamma =0$ degrees. The alternating force obtained during the post-processing step is also shown

Fig. 10 [Images not available. See PDF.]

Relative error of the sPGD model (nonlinear B–H relationship) to determine the real part of $B_n$ for different untrained ${\overrightarrow{z}}_i$. Error criteria of Eq. (23) is used

Fig. 11 [Images not available. See PDF.]

Relative error of the sPGD model (nonlinear B–H relationship) to determine the imaginary part of $B_n$ for different untrained ${\overrightarrow{z}}_i$. Error criteria of Eq. (23) is used

The main advantage is that the PGD model computes induction and force for a given ${\overrightarrow{z}}_i$ in less than 0.2 s indepentently of the computational cost of the Finite Element solutions used for the snapshots.

It is important to highlight that the computational cost of the PGD model is independent of the one of the FEM software used for the snapshots. Hence, if the computational cost of the FEM software was some days, the time needed for the PGD still would be less than 0.2 s.

Synchronous motor

The approach shown is based on a current based model. In addition, perfect sine wave current functions are supposed for the balanced three-phase system. In addition, a nonlinear B–H relationship is used for the materials.

The searched PGD function is:

25

$$\begin{array}{c}\hfill B_n(\alpha ,\gamma ,I_p,\tau ,\rho );B_n\in \mathbb{R},\end{array}$$where $\alpha$ refers to an angle pointing to a node in the air-gap midline, $\gamma$ is the relative rotor position in relation to stator, $I_p$ is the current peak value of the fundamental frequency of the source, $\tau$ is the torque angle, namely, the phase difference between rotor and stator magnetic fields and $\rho$ is a parameter defining the dynamic eccentricity of the machine ($\rho =\frac{D_{O_s,O_r}}{\delta }$; where $O_s$ is the stator symmetry center which in this eccentricity is equals to the rotor rotation center $O_w$, $O_r$ is the rotor symmetrical axis, $D_{O_s,O_r}$ is the distance between $O_s$ and $O_r$ and $\delta$ is the uniform air-gap length when there is no eccentricity). The above eccentricity parameters can be seen in Fig. 12.

Fig. 12 [Images not available. See PDF.]

Sketch to visualize the parameters $O_s,O_w,O_r$ and $D_{O_s,O_r}$

In this section, the problematic dimensions described in Sect. 3.3 are $\alpha$ and $\gamma$. Therefore, the RB procedure will be applied in these dimensions as explained in the reformulation of Sect. 3.3.

Furthermore, to sample the training set, different Latin Hypercubes (LHs) are taken using a grid composed of Chebyshev nodes with the exception of $\gamma$ (where it is preferred to use equidistant nodes for each parameter combination to complete a 180-degree turn). In addition, eight hours are used to obtain the training set in the offline stage.

In Figs. 13 and 14, the induction and the magnetic pressure are depicted for two points ${\overrightarrow{z}}_1,{\overrightarrow{z}}_2\notin \mathcal{T}$. The error associated according to Eq. (23) is respectively ${\text{err}}_1^{\mathrm{real}}=0.01$ and ${\text{err}}_2^{\mathrm{real}}=0.009$.

In Fig. 15, we can see the time evolution of the magnetic pressure of figure 14 supposing a rotor and synchronous speed of $N=2000$ rpm in steady state.

In Fig. 16, a comparison of the results between FEMM and the PGD are shown to analyze the error of the proposed ROM.

As it can be seen, the error analysis for this PGD solution has similar results to the other PGD functions obtained for the induction motor.

Finally, as in the other cases, the big advantage is that the PGD model computes induction and force for a given ${\overrightarrow{z}}_i$ in less than 0.2 s and with accuracy indepentently of the computational cost of the Finite Element solutions used for the snapshots.

Fig. 13 [Images not available. See PDF.]

Induction obtained with the sPGD model for the parameters: $\gamma =0$ degrees, $I_p=3$ A, $\tau =\pi /2$ rad, $\rho =0$

Fig. 14 [Images not available. See PDF.]

Radial force wave obtained with the sPGD model for a set of parameters: $\gamma =0$ degrees, $I_p=13$ A, $\tau =\pi /2$ rad, $\rho =0.5$

Fig. 15 [Images not available. See PDF.]

Time evolution of the magnetic pressure wave obtained with the sPGD supposing $N=2000$ rpm. Set of parameters: $\gamma =0$ degrees, $I_p=13$ A, $\tau =\pi /2$ rad, $\rho =0.5$

Fig. 16 [Images not available. See PDF.]

Comparison between sPGD and FEMM results. Set of the parameters (top plot): $\gamma =0$ degrees, $I_p=3$ A, $\tau =\pi /2$ rad, $\rho =0$. Set of the parameters (bottom plot): $\gamma =0$ degrees, $I_p=13$ A, $\tau =\pi /2$ rad, $\rho =0.5$

Sensitivity analysis of the parametric solutions

Sensitivity analysis is often interesting in parametric models because of:

• the need to characterize how sensitive the response is with respect to uncertainties in the input data; for example, manufacturing tolerances or material properties.

• the need to characterize how sensitive the response is in function of the operation point of the system.

• the need to make changes to improve the performance of a design and want to find out which changes that are most efficient for attaining the expected goals.

As a result, it seems evident that the importance of doing this type of analysis arises in industrial applications like the one treated in this paper. The most appealing point is that the ROMs constructed with the separated representation proposed in this work make easier performing this analysis. The reason is that computing partial derivatives in the separated representation is translated to compute the derivatives of one-dimensional functions. As the vast part of the one-dimensional functions used in the ROMs are polynomials, the computation remains quite simple and with low-computational cost enabling online real-time calculations.

To illustrate the procedure, imagine a PGD solution concerning 3 parameters:$$\begin{array}{c}\hfill {\displaystyle {\tilde{f}}^M=\sum _{\begin{array}{c}m=1\end{array}}^M{X}_m(x){\Theta }_m^1({\theta }_1){\Theta }_m^2({\theta }_2){\Theta }_m^3({\theta }_3),}\end{array}$$where we want to analyze the how sensitive is the solution around the nominal point ${\theta }_1=a_1,{\theta }_2=a_2,{\theta }_3=a_3$. Then, the partial derivatives $\frac{\partial {\tilde{f}}^M}{\partial {\theta }_1}(x,a_1,a_2,a_3)$, $\frac{\partial {\tilde{f}}^M}{\partial {\theta }_2}(x,a_1,a_2,a_3)$ and $\frac{\partial {\tilde{f}}^M}{\partial {\theta }_3}(x,a_1,a_2,a_3)$ are computed as:$$\begin{array}{cc}\hfill \frac{\partial {\tilde{f}}^M}{\partial {\theta }_1}(x,a_1,a_2,a_3)=& {\displaystyle \sum _{\begin{array}{c}m=1\end{array}}^M{X}_m(x)\frac{\partial {\Theta }_m^1}{\partial {\theta }_1}(a_1){\Theta }_m^2(a_2){\Theta }_m^3(a_3),}\hfill \\ \hfill \frac{\partial {\tilde{f}}^M}{\partial {\theta }_2}(x,a_1,a_2,a_3)=& {\displaystyle \sum _{\begin{array}{c}m=1\end{array}}^M{X}_m(x){\Theta }_m^1(a_1)\frac{\partial {\Theta }_m^2}{\partial {\theta }_2}(a_2){\Theta }_m^3(a_3),}\hfill \\ \hfill \frac{\partial {\tilde{f}}^M}{\partial {\theta }_3}(x,a_1,a_2,a_3)=& {\displaystyle \sum _{\begin{array}{c}m=1\end{array}}^M{X}_m(x){\Theta }_m^1(a_1){\Theta }_m^2(a_2)\frac{\partial {\Theta }_m^3}{\partial {\theta }_3}(a_3),}\hfill \end{array}$$where the magnitude of the partial derivative indicates the sensitivity of the solution for the given parameter in the analysed point (higher magnitude corresponds to higher sensitivity). In addition, the sign of the partial derivative indicates the direct or indirect relationship between the function value and the parameter.

Higher-order derivatives and also derivatives concerning different parameters are also possible. It is important to note that the sensitivities are a field, offering large information, beyond the behavior around the operating point.

Since the ROMs are based in the B-field, the sensitivity relation between magnetic pressure and the B-field for a given parameter ${\theta }_i$ is [using Eq. (7)]:

26

$$\begin{array}{c}\hfill \frac{\partial P_n}{\partial {\theta }_i}(\alpha )=\frac{1}{\mu }·B_n\frac{\partial B_n}{\partial {\theta }_i}(\alpha ).\end{array}$$Now, as an example of use, we show some analysis that can be carried out with the extracted ROMs in different regions of the domain. In the first example, we are going to focus on the synchronous machine. Here, the sensitivity of the solution for the parameters $I_p$ and $\rho$ is explored in an area of interest $\alpha \in [134.7,224.8]$ deg. under the conditions: high current and no eccentricity, high current and eccentricity, low current and no eccentricity, and low current and eccentricity.

In Figs. 17, 18, 19 and 20 the sensitivity is compared with two operating points changing its $I_p$ from 0.5 A to 20 A. Here, we can observe how much a little change in the eccentricity or in the current peak value affects the magnetic pressure as well as the B-field for both cases. In this operating points, adding eccentricities in the order of 0–0.25 mm does not change a lot the sensitivity behavior.

Fig. 17 [Images not available. See PDF.]

Sensitivity of $B_n$ and $P_n$ for small changes in $I_p$. Operating point: $\gamma =0$ degrees, $I_p=0.5$ A, $\tau =\pi /2$ rad, $\rho =0$

Fig. 18 [Images not available. See PDF.]

Sensitivity of $B_n$ and $P_n$ for small changes in $D_{O_s,O_r}$. Operating point: $\gamma =0$ degrees, $I_p=0.5$ A, $\tau =\pi /2$ rad, $\rho =0$

Fig. 19 [Images not available. See PDF.]

Sensitivity of $B_n$ and $P_n$ for small changes in $I_p$. Operating point: $\gamma =0$ degrees, $I_p=20$ A, $\tau =\pi /2$ rad, $\rho =0$

Fig. 20 [Images not available. See PDF.]

Sensitivity of $B_n$ and $P_n$ for small changes in $D_{O_s,O_r}$. Operating point: $\gamma =0$ degrees, $I_p=20$ A, $\tau =\pi /2$ rad, $\rho =0$

The second example is reported in Fig. 21 concerning the induction motor. Here, the sensitivity of the solution for five parameters is explored under a nominal operation point. Also, the sensitivity of the conductivity when its value is not well known was studied concluding that this uncertainty does not affect strongly the variables of interest.

Fig. 21 [Images not available. See PDF.]

Sensitivity of $P_n$ (the modulus of the alternating component) for small changes in the parameters shown in the plot. The operating point is $f=50$ Hz, $\gamma =0$ degrees, $I_p=3$ A, $s=2.5$ %, $\sigma =34.45$ MS/m, $\delta =0.375$ mm. The region explored is $\alpha \in [22.2,53.8]$ deg.

Conclusions

In this paper a ROM is developed combining both the sPGD and RB techniques. It can be observed that the results of the FEMM software are reproduced with a high accuracy using the proposed ROM model. A reduction in the computational time and resources needed to obtain parametric electromagnetic forces is achieved. In fact, the computational cost can be carried out by a standard laptop in less than 0.2 s.

The saving in computational time and resources opens a door for design, analysis, optimization and simulation of NVH in electric motors under this rationale in the electromagnetic step.

In addition, the proposed ROM facilitates the integration and coupling of the force computation in electric motors to other systems (such as the EV system) because of the simplicity of the obtained algebraic expression.

The extremely low computational cost of the proposed ROM is independent of the complexity of the model used to offline obtain the snapshots. In this paper the two-dimensional models presented in Sect. 2 are used through formulation available in the free software FEMM.

Furthermore, richer finite element models can be used to obtain the electromagnetic forces to construct the ROM without major difficulties. For example, a transient three-dimensional model taking into consideration motion.

Although it is true that using these models to obtain the snapshots for the offline stage is more time consuming, the computing time needed for the ROM once it is constructed still would not be affected (as it discussed in the results, the computing time of the ROM is less than half a second).

Our works in progress address more complex models such as a transient 3-dimensional model as well as the coupling of the numerical model with circuit equations. Thus, the computational cost reduction would be more drastic allowing to deal with the NVH problem with high-accurate models which can be computationally prohibitive.

In addition, while analysing these cases, it would be interesting to add other type of faults as a parameters in the sPGD model to see its effect under different circumstances.

Acknowledgements

We thank Professor Jose Roger-Folch (Universitat Politècnica de València) for useful discussion, comments and feedback about the present work.

Availability of data and material

Not applicable.

Code availability

Not applicable.

Compliance with ethical standards

Conflict of interest

The authors declare that they have no conflict of interest.

Publisher's Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

References

1. (2015) Paris declaration on electro-mobility and climate change & call to action. https://unfccc.int/media/521376/paris-electro-mobility-declaration.pdf. Accessed 24 Oct 2019

2. (2019) European energy research alliance - eu projects - batteries europe. https://www.eera-set.eu/eu-projects/batteries-europe/. Accessed 24 Oct 2019

3. Bai G (2018) Numerical simulation and optimization for electromagnetic noises of permanent magnet synchronous motors in vehicles. J Vibroeng. https://doi.org/10.21595/jve.2017.18720

4. Bianchi, N. Electrical machine analysis using finite elements; 2005; Boca Raton, CRC Press Taylor & Francis Group: [DOI: https://dx.doi.org/10.1201/9781315219295]

5. Chinesta, F; Keunings, R; Leygue, A. The proper generalized decomposition for advanced numerical simulations: a primer; 2013; New York, Springer:zbMath ID: 1287.65001

6. Cueto, E; Gonzalez, D; Alfaro, I. Proper generalized decompositions: an introduction to computer implementation with Matlab; 2016; 1 New York, Springer: [DOI: https://dx.doi.org/10.1007/978-3-319-29994-5]

7. Dupont JB, Bouvet P (2012) Multiphysics modelling to simulate the noise of an automotive electric motor. In: 7th international styrian noise, vibration & harshness congress: the European automotive noise conference. SAE International. https://doi.org/10.4271/2012-01-1520

8. Dupont JB, Bouvet P (2013) Noise radiated by an electrical powertrain: multiphysical simulation. In: CFM2013. CFM 2013

9. Dupont, JB; Saucy, H. Noise radiated by electric motors: simulation process and overview of the optimization approaches; 2019; Wiesbaden, Springer: pp. 107-121. [DOI: https://dx.doi.org/10.1007/978-3-658-20251-4_7]

10. Gieras, JF; Wang, CR; Lai, JCS. Noise of polyphase electric motors; 2006; Boca Raton, CRC Press Taylor & Francis Group: [DOI: https://dx.doi.org/10.1201/9781420027730]

11. Haas, S; Ellermann, K. Development and analysis of radial force waves in electrical rotating machines. Technische Mechanik; 2017; 37, 2–5 pp. 218-225. 10.24352/UB.OVGU-2017-098

12. Ibáñez Pinillo, R; Abisset-Chavanne, E; Ammar, A; González, D; Cueto, E; Huerta, A; Louis Duval, J; Chinesta, F. A multidimensional data-driven sparse identification technique: the sparse proper generalized decomposition. Complexity; 2018; 2018, pp. 1-11. [DOI: https://dx.doi.org/10.1155/2018/5608286]zbMath ID: 1407.93063

13. Jack, AG; Mecrow, BC. Methods for magnetically nonlinear problems involving significant hysteresis and eddy currents. IEEE Trans Magn; 1990; 26, 2 pp. 424-429. [DOI: https://dx.doi.org/10.1109/20.106344]

14. Kumar D, Kottalgi S, Sambharam T, Mandloi P (2017) A multiphysics optimization approach to design low noise and light weight electric powertrain noise, vibration and harshness (nvh) prediction of electric powertrain using finite element analysis (fea) and optimization. In: IECON 2017: 43rd Annual Conference of the IEEE Industrial Electronics Society, pp 1692–1697. https://doi.org/10.1109/IECON.2017.8216287

15. Mazgaonkar N, Chowdhury M, Fernandes LF (2019) Design of electric motor using coupled electromagnetic and structural analysis and optimization. In: SAE Technical Paper. SAE International. https://doi.org/10.4271/2019-01-0937

16. Meeker DC, Finite element method magnetics, version 4.2 (28feb2018 build) edn. http://www.femm.info

17. Pile, R; Devillers, E; Le Besnerais, J. Comparison of main magnetic force computation methods for noise and vibration assessment in electrical machines. IEEE Trans Magn; 2018; 54, 7 pp. 1-13. [DOI: https://dx.doi.org/10.1109/TMAG.2018.2828388]

18. Salon, SJ. Finite element analysis of electrical machines; 1995; New York, Springer: [DOI: https://dx.doi.org/10.1007/978-1-4615-2349-9]

19. Sancarlos, A; Cameron, M; Abel, A; Cueto Elias Duval, J; Chinesta, F. From rom of electrochemistry to ai-based battery digital and hybrid twin. Arch Comput Methods Eng; 2020; [DOI: https://dx.doi.org/10.1007/s11831-020-09404-6]

20. Sancarlos-González A, Pineda SM, Puche-Panadero R, Sapena-Bañó A, Riera-Guasp M, Martinez-Roman J, Perez-Cruz J, Roger-Folch J (2017) Application of the parametric proper generalized decomposition to the frequency-dependent calculation of the impedance of an ac line with rectangular conductors. Open Phys. https://doi.org/10.1515/phys-2017-0113

21. Sathyan S, Aydin U, Belahcen A (2020) Acoustic noise computation of electrical motors using the boundary element method. Energies 13(1). https://doi.org/10.3390/en13010245. https://www.mdpi.com/1996-1073/13/1/245

22. Scholkopf, B; Smola, A; Muller, K. Burges, CJC; Schölkopf, B; Smola, AJ. Kernel principal component analysis. Advances in kernel methods: support vector learning; 1999; Cambridge, MIT Press: pp. 327-352.

23. Wasserman, L. Topological data analysis. Ann Rev Stat Appl; 2018; 5, 1 pp. 501-532.MathSciNet ID: 3774757[DOI: https://dx.doi.org/10.1146/annurev-statistics-031017-100045]

24. Xu, X; Han, Q; Chu, F. Review of electromagnetic vibration in electrical machines. Energies; 2018; 11, 1779. [DOI: https://dx.doi.org/10.3390/en11071779]