Magnetic Equivalent Circuit and Lagrange Interpolation Function Modeling of Induction Machines Under Broken Bar Faults

This article introduces a mesh-based magnetic equivalent circuit (MEC) modeling technique for induction machines (IMs) in healthy and broken rotor bars conditions. The MEC model is presented as a highly accurate and computationally efficient alternative to finite element (FE) models. By incorporating modifications to the air gap coupling method, including a new Lagrange interpolation function, and utilizing a harmonic MEC model, the accuracy of the solution is improved while reducing electrical and mechanical transients. Compared to experiments and 2-D FE models, this model achieves precise results for electromagnetic torque, rotational speed, and forces across various conditions. The Lagrange interpolation function forms the basis for the air gap coupling between stator and rotor flux densities. The results demonstrate the MEC model’s exceptional accuracy in predicting speed oscillations, calculating forces, and analyzing current harmonics in faulty IMs. Furthermore, the MEC model performs over 30 times faster than the 2-D FE models.


I. INTRODUCTION
I NDUCTION machines (IMs) are extensively employed in various applications owing to their reliability.As a result, multiple modeling techniques have been developed over time.While finite element (FE) models offer high accuracy [1], they may not be the most efficient tool for quick results during initial machine design or real-time simulation integration.Consequently, analytical tools have been devised [1], [2].A comprehensive overview of these models, covering various types of faults in IMs, can be found in [1].
The magnetic equivalent circuit (MEC) model [3], [4] offers a valuable alternative in electric machine modeling.It strikes a balance between accuracy and computational efficiency by providing higher accuracy than analytical models while requiring fewer computational resources compared to FE models.However, the study presented in [4] indicates that while the computational times are low, the accuracy of torque ripple prediction falls short compared to 2-D FE models.Additionally, the forces generated in different fault scenarios were not included in any of the models discussed.These forces can be crucial for conducting further vibration analysis of electric machines, as highlighted in [5].
The Lagrange interpolation function described in [6] for the MEC model is used for coupling the stator and rotor nodes.This technique is applied to the nodal-based MEC for permanent magnet machines.Raising the interpolation order in the Lagrange function has been observed to result in instability in the solution and necessitates a greater number of iterations to attain the desired tolerance.This article presents modifications made to the air gap coupling method described in [6] for the mesh-based MEC.Instead of considering only the central angle in the Lagrange function, the approach now involves dividing each block into multiple angles and taking the average of the Lagrange function.This modification ensures a highly accurate solution.This approach is applied in this study to analyze various broken bar faults in IMs.The article provides a description of the air gap coupling methodology and the stator and rotor electric circuit equations.Furthermore, a harmonic MEC model is developed to initialize the time-stepping MEC and minimize mechanical and electrical transients.Non-linear iterations are handled using the Newton-Raphson technique.

II. MEC MODEL INTRODUCTION
The MEC models for a stator and a rotor tooth are depicted in Fig. 1 The model incorporates the non-linear behavior of the iron material using the Newton-Raphson technique.The non-linear elements are represented by the red-colored reluctances [7].In addition, the rotor skewing effect is modeled using the multi-slice 2-D modeling approach [7].To account for the impact of broken rotor bar conditions on mechanical transients, the model takes the mechanical torque as input and adjusts the rotational speed based on the electromagnetic output torque of the model and the system's mechanical inertia.

