# Mathematical Model of a Three-Phase Induction Machine in a Naturalabc Reference Frame Utilizing the Method of Numerical Integration of Average Voltages at the Integration Step and Its Application to the Analysis of Electromechanical Systems

Authors: Oleksiy Kuznyetsov
Journal: Mathematical Problems in Engineering
 (2019)

---

## Abstract

Recent advances in the real-time simulation of electric machines are linked with the increase in the operation speed of the numerical models retaining the calculation accuracy. We propose utilizing the method of average voltages at the integration step (AVIS) for the design of a three-phase induction machine’s model in its naturalabc reference frame. The method allows for avoiding rotational e.m.f. calculation at every step; in turn, the electromagnetic energy conversion is accounted by the change of flux-linkage. The model is integrated into the object-oriented environment in C++ for designing the computer models of electromechanical systems. The design of the model of an electromechanical system utilizing the proposed approach is explained in an example. The behavior of the numerical models of a three-phase IM has been compared for the set of conventional numerical methods as well as first- and second-order AVIS. It has been demonstrated that both first- and second-order AVIS methods are suitable tools for high-speed applications, namely, AVIS provides higher maximum possible integration step (e.g., first-order AVIS provides 4 times higher than the second-order Runge–Kutta method, and the second-order AVIS provides 2.5 times higher than the first-order method). Therefore, we consider the most preferable order of the AVIS method for the high-speed applications is the second order, while the first order may be a suitable alternative to increase the calculation speed by 30% with the acceptable decrease in the accuracy.


---

## Body

## 1. Introduction

Plenty of approaches are proposed to develop the models of rotating electric machines; many of them are examined in a review [1]. It is obvious that the models based on a finite element method provide the deepest understanding of the phenomena inside the machine, although this deals with the increase in the simulation duration. The finite element method is sometimes the only possible tool to analyze internal faults in a machine [2]. However, most problems related to the design and analysis of electromechanical systems (EMS) require less accurate representation of the machine; e.g., in [3], the simplified IM model is exploited for power system dynamic analysis (balancing of active and reactive powers).We focus on the rotating machines’ models for the interaction of physical equipment with mathematical models in the real-time mode [4]. The real-time techniques develop together with the evolution of modeling methods and simulation techniques related to electric machines [5, 6] to provide their ability of real-time operation.Commercial simulation environments with the real-time-operation ability utilize the equations based on the nodal (or the modified nodal) or state variable equations [7]. Thus, [8] proposes a set of rotating machines’ models realized using the Simulink block using the state-space-nodal theory that allows their application in an eMEGAsim digital real-time stimulator with increased stability.Before selecting a suitable model, we should first distinguish the area of its application. Thus, for offline computation, the most critical parameter is the accuracy of the results, while the real-time operation requires both accuracy and the operation speed. Focusing on an induction machine (IM), two basic families of modeling approaches are commonly used. The first one represents an IM in its naturalabc reference frame, while the second utilizes orthogonal reference frames. In general, a traditional representation of an IM in a natural abc reference frame [9] (based on voltage balances for stator and rotor circuits) is not far from the first attempt in the 1930s [10]. The reference frames belonging to the second family are stationary αβγ representing Clarke transform and a rotating dq0 reference frame representing Park transform [9, 11]. In [9], both are examined in detail, but simply, the use of orthogonal reference frames reduces the number of equations and allows avoiding interdependences between axes.Nevertheless, as it was mentioned in [12], the standard dq0 IM model calculates the current after the fundamental harmonic and the average torque—only near the operating point. However, reducing the number of equations is critical for the high-speed applications; therefore, the representation of an IM in an orthogonal reference frame is widely utilized for the real-time applications, e.g., in [13], for the real-time emulation of induction motor loaded with a centrifugal compressor. The model of an IM in the voltage-behind-reactance formulation [14] combines both abc and dq0 representations, so as slower subsystem represents rotor dynamics and faster subsystem represents stator dynamics [15]. In [16], energy-based models of electric motors are utilized to simply take into account nonlinearity and nonsinusoidality in the machine. In [17], internal faults in the machine for the real-time simulation are considered utilizing a magnetic equivalent circuit.In the current work, we limit ourselves to the model of an IM in a naturalabc reference frame as far as it provides the ability of simply modeling plenty of possible machine’s modes of operation, thus remaining simple for the calculation, e.g., in comparison with finite element models. The present work focuses on the increase in the speed of processing a model of a three-phase IM in the abc reference frame utilizing the method of numerical integration of average voltages at the integration step (AVIS) [18]. The method is based on a reformulation of Kirchhoff’s equations to increase the simulation speed [18]. Thus, it has been successfully applied for the real-time emulation of a synchronous machine [19, 20]. In the current work, we summarize the results reported earlier by us in [21–23] as well as deepen the mathematical background of our previous works.


## 2. Representation of an IM in theabc-Reference Frame

### 2.1. Basic Equations of the IM

The analysis of an IM in a naturalabc reference frame is based upon the assumptions as follows:(i)
The sinusoidal magnetic flux distribution in the air gap of an IM(ii)
The constancy of the parameters of windings while the operation(iii)
Neglecting the effects of hysteresis and eddy currents (considering no iron losses)(iv)
Neglecting the mechanical losses (considering the equality of mechanical and electromechanical torques)(v)
Neglecting the skin effectThese assumptions are common for the representation of an IM in the form of an equivalent circuit. Therefore, in theabc reference frame, the electrical-balance equation for a stator phase of an IM is depicted as(1)uA=rAiA+dψAdt.Since the flux-linkage of a stator phase is a function of the stator (iA, iB, iC) and the rotor (ia, ib, ic) currents and of the angle between the corresponding axes of stator and rotor ϑ, the flux-linkage derivative in (1) is(2)dψAdt=∂ψA∂iA⋅diAdt+∂ψA∂iB⋅diBdt+∂ψA∂iC⋅diCdt+∂ψA∂ia⋅diadt+∂ψA∂ib⋅dibdt+∂ψA∂ic⋅dicdt+∂ψA∂ϑ⋅dϑdt.The partial derivatives with respect to currents in (2) are the dynamical inductances. Neglecting nonlinearities in a magnetic path of an IM, they are represented as static inductances: the self-inductance of a stator A-phase LAA, mutual inductances between the stator phases LAB and LAC, and mutual inductances between the stator and rotor LAa, LAb, and LAc, respectively. The item,(3)∂ψA∂ϑ⋅dϑdt=∂ψA∂ϑzpω,is the rotational e.m.f.The self-inductance of a statorA-phase LAA includes the components from the linkage magnetic flux and the leakage magnetic flux. Taking into account that the inductance from the main magnetic flux for a stator phase Lm and the inductance of a single phase Lms(4)Lm=32Lms,the similar ratio for the leakage inductance Lσs [9](5)LAA=23Lm+Lσs,where Lm and Lσs are the parameters of a traditional equivalent circuit of an IM. It is obvious that an IM is represented by an equivalent circuit considering the symmetrical machine and the symmetrical power supply. The parameters of the equivalent circuit are conveniently exploited in the modeling of an IM, and the methodology of determining those parameters is obvious.The equation for a mutual inductance between stator windings is calculated under the assumption of their symmetry:(6)LAB=LAC=LAAcos2π3=−13Lm+Lσs.The mutual inductances between the stator and the rotorLAa, LAb, and LAc are the functions of the angle ϑ, and the mutual inductance Lm between the stator and the rotor being calculated for ϑ=0 is(7)LAa=23kLmcosϑ,LAb=23kLmcosϑ+2π3,LAc=23kLmcosϑ−2π3,where the rotor winding is referred to the stator winding with the turn ratio k=Nr/Ns, where Ns and Nr are the numbers of turns in the stator and the rotor windings.The equation for a rotora phase (referred to a stator winding) is similar to (1):(8)ua=raia+dψadt,where the equation for the flux-linkage derivative dψa/dt is similar to (2), including a rotor phase’s self-inductance Laa, mutual inductances between the rotor windings Lab, Lac, and mutual inductances between the rotor and the stator LaA, LaB, and LaC.The self-inductance of a rotor phase is calculated similarly to (5) considering the turn ratio k as(9)Laa=23k2Lm+Lσr.The mutual inductance between the rotor windings considering their symmetry is(10)Lab=Lac=Laacos2π3=−13k2Lm+Lσr.The mutual inductances between the rotor and the stator are calculated similarly to (7).Other phases of an IM are depicted by the equations similar to (1), (2), and (5)–(10).