A. Air Gap Coupling
The coupling between the stator and rotor nodes is done using the Lagrange interpolation function [6].This interpolation function is a high-order polynomial of order n that couples the radial flux density variation in the air gap region and also the scalar potential in the stator V s1 -V s4 and the rotor V r1 -V r4 , shown in Fig. 1(a) and (b) respectively.To obtain a complete coupling between the stator and rotor, additional equations are included in the system matrix.In order to improve the accuracy of the interpolation technique for each MEC block, instead of using a single central angle for each stator and rotor block, k child angles for each block are used.This leads to a total of k 2 combinations between each angle.The average of all combinations is then obtained by dividing the sum of all combinations by k 2 .This approach results in a much more accurate interpolation technique, thus improving the overall performance of the MEC model.The stator's scalar potential at a primary stator block j, denoted by V s j , based on the rotor's scalar potential at a primary block w, denoted by V rw is expressed as where the stator angle at a primary block j and a child angle division u is denoted by θ s j u .In contrast, the rotor angle at a primary block b and a child angle division v is denoted by θ rb v , and the rotor angle at a primary block w and a child angle division v is denoted by θ rw v .The third summation and the product indicate the Lagrange interpolation order n for the adjacent main angles.Similarly, the rotor flux density at each block w denoted as B rw can be expressed as a function of the stator flux density at a primary block j, denoted by V s j .The total vector φ tot to be solved can be expressed as where φ s is the flux loop vector in the stator excluding the flux loops in the stator's air gap φ sg , shown in Fig. 1(a).Similarly, φ r is the loop flux vector in the rotor, excluding the ones in the rotor's air gap φ rg , shown in Fig. 1(b).V rg is the rotor scalar potential vector defined in (1) and V bo is a scalar potential vector to describe the symmetry conditions.In this case, the residual function r can be expressed as where L is a connection matrix of each branch flux and the loop flux [7].R totd is a diagonal matrix comprising of all reluctance elements.G is a matrix that relates the air gap stator loop fluxes φ sg and rotor loop fluxes φ rg with the rotor scalar potentials V rg and ensures the equality of the air gap flux densities.In addition, G includes the relation between the symmetry scalar potential vector V bo and the total flux loop vectors.M s and M r link each loop flux equation with the corresponding stator phase current vector I s and rotor loop current vector I r , respectively.I s and I r have sizes of (m × 1) and (N r × 1), respectively, where m and N r are the number of phases and the number of rotor bars, respectively.

B. Stator and Rotor Electric Circuit Equations
To include the stator electric circuit, additional equations are required to input the stator voltage at each time step.The model can consider either a star or a delta-connected machine but only a delta-connected one is presented here.For a delta-connected machine, the stator line voltages are applied directly to the terminals of the machine.They are organized into a vector V s of a size (m × 1).The phase resistances and end winding inductances are formed in diagonal matrices of R s and L s , respectively.The rotor electric circuit consists of N r loop current equations that form vector I N r .In each loop, the bar resistance and the end ring resistance are considered to form the rotor resistance matrix R r .Similarly, the end ring inductance in each loop is used to form the rotor inductance matrix L r .To consider the broken rotor bars, the conductivity of R b of a certain bar is controlled to present the fault.
The utilization of the Crank-Nicolson method transforms the stator and rotor time differential equations into finite difference form, which leads to the incorporation of a new residual function r where q s represents the number of symmetries considered.t is the time step considered between time instants k t and k t − 1.To minimize mechanical and electrical transients, the suggested approach is to employ a time-harmonic steadystate MEC model.In this model, the time derivatives in the stator electric circuit equations are converted to −iω s , where ω s represents the supply angular frequency.While the time derivatives in the rotor electric circuit equations are converted to −isω s , with s denoting the slip of the induction machine (IM).Furthermore, the slip is adjusted until the electromagnetic torque matches the desired input mechanical torque.The output flux loops and currents output from the time-harmonic model are the solution at the first time step for the time-stepping MEC model indicated in (4) for efficient transient reduction.