### 2.2. Representing the Equations of a Three-Phase IM in a Matrix Form

The equations of the electrical balance for stator and rotor circuits, (1), (2), and (5)–(10), are simply represented in a matrix form:(11)us=Rsis+ddtLsis+Lsrir+∂ψs∂ϑzpω,ur=Rrir+ddtLrsis+Lrir+∂ψr∂ϑzpω,where us=uAuBuCT, ur=uaubucT, is=iAiBiCT, ir=iaibicT, ψs=ψAψBψCT, and ψr=ψaψbψcT are the vectors of voltages, currents, and flux-linkages on the terminals of the stator and rotor winding; Rs=diagrA,rB,rC and Rr=diagra,rb,rc are the diagonal resistance matrices of stator and rotor windings. The elements of stator and rotor self-inductance matrices Ls and Lr and matrices of mutual inductances between the stator and the rotor Lsr and Lrs in (11) assuming the symmetrical IM are defined after (5)–(7) and (9)-(10).Thus, the matrices of mutual inductances between the stator and the rotor are(12)Lsr=LrsT=LAaLAbLAcLBaLBbLBcLCaLCbLCc,and self-inductance matrices Ls and Lr in the case of symmetrical IM are(13)Ls=23Lm+Lσs1−12−12−121−12−12−121,Lr=23k2Lm+Lσr1−12−12−121−12−12−121.Note that for a symmetrical IM, the matricesLs and Lr in (13) are singular. Thereafter, for a symmetrical IM, the calculation of the inverse matrix (related to the operation of solving the matrix equation) is impossible. We propose avoiding this issue as it is demonstrated in Appendix A.


## 2.1. Basic Equations of the IM

The analysis of an IM in a naturalabc reference frame is based upon the assumptions as follows:(i)
The sinusoidal magnetic flux distribution in the air gap of an IM(ii)
The constancy of the parameters of windings while the operation(iii)
Neglecting the effects of hysteresis and eddy currents (considering no iron losses)(iv)
Neglecting the mechanical losses (considering the equality of mechanical and electromechanical torques)(v)
Neglecting the skin effectThese assumptions are common for the representation of an IM in the form of an equivalent circuit. Therefore, in theabc reference frame, the electrical-balance equation for a stator phase of an IM is depicted as(1)uA=rAiA+dψAdt.Since the flux-linkage of a stator phase is a function of the stator (iA, iB, iC) and the rotor (ia, ib, ic) currents and of the angle between the corresponding axes of stator and rotor ϑ, the flux-linkage derivative in (1) is(2)dψAdt=∂ψA∂iA⋅diAdt+∂ψA∂iB⋅diBdt+∂ψA∂iC⋅diCdt+∂ψA∂ia⋅diadt+∂ψA∂ib⋅dibdt+∂ψA∂ic⋅dicdt+∂ψA∂ϑ⋅dϑdt.The partial derivatives with respect to currents in (2) are the dynamical inductances. Neglecting nonlinearities in a magnetic path of an IM, they are represented as static inductances: the self-inductance of a stator A-phase LAA, mutual inductances between the stator phases LAB and LAC, and mutual inductances between the stator and rotor LAa, LAb, and LAc, respectively. The item,(3)∂ψA∂ϑ⋅dϑdt=∂ψA∂ϑzpω,is the rotational e.m.f.The self-inductance of a statorA-phase LAA includes the components from the linkage magnetic flux and the leakage magnetic flux. Taking into account that the inductance from the main magnetic flux for a stator phase Lm and the inductance of a single phase Lms(4)Lm=32Lms,the similar ratio for the leakage inductance Lσs [9](5)LAA=23Lm+Lσs,where Lm and Lσs are the parameters of a traditional equivalent circuit of an IM. It is obvious that an IM is represented by an equivalent circuit considering the symmetrical machine and the symmetrical power supply. The parameters of the equivalent circuit are conveniently exploited in the modeling of an IM, and the methodology of determining those parameters is obvious.The equation for a mutual inductance between stator windings is calculated under the assumption of their symmetry:(6)LAB=LAC=LAAcos2π3=−13Lm+Lσs.The mutual inductances between the stator and the rotorLAa, LAb, and LAc are the functions of the angle ϑ, and the mutual inductance Lm between the stator and the rotor being calculated for ϑ=0 is(7)LAa=23kLmcosϑ,LAb=23kLmcosϑ+2π3,LAc=23kLmcosϑ−2π3,where the rotor winding is referred to the stator winding with the turn ratio k=Nr/Ns, where Ns and Nr are the numbers of turns in the stator and the rotor windings.The equation for a rotora phase (referred to a stator winding) is similar to (1):(8)ua=raia+dψadt,where the equation for the flux-linkage derivative dψa/dt is similar to (2), including a rotor phase’s self-inductance Laa, mutual inductances between the rotor windings Lab, Lac, and mutual inductances between the rotor and the stator LaA, LaB, and LaC.The self-inductance of a rotor phase is calculated similarly to (5) considering the turn ratio k as(9)Laa=23k2Lm+Lσr.The mutual inductance between the rotor windings considering their symmetry is(10)Lab=Lac=Laacos2π3=−13k2Lm+Lσr.The mutual inductances between the rotor and the stator are calculated similarly to (7).Other phases of an IM are depicted by the equations similar to (1), (2), and (5)–(10).


## 2.2. Representing the Equations of a Three-Phase IM in a Matrix Form