III. SIMULATION AND EXPERIMENTAL RESULTS
The experimental setup for data collection is illustrated in Fig. 2. In this setup, a four-pole, 50 Hz IM with a power rating of 7.5 kW is employed as the test machine.For every experimental trial, the rotor featuring a broken rotor bar is substituted.Various parameters including speed, torque, and three-phase currents are measured at a sampling rate of 20 kHz.
All 2-D FE and MEC simulations are conducted on a 11th Gen Intel1 Core2 i5-11600 with 32 GB installed RAM.Two radial slices are considered to model the skewing effect of the rotor in both the 2-D FE and MEC models.The simulations are conducted for 100 cycles, each cycle with 400-time steps.The 2-D multi-slice FE models are developed using COMSOL software.The number of triangles in the 2-D FE model is 38 338 with linear shape function elements for the two axial slices.The MEC model has 1380 mesh elements for the stator and 1244 for the rotor for the two axial slices.Various simulation scenarios are performed to test the resilience of the MEC model.These simulations include different conditions, such as a healthy state, one broken rotor bar (1BRB), two broken rotor bars (2BRBs), and three broken rotor bars (3BRBs).In order to examine the impact of child number divisions k in the Lagrange interpolation, diverse loading conditions are applied to both the experimental setup and the MEC models, specifically for values of k equal to 1 and 3.
The rotor speed of the MEC model, shown in Fig. 3, is compared with the experimental and 2-D FE model results for different broken rotor bars at rated loading conditions.The results show that the MEC model can obtain similar speed oscillations and average speed with great accuracy compared  to experiments and the 2-D FE model.The amplitudes of the harmonic spectrum of the current waveform Î a at the steady state of the MEC model, shown in Fig. 4, is compared with the experiments and the 2-D FE model for a different number of broken rotor bars at rated loading conditions.It is clear that the MEC and 2-D FE models can capture the lower sidebands (LSBs) and upper sidebands (USBs) amplitudes and frequencies of the faults with great accuracy.The horizontal forces F x are also compared in Fig. 5 for both models for different broken bars at rated loading conditions.The forces are computed using Maxwell's stress tensor [8].The results show the robustness of the MEC model.To assess the efficacy of the proposed model, varying loading levels are tested in the MEC model, considering different child division numbers, namely k equal 1 and 3.The investigation also includes one, two, and three broken rotor bars.The amplitudes of the LSB and USB currents denoted as Î LSB and Î USB , are observed in Fig. 6 for different broken rotor bar conditions.As the number of broken rotor bar increase, the MEC model with k equal to 3 exhibits higher accuracy compared to the case of k equal to 1 for detecting the LSB and USB amplitudes.
Table I shows the CPU time comparison between the 2-D FE model and the MEC model for both cases of child number discretizations of k equal to 1 and 3.The table shows that the MEC model is 30 times faster than the 2-D FE model without a significant difference between k equal to 1 and 3.

IV. CONCLUSION
The MEC model efficiently represents the electromagnetic properties of IMs, both in healthy and broken bar conditions.Compared to 2-D FE models, the MEC model offers high accuracy while requiring shorter computational times.It employs a novel Lagrange interpolation technique, enabling it to effectively capture harmonic sidebands and accurately compute electromagnetic forces during faults.Moreover, the MEC model can serve as valuable input for future vibration analysis purposes and can be extended for various types of electric machines.

Fig. 1 .
Fig. 1.MEC models for the stator and the rotor.The MMFs F s (d, 1−2), F s (d + 1, 1−2), F r (c, 1−3), and F r (c + 1, 1−3) are associated with the stator slots d and d + 1, as well as the rotor slots c and c + 1, respectively.The vector φ s represents the flux loop vector in the stator, excluding the flux loops present in the stator's air gap denoted as φ sg .Similarly, φ r corresponds to the loop flux vector in the rotor, excluding those in the rotor's air gap represented as φ rg .Moreover, V sg and V rg indicate the stator and rotor scalar potential vectors in the air gap region, respectively.(a) Stator MEC.(b) Rotor MEC.

Fig. 4 .
Fig. 4. Harmonic spectrum for the stator phase current for (a) 1BRB, (b) 2BRBs, and (c) 3BRBs for the experiments, 2-D FE, and MEC models at rated loading condition.