The equations of the electrical balance for stator and rotor circuits, (1), (2), and (5)–(10), are simply represented in a matrix form:(11)us=Rsis+ddtLsis+Lsrir+∂ψs∂ϑzpω,ur=Rrir+ddtLrsis+Lrir+∂ψr∂ϑzpω,where us=uAuBuCT, ur=uaubucT, is=iAiBiCT, ir=iaibicT, ψs=ψAψBψCT, and ψr=ψaψbψcT are the vectors of voltages, currents, and flux-linkages on the terminals of the stator and rotor winding; Rs=diagrA,rB,rC and Rr=diagra,rb,rc are the diagonal resistance matrices of stator and rotor windings. The elements of stator and rotor self-inductance matrices Ls and Lr and matrices of mutual inductances between the stator and the rotor Lsr and Lrs in (11) assuming the symmetrical IM are defined after (5)–(7) and (9)-(10).Thus, the matrices of mutual inductances between the stator and the rotor are(12)Lsr=LrsT=LAaLAbLAcLBaLBbLBcLCaLCbLCc,and self-inductance matrices Ls and Lr in the case of symmetrical IM are(13)Ls=23Lm+Lσs1−12−12−121−12−12−121,Lr=23k2Lm+Lσr1−12−12−121−12−12−121.Note that for a symmetrical IM, the matricesLs and Lr in (13) are singular. Thereafter, for a symmetrical IM, the calculation of the inverse matrix (related to the operation of solving the matrix equation) is impossible. We propose avoiding this issue as it is demonstrated in Appendix A.


## 3. IM Model Utilizing the Method of AVIS

### 3.1. General Description of the Method

Assume a branch of an electric circuit containing serial-connected e.m.f.e, resistance R, inductance L, and capacitance C with voltage u on its terminals (Figure 1(a)). The method of AVIS is based on the representation of electrical balance at the branch in the form of the balance of average voltages at the integration step Δt starting with t0 [18]:(14)1Δt∫t0t0+Δtu−e−uR−uC−uLdt=0,where the integrals from u, e, uR, uC, and uL define the average values of terminal voltages and voltages at the elements of the branch at the integration step.Figure 1
The branch of an electric circuit (a) and its equivalent circuit at the integration step, according to the method of AVIS (b).According to the method, equation (14) is algebraized, considering the current waveform at the integration step between i0 and i1 as an m-order polynomial. Therefore, (14) following [18] is transformed into(15)U+E−uR0−uC0+Rm+1+ΔtC⋅2−m+1m+22m+1m+2+L0Δti0−∑k=1m−1RΔtkk+1!⋅m−km+1+Δtk+1Ck+2!⋅m+1m+2−k+1k+2m+1m+2⋅dki0dtk−Rm+1+ΔtCm+1m+2+L1Δti1=0,where i0, uR0, and uC0 are the currents and the voltages at the resistance and the capacitance at the beginning of the step; L0 and L1 are the branch inductances at the beginning and the end of the step; and U and E are the average values at the step of the applied voltage and electromotive force: U=1/Δt∫t0t0+Δtudt and E=1/Δt∫t0t0+Δtedt. The unknown quantities in (15) are the currents at the end of the integration step i1 and the average value of a terminal voltage.Consequently, if the values ofEs and Rs are calculated at each integration step after (16) and (17) as [18],(16)Es=U+E−uR0−uC0+Rm+1+ΔtC⋅2−m+1m+22m+1m+2+L0Δti0−∑k=1m−1RΔtkk+1!⋅m−km+1+Δtk+1Ck+2!⋅m+1m+2−k+1k+2m+1m+2⋅dki0dtk,(17)Rs=Rm+1+ΔtCm+1m+2+L1Δt,then the branch in Figure 1(a) is represented at the integration step by an equivalent electromotive force Es and an equivalent resistance Rs connected in series with a current i1 flowing through the branch. That assumption allows representing the branch in Figure 1(a) in the form of an equivalent circuit at the integration step as in Figure 1(b).


### 3.2. Reformulation of the Matrix Equations of an IM with the Use of the Method of AVIS

Consider IM as a multiport network (Figure2(a)) with its nodes as IM winding terminals (in a squirrel-cage IM, the rotor nodes are short-circuited). Therefore, the matrix equations of an IM (11) concerning (15) are written as(18)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt+Rmim0−Rmm+1im0−∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim1dtk−Rmm+1Rm1−1Δtψm1−ψm0=0,where v1=vA1vB1vC1va1vb1vc1T and v2=vA2vB2vC2va2vb2vc2T are vectors of electrical potentials on IM terminals; im0=iA0iB0iC0ia0ib0ic0T is the vector of stator and rotor currents at the beginning of the integration step, and im1=iA1iB1iC1ia1ib1ic1T is the vector of the same currents at the end of the integration step; Rm=diagrA,rB,rC,ra,rb,rc is the diagonal matrix of resistances of the stator and rotor winding; ψm0=ψA0ψB0ψC0ψa0ψb0ψc0T and ψm1=ψA1ψB1ψC1ψa1ψb1ψc1T are the vectors of flux-linkages of the stator and rotor winding at the beginning and the end of the integration step, respectively.Figure 2
Representation of an IM as an electric multiport network (a) and an equivalent multiport network at the integration step according to the method of AVIS (b).
(a)
(b)Assumingψm0=Lmim0 and ψm1=Lmim1, where Lm is the matrix of self- and mutual inductances of stator and rotor windings of an IM(19)Lm=LsLsrLrsLr.(18) is transformed into(20)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt+Rmim0−Rmm+1im0−∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim1dtk−Rmm+1Rm1−1ΔtLmim1+1ΔtLmim0=0,where the unknown quantity is the vector of currents at the end of the integration step im1.Consequently, the representation (20) of the IM equations contains no components associated with the rotational e.m.f., as far as electromagnetic energy conversion is accounted with the difference in flux-linkage.Following the method of AVIS, each branch of a multiport network is represented by its equivalent circuit at the integration step as it is demonstrated in Figure1. Thus, a multiport network representing an IM (Figure 2(a)) is transformed into the equivalent multiport network at the integration step in Figure 2(b). Thus, (18) for the network in Figure 2(b) is transformed into(21)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt−Rsim1−es=0,where the vector of equivalent e.m.f. es and the matrix of equivalent branch resistances Rs are defined as(22)es=m⋅Rmm+1im0+∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim0dtk−1ΔtLmim0,(23)Rs=Rmm+1−1ΔtLm.Besides (21), the IM model requires the mechanical component(24)ωm1=ωm0+1JTm−TnΔt,resulting from the equation of the mechanical state of IM and(25)ϑm1=ϑm0+zpωm0Δt,for the rotor rotation angle, where Tm and Tn are the average values at the integration step of electromagnetic torque and loading torque, J is the rotor inertia, and zp is the number of pole pairs.The electromagnetic torque of the symmetrical machine is defined as(26)Tm=32p0Lmirβisα−irαisβ,where isα, isβ, irα, and irβ are stator and rotor currents represented in a stationary αβ reference frame with α-axis aligned with the A phase axis and β-axis perpendicular to it:(27)isα=23iAcos0+iBcos2π3+iCcos−2π3,isβ=23iAsin0+iBsin2π3+iCsin−2π3,irα=23iacosϑ+ibcosϑ−2π3+iccosϑ+2π3,irβ=23iasinϑ+ibsinϑ−2π3+icsinϑ+2π3.Therefore, a mathematical model of an IM according to the method of AVIS contains equations (21) and (24)–(27).


### 3.3. Algorithm of Solving the Equations of an IM at the Integration Step

The equations of an IM following the method of AVIS are computed after the algorithm:(1)
For a new rotor position (rotation angleϑ), the matrices of mutual inductances Lsr and Lsr are calculated.(2)
Derivatives of the currents at the beginning of the integration stepdki0/dtk are defined. In general, the m-order method requires derivatives up to m−1 order inclusive.(3)
The vector of equivalent electromotive forceses is defined after (22).(4)
The matrix of equivalent resistancesRs is defined after (23).(5)
Average values of potentials of IM outer terminals at the integration step (vectorsv1 and v2) are defined.(6)
After (21), the currents at the end of the integration step are defined.(7)
Other variables are calculated after (24)–(27).Here, the input quantities for the algorithm are the values of electric potentials at the terminals of an IMv1 and v2, currents im0, and loading torque Tn at the beginning of the integration step.


## 3.1. General Description of the Method

Assume a branch of an electric circuit containing serial-connected e.m.f.e, resistance R, inductance L, and capacitance C with voltage u on its terminals (Figure 1(a)). The method of AVIS is based on the representation of electrical balance at the branch in the form of the balance of average voltages at the integration step Δt starting with t0 [18]:(14)1Δt∫t0t0+Δtu−e−uR−uC−uLdt=0,where the integrals from u, e, uR, uC, and uL define the average values of terminal voltages and voltages at the elements of the branch at the integration step.Figure 1
The branch of an electric circuit (a) and its equivalent circuit at the integration step, according to the method of AVIS (b).According to the method, equation (14) is algebraized, considering the current waveform at the integration step between i0 and i1 as an m-order polynomial. Therefore, (14) following [18] is transformed into(15)U+E−uR0−uC0+Rm+1+ΔtC⋅2−m+1m+22m+1m+2+L0Δti0−∑k=1m−1RΔtkk+1!⋅m−km+1+Δtk+1Ck+2!⋅m+1m+2−k+1k+2m+1m+2⋅dki0dtk−Rm+1+ΔtCm+1m+2+L1Δti1=0,where i0, uR0, and uC0 are the currents and the voltages at the resistance and the capacitance at the beginning of the step; L0 and L1 are the branch inductances at the beginning and the end of the step; and U and E are the average values at the step of the applied voltage and electromotive force: U=1/Δt∫t0t0+Δtudt and E=1/Δt∫t0t0+Δtedt. The unknown quantities in (15) are the currents at the end of the integration step i1 and the average value of a terminal voltage.Consequently, if the values ofEs and Rs are calculated at each integration step after (16) and (17) as [18],(16)Es=U+E−uR0−uC0+Rm+1+ΔtC⋅2−m+1m+22m+1m+2+L0Δti0−∑k=1m−1RΔtkk+1!⋅m−km+1+Δtk+1Ck+2!⋅m+1m+2−k+1k+2m+1m+2⋅dki0dtk,(17)Rs=Rm+1+ΔtCm+1m+2+L1Δt,then the branch in Figure 1(a) is represented at the integration step by an equivalent electromotive force Es and an equivalent resistance Rs connected in series with a current i1 flowing through the branch. That assumption allows representing the branch in Figure 1(a) in the form of an equivalent circuit at the integration step as in Figure 1(b).


## 3.2. Reformulation of the Matrix Equations of an IM with the Use of the Method of AVIS

Consider IM as a multiport network (Figure2(a)) with its nodes as IM winding terminals (in a squirrel-cage IM, the rotor nodes are short-circuited). Therefore, the matrix equations of an IM (11) concerning (15) are written as(18)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt+Rmim0−Rmm+1im0−∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim1dtk−Rmm+1Rm1−1Δtψm1−ψm0=0,where v1=vA1vB1vC1va1vb1vc1T and v2=vA2vB2vC2va2vb2vc2T are vectors of electrical potentials on IM terminals; im0=iA0iB0iC0ia0ib0ic0T is the vector of stator and rotor currents at the beginning of the integration step, and im1=iA1iB1iC1ia1ib1ic1T is the vector of the same currents at the end of the integration step; Rm=diagrA,rB,rC,ra,rb,rc is the diagonal matrix of resistances of the stator and rotor winding; ψm0=ψA0ψB0ψC0ψa0ψb0ψc0T and ψm1=ψA1ψB1ψC1ψa1ψb1ψc1T are the vectors of flux-linkages of the stator and rotor winding at the beginning and the end of the integration step, respectively.Figure 2
Representation of an IM as an electric multiport network (a) and an equivalent multiport network at the integration step according to the method of AVIS (b).
(a)
(b)Assumingψm0=Lmim0 and ψm1=Lmim1, where Lm is the matrix of self- and mutual inductances of stator and rotor windings of an IM(19)Lm=LsLsrLrsLr.(18) is transformed into(20)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt+Rmim0−Rmm+1im0−∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim1dtk−Rmm+1Rm1−1ΔtLmim1+1ΔtLmim0=0,where the unknown quantity is the vector of currents at the end of the integration step im1.Consequently, the representation (20) of the IM equations contains no components associated with the rotational e.m.f., as far as electromagnetic energy conversion is accounted with the difference in flux-linkage.Following the method of AVIS, each branch of a multiport network is represented by its equivalent circuit at the integration step as it is demonstrated in Figure1. Thus, a multiport network representing an IM (Figure 2(a)) is transformed into the equivalent multiport network at the integration step in Figure 2(b). Thus, (18) for the network in Figure 2(b) is transformed into(21)1Δt∫t0t0+Δtv1dt−1Δt∫t0t0+Δtv2dt−Rsim1−es=0,where the vector of equivalent e.m.f. es and the matrix of equivalent branch resistances Rs are defined as(22)es=m⋅Rmm+1im0+∑k=1m−1RmΔtkk+1!⋅m−km+1⋅dkim0dtk−1ΔtLmim0,(23)Rs=Rmm+1−1ΔtLm.Besides (21), the IM model requires the mechanical component(24)ωm1=ωm0+1JTm−TnΔt,resulting from the equation of the mechanical state of IM and(25)ϑm1=ϑm0+zpωm0Δt,for the rotor rotation angle, where Tm and Tn are the average values at the integration step of electromagnetic torque and loading torque, J is the rotor inertia, and zp is the number of pole pairs.The electromagnetic torque of the symmetrical machine is defined as(26)Tm=32p0Lmirβisα−irαisβ,where isα, isβ, irα, and irβ are stator and rotor currents represented in a stationary αβ reference frame with α-axis aligned with the A phase axis and β-axis perpendicular to it:(27)isα=23iAcos0+iBcos2π3+iCcos−2π3,isβ=23iAsin0+iBsin2π3+iCsin−2π3,irα=23iacosϑ+ibcosϑ−2π3+iccosϑ+2π3,irβ=23iasinϑ+ibsinϑ−2π3+icsinϑ+2π3.Therefore, a mathematical model of an IM according to the method of AVIS contains equations (21) and (24)–(27).


## 3.3. Algorithm of Solving the Equations of an IM at the Integration Step

The equations of an IM following the method of AVIS are computed after the algorithm:(1)
For a new rotor position (rotation angleϑ), the matrices of mutual inductances Lsr and Lsr are calculated.(2)
Derivatives of the currents at the beginning of the integration stepdki0/dtk are defined. In general, the m-order method requires derivatives up to m−1 order inclusive.(3)
The vector of equivalent electromotive forceses is defined after (22).(4)
The matrix of equivalent resistancesRs is defined after (23).(5)
Average values of potentials of IM outer terminals at the integration step (vectorsv1 and v2) are defined.(6)
After (21), the currents at the end of the integration step are defined.(7)
Other variables are calculated after (24)–(27).Here, the input quantities for the algorithm are the values of electric potentials at the terminals of an IMv1 and v2, currents im0, and loading torque Tn at the beginning of the integration step.


## 4. Example of the Application of a Designed Model for the Analysis of EMS

In this section, we would like to demonstrate the design of a model of an EMS on an example from our practice. The model is based on a model of a three-phase IM utilizing the method of AVIS.


### 4.1. Case Study: Doubly Fed Induction Machine with the Full Loss of Operation of a Grid-Side Converter

Figure3(a) demonstrates the topology of a DFIM with a back-to-back converter (containing the rotor-side converter (RSC) and the grid-side converter (GSC)). In our practice, we had to analyze the operation of a DFIM with the full loss of operation of a GSC when no power supply is provided to the IM rotor circuit. The equivalent topology for the mode is demonstrated in Figure 3(b). The objective of the study was to analyze the operation of a DFIM without additional power supply in the rotor circuit and the limitations of the operation. Although the detailed description of the study comes outside the scope of the present work, we would like to demonstrate the numerical and computer models’ design algorithm for this case study. The analytical expressions and experiments related to the case study were provided by us in [23].Figure 3
Topology of a DFIM, where RSC is a rotor-side converter and GSC is a grid-side converter of a back-to-back converter (a); the equivalent topology of the EMS with the full loss of operation in GSC (b).


### 4.2. Numerical and Computer Models of the EMS

The mathematical model of the analyzed EMS in Figure3(b) is based on the ideas described in [22]. Thus, the model of a complicated EMS is constructed from the models of the structural elements taken separately; each of them is represented as a multiport network (power network, IM, switch groups, and DC link), as it is demonstrated in Figure 4(a). Switch groups contain branches of individual semiconductor switches; each of them is represented by the serial-connected resistance and inductance. In the “on” state, those resistance and inductance have the low values, and in the “off” state—the high values. The details of the representation and the procedure of defining the switches’ turn-on and turn-off instants have been previously reported by us in [22]. Note that the group of IGBT switches is represented by two individual groups, one of them representing the switches themselves and the other—antiparallel diodes included in an IGBT package (in Figure 4(b), these are switch groups 1 and 3 and switch groups 2 and 4).Figure 4
Computational equivalent diagram of the topology in Figure3(b) as a combination of multiport networks (a) and computational equivalent diagram at the integration step following the AVIS (b).
(a)
(b)Following the method of AVIS, each branch of the networks is represented at the integration step by equivalent e.m.f. and resistance connected in series. Thus, the equivalent circuit diagram of the EMS at the integration step is represented in a form as in Figure4(b).The description of the equations required for the numerical modeling of EMS is provided by [22]. The algorithm of solving the equations is similar to the modified nodal analysis; however, the average values of the independent node potentials are utilized instead of the independent node potentials.Thus, each multiport network is depicted by the terminal matrix equation(28)ie+Gse1Δt∫t0t0+Δtvedt+cse=0,where ve is a vector of potentials of outer terminals; ie is a vector of currents in external branches of a multiport network; and Gse and cse are a matrix of coefficients and a vector of absolute terms. This equation is written after (15) for each element.The design of a model of EMS is reduced to the coupling of single multiport networks at their nodes following the modified nodal analysis. Incidence matricesΠ are used to define the way the networks are interconnected within the EMS according to the relation ve=ΠTvc, where ve is a vector of the potentials of external terminals of the network and vc is a vector of the potentials of independent nodes of the whole EMS. The same relation is written for average values at the integration step:(29)1Δt∫t0t0+Δtvedt=ΠΤ1Δt∫t0t0+Δtvcdt.The average values of the potentials of independent nodes1/Δt∫t0t0+Δtvcdt at the integration step are defined after the algebraic matrix equation(30)Gsc1Δt∫t0t0+Δtvcdt+csc=0.In (30), the relations between Gsc and Gse and between csc and cse for each multiport network are defined as(31)Gsc=∑j=1LΠjGsejΠjT,csc=∑j=1LΠjcsej,where L is the quantity of single multiport networks within the EMS.The algorithm of solving the equations for EMS is as follows:(i)
First, (29) is used to calculate the average values of the potentials of independent nodes(ii)
Afterward, the average values of terminal potentials for each network are defined after (30)(iii)
And finally, the currents in external branches of each network at the end of the step are defined after (28)The computer model of the EMS is designed in C++ utilizing the object-oriented environment developed under the supervising of Prof. Plakhtyna with the participation of the author [22]. Here, we explain the general idea of the environment using the definitions and terms in the form they provided in [22].Thus, constructing the model of EMS reduces to the declaration of the developed ready-to-use objects (the models of the standard elements) and to setting the communication between them. In [22], we describe the developed hierarchical class structure: the basic superclass is the basis for the design of the classes representing the standard elements; the basic superclass includes the methods required for the functioning of all objects (models of structural elements). In particular, these are the methods realizing the matrix operations (28)–(31). The derived classes represent the standard elements (in the analyzed EMS, those are the power network, the IM, the switch groups, and the DC link). Those are inherited from the basic superclass. The separate classes are designed for a power circuit as a whole and a control system.The detailed description of the ideology and the principles of the communication between the objects within the developed object-oriented environment are offered in [22]. Thus, Figure 5 depicts the structure of the objects and their communication while functioning. The communication between the objects is provided by sending the messages, activating the appropriate methods (in Figure 5, those are the arrows marked as CNTR, DIFF, SAR, DY, and GsC). The numbers near the arrows represent the order of sending messages. For the detailed description, we recommend referring to [22].Figure 5
The diagram of objects within an EMS after [22].


## 4.1. Case Study: Doubly Fed Induction Machine with the Full Loss of Operation of a Grid-Side Converter

Figure3(a) demonstrates the topology of a DFIM with a back-to-back converter (containing the rotor-side converter (RSC) and the grid-side converter (GSC)). In our practice, we had to analyze the operation of a DFIM with the full loss of operation of a GSC when no power supply is provided to the IM rotor circuit. The equivalent topology for the mode is demonstrated in Figure 3(b). The objective of the study was to analyze the operation of a DFIM without additional power supply in the rotor circuit and the limitations of the operation. Although the detailed description of the study comes outside the scope of the present work, we would like to demonstrate the numerical and computer models’ design algorithm for this case study. The analytical expressions and experiments related to the case study were provided by us in [23].Figure 3
Topology of a DFIM, where RSC is a rotor-side converter and GSC is a grid-side converter of a back-to-back converter (a); the equivalent topology of the EMS with the full loss of operation in GSC (b).


## 4.2. Numerical and Computer Models of the EMS

The mathematical model of the analyzed EMS in Figure3(b) is based on the ideas described in [22]. Thus, the model of a complicated EMS is constructed from the models of the structural elements taken separately; each of them is represented as a multiport network (power network, IM, switch groups, and DC link), as it is demonstrated in Figure 4(a). Switch groups contain branches of individual semiconductor switches; each of them is represented by the serial-connected resistance and inductance. In the “on” state, those resistance and inductance have the low values, and in the “off” state—the high values. The details of the representation and the procedure of defining the switches’ turn-on and turn-off instants have been previously reported by us in [22]. Note that the group of IGBT switches is represented by two individual groups, one of them representing the switches themselves and the other—antiparallel diodes included in an IGBT package (in Figure 4(b), these are switch groups 1 and 3 and switch groups 2 and 4).Figure 4
Computational equivalent diagram of the topology in Figure3(b) as a combination of multiport networks (a) and computational equivalent diagram at the integration step following the AVIS (b).
(a)
(b)Following the method of AVIS, each branch of the networks is represented at the integration step by equivalent e.m.f. and resistance connected in series. Thus, the equivalent circuit diagram of the EMS at the integration step is represented in a form as in Figure4(b).The description of the equations required for the numerical modeling of EMS is provided by [22]. The algorithm of solving the equations is similar to the modified nodal analysis; however, the average values of the independent node potentials are utilized instead of the independent node potentials.Thus, each multiport network is depicted by the terminal matrix equation(28)ie+Gse1Δt∫t0t0+Δtvedt+cse=0,where ve is a vector of potentials of outer terminals; ie is a vector of currents in external branches of a multiport network; and Gse and cse are a matrix of coefficients and a vector of absolute terms. This equation is written after (15) for each element.The design of a model of EMS is reduced to the coupling of single multiport networks at their nodes following the modified nodal analysis. Incidence matricesΠ are used to define the way the networks are interconnected within the EMS according to the relation ve=ΠTvc, where ve is a vector of the potentials of external terminals of the network and vc is a vector of the potentials of independent nodes of the whole EMS. The same relation is written for average values at the integration step:(29)1Δt∫t0t0+Δtvedt=ΠΤ1Δt∫t0t0+Δtvcdt.The average values of the potentials of independent nodes1/Δt∫t0t0+Δtvcdt at the integration step are defined after the algebraic matrix equation(30)Gsc1Δt∫t0t0+Δtvcdt+csc=0.In (30), the relations between Gsc and Gse and between csc and cse for each multiport network are defined as(31)Gsc=∑j=1LΠjGsejΠjT,csc=∑j=1LΠjcsej,where L is the quantity of single multiport networks within the EMS.The algorithm of solving the equations for EMS is as follows:(i)
First, (29) is used to calculate the average values of the potentials of independent nodes(ii)
Afterward, the average values of terminal potentials for each network are defined after (30)(iii)
And finally, the currents in external branches of each network at the end of the step are defined after (28)The computer model of the EMS is designed in C++ utilizing the object-oriented environment developed under the supervising of Prof. Plakhtyna with the participation of the author [22]. Here, we explain the general idea of the environment using the definitions and terms in the form they provided in [22].Thus, constructing the model of EMS reduces to the declaration of the developed ready-to-use objects (the models of the standard elements) and to setting the communication between them. In [22], we describe the developed hierarchical class structure: the basic superclass is the basis for the design of the classes representing the standard elements; the basic superclass includes the methods required for the functioning of all objects (models of structural elements). In particular, these are the methods realizing the matrix operations (28)–(31). The derived classes represent the standard elements (in the analyzed EMS, those are the power network, the IM, the switch groups, and the DC link). Those are inherited from the basic superclass. The separate classes are designed for a power circuit as a whole and a control system.The detailed description of the ideology and the principles of the communication between the objects within the developed object-oriented environment are offered in [22]. Thus, Figure 5 depicts the structure of the objects and their communication while functioning. The communication between the objects is provided by sending the messages, activating the appropriate methods (in Figure 5, those are the arrows marked as CNTR, DIFF, SAR, DY, and GsC). The numbers near the arrows represent the order of sending messages. For the detailed description, we recommend referring to [22].Figure 5
The diagram of objects within an EMS after [22].


## 5. Comparison of the Behavior of the IM Models Developed with the Use of Different Numerical Methods

### 5.1. Set of the Numerical Methods for the Comparison

Numerical models of an IM directly connected to a power network have been designed utilizing the set of the methods of numerical integration [21]:(1)
Second-order Runge–Kutta method (RK2) because of the sufficient accuracy of the results provided by a small number of calculations(2)
Fourth-order Adams–Bashforth (AB4) and fourth-order Adams–Moulton (AM4) methods, as far as they provide the best operation speed and accuracy for the modeling of electric circuits, according to the comparison in [18](3)
First- and second-order AVIS (AVIS1 and AVIS2)The set of the computer models has been used for the simulation of a direct starting mode of an IM neglecting the saturation (Lm=const) with the different values of the integration step. The parameters of the machine applied in those models are presented in Appendix B.


### 5.2. Determining the “Standard” Values and the Numerical Validation of the Developed Models of an IM

The numerical validation of a method of numerical integration is often held by the comparison of the result obtained with the use of the method and an exact analytical solution of the equation. The same technique has been used in [18], where the method of AVIS has been compared with conventional numerical methods for the analysis of the simple DC and AC electric circuits.Due to the higher complexity of the IM equations, we propose another technique for the validation of the models [21]. The first stage of the numerical experiment was selecting the standard simulation results (here, the standard waveforms of speed, electromagnetic torque, and stator and rotor currents); the standard results were considered as the fictional exact analytical solutions of the IM equations. We assume that with a very low value of the integration step, all methods provide the same results (that means the results for the numerical methods are characterized by the minimal difference). That means, for each calculated point, the results, obtained with the use of different numerical methods, were compared in the terms of relative error. We also consider the fact that the standard result exists for the set of the methods to be the validation of the numerical models of an IM.Thus, the simulation results for the integration step of10−6s provided the maximum relative error for all quantities below 10−3%, while for most compared quantities, it was below 10−6%. Those results were further considered as standard results. The further simulation results obtained with the higher values of the integration step have been compared with the standard results to obtain a numerical evaluation of the results in terms of relative error [21].


### 5.3. Comparison of the Behavior of Different Numerical Methods for the Analysis of an IM

Afterward, we proceed with higher values of the integration step. The simulated results were compared with thestandard values for each calculated point, so as each quantity x is characterized by an average relative error referred to the steady-state value [21]:(32)δx=1N⋅∑i=1Nxstdti−xcalctixss⋅100%,where xstdti is a standard value at ti (at the i-th point of calculation), xcalcti is a calculated value of the quantity obtained with the use of one of the methods, xss is a value of x in a steady-state or a steady-state amplitude for the stator and rotor currents, and N is the number of the calculation points.A significant feature of the rotor current’s error estimation is caused by the dependence of its frequency on the slip. Thus, an insignificant difference in the rotational speed may cause higher values of average relative errors after (32). Here, in Figures 6 and 7, we demonstrate an example in the case of 1.5×10−3s integration step.Figure 6
Simulated quantities of the IM for the integration step of1.5×10−3s for a direct starting mode.Figure 7
Enlarged graph of the rotor current waveforms from Figure6.Both Adams methods are unstable for that value of an integration step. The results provided by the other methods are satisfactory for the engineering practice, while the worst of them are the results obtained using the RK2 method. However, while the results obtained with the use of AVIS1 are closer to thestandard waveforms (according to Figures 6 and 7); the numerical evaluation of the errors in calculating the rotor currents after (32) is 140% for AVIS1 and 110% for RK2.Therefore, the behavior of the rotor current is unclearly evaluated using (32). On the contrary, the parameters of the EMS elements are known with the accuracies of 1,…,5% [24]. Consequently, for the engineering applications, we recommend considering the general behavior of the periodical quantities of IM by the comparison of waveforms and envelopes instead of instantaneous values [21]; thus, integrals for each quantity xcalcti are calculated as(33)Sx=∑i=1NxcalctiΔt.The values ofSx were afterward compared to the value of Sxstd for a standard waveform xstdti to obtain an “integral assessment” of a quantity; it has been used in addition to the error for the instantaneous values being calculated after (32). We consider the integral assessment as sufficient for the engineering application characteristic of the quantities of an IM.


### 5.4. Numerical Comparison of the Methods of Numerical Integration for the Analysis of an IM

The further experiment was held as follows: starting with the integration step of10−6s, we calculated the maximum relative errors for the instantaneous values and the integral assessment [21]. In Figure 8, the step sizes providing a maximum relative error of 1%, 5%, and 10% for both criteria are summarized.Figure 8
Comparison of the steps providing maximum relative errors of1% (lilac), 5% (orange), and 10% (red) for different numerical methods for instantaneous values (a) and for integral assessment (b).
(a)
(b)We also compared the simulation duration for different numerical methods for an integration step of10−6s and maximum integration step size providing the stable results. In Figures 9 and 10, those results are demonstrated. Thus, the difference between the simulation duration for Adams methods and AVIS2 was less than 2% (150 s and 149 s for AB4 and AM4, respectively, and 153 s for AVIS2), while RK2-based simulation lasted 158 s. As a result, the difference between all those values excluding AVIS1 was below 8%, while AVIS1 itself provided an increase in the operation speed of 30% (up to 90 s).Figure 9
Comparison of the simulation duration for different numerical methods for an integration step of10−6s.Figure 10
Comparison of maximum integration step sizes, providing stable results.The result for both AVIS methods is explained as follows: according to (13), the m-order method requires the information about the derivatives up to m−1 order. Thereby, AVIS1 requires no derivative calculation, AVIS2 requires the calculation of only first-order derivative, AVIS3 requires the calculation of first- and second-order derivatives, and so on. In other words, the increase in the order of the method results in the increase in the number of calculations and in the simulation duration.In Figure10, the maximum integration step sizes, providing stable results for all methods are compared. The figure demonstrates the most significant feature of the method—its ability to operate with high values of the integration step providing a stable calculation [21].


### 5.5. Discussion on the Behavior of the Methods of Numerical Integration under Different Values of the Integration Step

The behavior of the models developed with the use of different numerical methods demonstrated earlier in this section circumscribes the field where the proposed numerical models of an IM are perspective.The distinguishing feature of the Adams methods is their nonaccumulation of error [25]. On the contrary, AM4, as well as AB4, are characterized by the narrowest stability regions among the compared methods (Figure 10). Therefore, while they are well-suited for offline computation, real-time applications require other models. Thus, AM4 is characterized by the highest accuracy for instantaneous values considering the same integration step (Figure 8(a)), while the calculation time is within the same range as most of the compared methods excluding AVIS1 (Figure 9).However, it is important to consider the fact that is demonstrated in Figures6 and 7—the errors for the instantaneous values unclearly characterize the accuracy of the results. With the integral assessment (we consider it as a sufficient characteristic for most engineering applications), AVIS2 provides the ability to increase the integration step 7.5 times as compared to AM4 (for the same 10% error, Figure 8(b)).On the contrary, AVIS1 provides an increase in the operation speed of 30% (Figure9), while it provides 4 times increase in the maximum integration step as compared to RK2. However, AVIS2 provides 2.5 times higher maximum integration step.Summarizing the results, we consider both first- and second-order AVIS as suitable tools for the high-speed applications. Among them, the second-order method is preferable, while the first-order method may be an alternative to increase the calculation speed with the acceptable decrease in the accuracy.


## 5.1. Set of the Numerical Methods for the Comparison

Numerical models of an IM directly connected to a power network have been designed utilizing the set of the methods of numerical integration [21]:(1)
Second-order Runge–Kutta method (RK2) because of the sufficient accuracy of the results provided by a small number of calculations(2)
Fourth-order Adams–Bashforth (AB4) and fourth-order Adams–Moulton (AM4) methods, as far as they provide the best operation speed and accuracy for the modeling of electric circuits, according to the comparison in [18](3)
First- and second-order AVIS (AVIS1 and AVIS2)The set of the computer models has been used for the simulation of a direct starting mode of an IM neglecting the saturation (Lm=const) with the different values of the integration step. The parameters of the machine applied in those models are presented in Appendix B.


## 5.2. Determining the “Standard” Values and the Numerical Validation of the Developed Models of an IM

The numerical validation of a method of numerical integration is often held by the comparison of the result obtained with the use of the method and an exact analytical solution of the equation. The same technique has been used in [18], where the method of AVIS has been compared with conventional numerical methods for the analysis of the simple DC and AC electric circuits.Due to the higher complexity of the IM equations, we propose another technique for the validation of the models [21]. The first stage of the numerical experiment was selecting the standard simulation results (here, the standard waveforms of speed, electromagnetic torque, and stator and rotor currents); the standard results were considered as the fictional exact analytical solutions of the IM equations. We assume that with a very low value of the integration step, all methods provide the same results (that means the results for the numerical methods are characterized by the minimal difference). That means, for each calculated point, the results, obtained with the use of different numerical methods, were compared in the terms of relative error. We also consider the fact that the standard result exists for the set of the methods to be the validation of the numerical models of an IM.Thus, the simulation results for the integration step of10−6s provided the maximum relative error for all quantities below 10−3%, while for most compared quantities, it was below 10−6%. Those results were further considered as standard results. The further simulation results obtained with the higher values of the integration step have been compared with the standard results to obtain a numerical evaluation of the results in terms of relative error [21].


## 5.3. Comparison of the Behavior of Different Numerical Methods for the Analysis of an IM

Afterward, we proceed with higher values of the integration step. The simulated results were compared with thestandard values for each calculated point, so as each quantity x is characterized by an average relative error referred to the steady-state value [21]:(32)δx=1N⋅∑i=1Nxstdti−xcalctixss⋅100%,where xstdti is a standard value at ti (at the i-th point of calculation), xcalcti is a calculated value of the quantity obtained with the use of one of the methods, xss is a value of x in a steady-state or a steady-state amplitude for the stator and rotor currents, and N is the number of the calculation points.A significant feature of the rotor current’s error estimation is caused by the dependence of its frequency on the slip. Thus, an insignificant difference in the rotational speed may cause higher values of average relative errors after (32). Here, in Figures 6 and 7, we demonstrate an example in the case of 1.5×10−3s integration step.Figure 6
Simulated quantities of the IM for the integration step of1.5×10−3s for a direct starting mode.Figure 7
Enlarged graph of the rotor current waveforms from Figure6.Both Adams methods are unstable for that value of an integration step. The results provided by the other methods are satisfactory for the engineering practice, while the worst of them are the results obtained using the RK2 method. However, while the results obtained with the use of AVIS1 are closer to thestandard waveforms (according to Figures 6 and 7); the numerical evaluation of the errors in calculating the rotor currents after (32) is 140% for AVIS1 and 110% for RK2.Therefore, the behavior of the rotor current is unclearly evaluated using (32). On the contrary, the parameters of the EMS elements are known with the accuracies of 1,…,5% [24]. Consequently, for the engineering applications, we recommend considering the general behavior of the periodical quantities of IM by the comparison of waveforms and envelopes instead of instantaneous values [21]; thus, integrals for each quantity xcalcti are calculated as(33)Sx=∑i=1NxcalctiΔt.The values ofSx were afterward compared to the value of Sxstd for a standard waveform xstdti to obtain an “integral assessment” of a quantity; it has been used in addition to the error for the instantaneous values being calculated after (32). We consider the integral assessment as sufficient for the engineering application characteristic of the quantities of an IM.


## 5.4. Numerical Comparison of the Methods of Numerical Integration for the Analysis of an IM

The further experiment was held as follows: starting with the integration step of10−6s, we calculated the maximum relative errors for the instantaneous values and the integral assessment [21]. In Figure 8, the step sizes providing a maximum relative error of 1%, 5%, and 10% for both criteria are summarized.Figure 8
Comparison of the steps providing maximum relative errors of1% (lilac), 5% (orange), and 10% (red) for different numerical methods for instantaneous values (a) and for integral assessment (b).
(a)
(b)We also compared the simulation duration for different numerical methods for an integration step of10−6s and maximum integration step size providing the stable results. In Figures 9 and 10, those results are demonstrated. Thus, the difference between the simulation duration for Adams methods and AVIS2 was less than 2% (150 s and 149 s for AB4 and AM4, respectively, and 153 s for AVIS2), while RK2-based simulation lasted 158 s. As a result, the difference between all those values excluding AVIS1 was below 8%, while AVIS1 itself provided an increase in the operation speed of 30% (up to 90 s).Figure 9
Comparison of the simulation duration for different numerical methods for an integration step of10−6s.Figure 10
Comparison of maximum integration step sizes, providing stable results.The result for both AVIS methods is explained as follows: according to (13), the m-order method requires the information about the derivatives up to m−1 order. Thereby, AVIS1 requires no derivative calculation, AVIS2 requires the calculation of only first-order derivative, AVIS3 requires the calculation of first- and second-order derivatives, and so on. In other words, the increase in the order of the method results in the increase in the number of calculations and in the simulation duration.In Figure10, the maximum integration step sizes, providing stable results for all methods are compared. The figure demonstrates the most significant feature of the method—its ability to operate with high values of the integration step providing a stable calculation [21].


## 5.5. Discussion on the Behavior of the Methods of Numerical Integration under Different Values of the Integration Step

The behavior of the models developed with the use of different numerical methods demonstrated earlier in this section circumscribes the field where the proposed numerical models of an IM are perspective.The distinguishing feature of the Adams methods is their nonaccumulation of error [25]. On the contrary, AM4, as well as AB4, are characterized by the narrowest stability regions among the compared methods (Figure 10). Therefore, while they are well-suited for offline computation, real-time applications require other models. Thus, AM4 is characterized by the highest accuracy for instantaneous values considering the same integration step (Figure 8(a)), while the calculation time is within the same range as most of the compared methods excluding AVIS1 (Figure 9).However, it is important to consider the fact that is demonstrated in Figures6 and 7—the errors for the instantaneous values unclearly characterize the accuracy of the results. With the integral assessment (we consider it as a sufficient characteristic for most engineering applications), AVIS2 provides the ability to increase the integration step 7.5 times as compared to AM4 (for the same 10% error, Figure 8(b)).On the contrary, AVIS1 provides an increase in the operation speed of 30% (Figure9), while it provides 4 times increase in the maximum integration step as compared to RK2. However, AVIS2 provides 2.5 times higher maximum integration step.Summarizing the results, we consider both first- and second-order AVIS as suitable tools for the high-speed applications. Among them, the second-order method is preferable, while the first-order method may be an alternative to increase the calculation speed with the acceptable decrease in the accuracy.


## 6. Conclusions

The equations of a three-phase induction machine in the naturalabc reference frame are reformulated utilizing the method of numerical integration of average voltages at the integration step (AVIS). However, the model of an IM in the abc reference frame requires the calculation of rotational e.m.f. at every integration step; utilizing the method of AVIS allows avoiding this calculation, while electromagnetic energy conversion is accounted by the change of flux-linkage.The designed model of an IM is considered as a multiport network and implemented in the object-oriented environment in C++ to design computer models of electromechanical systems (EMS) as it has been demonstrated in the work. The environment is proposed as a tool for designing the computer models of EMS based on an IM utilizing AVIS.The behavior of the designed numerical model of IM (utilizing the first- and the second-order methods) has been compared with the behavior of numerical models designed utilizing the conventional methods of numerical integration (the second-order Runge–Kutta method and the fourth-order Adams–Bashforth and Adams–Moulton methods). The models designed utilizing the first- and the second-order AVIS have demonstrated definite advantages over the conventional methods for the high-speed applications.



---

*Source: 4581769-2019-12-14.xml*