Accounting for Magnetic Saturation Effects in Complex Multi-harmonic Model of Induction Machine

Abstract

Computations of quasi-dynamic electromagnetic field of induction machines using the complex magnetic vector potential require the use of the so-called effective magnetization curves, i.e., such in which the magnetic permeability is proportional to the amplitudes of magnetic flux density B or magnetic field strength H, not their instantaneous values. There are several definitions of that parameter mentioned in the literature provided for the case when B or H are monoharmonic. In this paper, seven different methods of determining the effective magnetization curves are compared in relation to the use of a field-circuit multi-harmonic model of an induction machine. The accuracy of each method was assessed by computing the performance characteristics of a solid-rotor induction machine. One new definition of the effective permeability was also introduced, being a function of multiple variables dependent on amplitudes of all the harmonics considered. The analyses demonstrated that the best practical approach, even for the multi-harmonic case, is to express the effective magnetic permeability as the ratio of the amplitudes of the fundamental harmonics of the magnetic flux density and the magnetic field strength, and assuming the sinusoidal variation of the latter.

1. Introduction

Numerical models based on the use of the finite element method have become an indispensable tool in the design of electrical machines. As is well known, the most accurate approach for overall consideration of physical phenomena in an analyzed converter is the use of time-domain models. Unfortunately, this approach usually involves significant computational cost, sometimes making its practical application impossible. For that reason, it is justified to explore methods which allow, particularly in the analysis of steady states, to achieve similar results, but in a much shorter time. While analyzing induction machines, one of such methods is the use of the multi-harmonic field-circuit model [1,2,3,4,5,6,7,8,9]. The results presented in [1,2,9] prove that this type of the field-circuit model can be particularly useful when analyzing steady states of the special constructions of induction machines (single-phase machines, high-speed solid-rotor machines). As demonstrated by Garbiec et al., with an appropriate modification of the coupling scheme between individual sub-models, the multi-harmonic field-circuit model makes it possible to achieve comparable calculation results with those obtained using the time-domain model, at a calculation cost comparable to the monoharmonic model. The multi-harmonic model in this case is understood as the De-Gersem type model formulated using the magnetic vector potential, where permeance slot harmonics of the magnetic field distribution in the air gap are extracted via Fourier transform of an air-gap magnetic field. These harmonics can affect the power losses in the rotor and the electromagnetic torque developed by the machine. The magnetic field in the stator in these models is considered as a function of the fundamental time harmonic only.

The application of the above-mentioned approach essentially involves the problem of modeling the nonlinearity of the magnetic materials. Taking this into account in the monoharmonic field models based on the use of the complex magnetic vector potential method is not a new task and it has been described to date in several works, including [10,11,12,13,14,15,16]. They propose several different methods of defining the so-called effective magnetic permeability based, among others, on the development of nonlinear waveforms in the Fourier series [14,16], averaging the magnetic reluctivity [10,15,16] or equivalence of energy stored in ferromagnetic components [10,11,12,13]. However, analyzing the conclusions drawn in the mentioned works, it is not possible to clearly state which method of defining the effective magnetic permeability will be most appropriate when using the multi-harmonic model with a strong coupling, as proposed by Garbiec et al. Some of them assume the sinusoidal variation of the magnetic flux density in calculating the effective magnetic permeability (which facilitates the calculations using the magnetic vector potential) [12,13,14,15], while others assume that it should be determined assuming sinusoidal variation of the magnetic field strength, or that it is insignificant for the accuracy of the calculations [10,11,14]. The latter approach was adopted in the authors’ previous work, and satisfactory results were achieved [9]. It was based on the results of the earlier work of the authors, where the effective magnetic permeability defined with the assumption of the sinusoidal variation of the magnetic field strength provided correct results, modeling a solid-rotor induction machine with the use of the multi-harmonic field-circuit model with a weak coupling [17]. To the best of authors’ knowledge, the research providing an analysis of the accuracy of definitions of the effective magnetic permeability on the results of computations using multi-harmonic field-circuit models was not reported.

This paper is a continuation of works by Garbiec et al. that aim at an in-depth verification of the De-Gersem type multi-harmonic complex model of induction machines from the point of view of their application in designing and optimization of real machines. The authors focus on the analysis of the impact of formulating the effective magnetic permeability. For that purpose, the basic operating characteristics of an exemplary machine were calculated for seven different methods of formulating the effective magnetic permeability. The first six methods use existing definitions derived in literature. The available formulas use averaging of magnetic reluctivity in terms of mean or RMS value. Also, the Fourier series truncated to the fundamental component and assuming sinusoidal variation of the magnetic flux density or magnetic field strength is used [14,15,16]. The seventh method considered in this paper is a new approach proposed by the authors that relies upon formulating the so-called multidimensional effective magnetic permeability. The results of calculations were compared with those obtained with the time-domain model and with the results of measurements performed on the physical model under one of the previous works where the efficient computational method for fast determination of the performance characteristics of solid-rotor induction machines was presented [17].

2. Mathematical Model

2.1. Effective Magnetic Permeability

Definition of the effective magnetic permeability, as used in the previous works of the authors [9,14,16,17], is the reference point for the considerations presented herein:

$${\mu }_{\mathrm{I}}\left(H_{pk}\right)=\frac{2}{\pi H_{pk}}{{\displaystyle \int }}_0^{\pi }\left({\mu }_{DC}\left(H_{pk}\sin \alpha \right)H_{pk}\sin \alpha \right)\sin \alpha d\alpha ,$$

(1)

where $H_{pk}$—amplitude of the sinusoidal magnetic field strength H, and ${\mu }_{DC}$—DC magnetic permeability. The above definition corresponds to the ratio of the fundamental harmonic amplitude of the magnetic flux density waveform (obtained on the assumption of the sinusoidal variation of the magnetic field strength) and $H_{pk}$. This approach requires a successive determination of the magnetic field strength in a given iteration of the algorithm based on the knowledge of its value determined in the previous iteration.

Assuming the sinusoidal variation of the magnetic field strength, two further definitions can be introduced [16]:

$${\mu }_{\mathrm{II}}\left(H_{pk}\right)=\frac{\pi }{{{\displaystyle \int }}_0^{\pi }{\nu }_{DC}\left(H_{pk}\sin \alpha \right)d\alpha },$$

(2)

$${\mu }_{\mathrm{III}}\left(H_{pk}\right)=\frac{1}{\sqrt{\frac{1}{\pi }{{\displaystyle \int }}_0^{\pi }{\nu }_{DC}^2\left(H_{pk}\sin \alpha \right)d\alpha }},$$

(3)

where ${\nu }_{DC}=\frac{1}{{\mu }_{DC}}$. The expressions (2) and (3) correspond to the inverse of the mean value and the inverse of the RMS value of the magnetic reluctivity obtained for the half of the period.

Similarly, the definitions of the effective magnetic permeability can be derived, assuming the sinusoidal variation of the magnetic flux density B with the amplitude $B_{pk}$ [15,16]:

$${\mu }_{\mathrm{IV}}\left(B_{pk}\right)=\frac{\pi B_{pk}}{2{{\displaystyle \int }}_0^{\pi }\left({\nu }_{DC}\left(B_{pk}\sin \alpha \right)B_{pk}\sin \alpha \right)\sin \alpha d\alpha },$$

(4)

$${\mu }_{\mathrm{V}}\left(B_{pk}\right)=\frac{\pi }{{{\displaystyle \int }}_0^{\pi }{\nu }_{DC}\left(B_{pk}\sin \alpha \right)d\alpha },$$

(5)

$${\mu }_{\mathrm{VI}}\left(B_{pk}\right)=\frac{1}{\sqrt{\frac{1}{\pi }{{\displaystyle \int }}_0^{\pi }{\nu }_{DC}^2\left(B_{pk}\sin \alpha \right)d\alpha }}.$$

(6)

Each of those definitions assumes that the magnetic saturation coefficient of the ferromagnetic material is a function of the fundamental harmonic of the magnetic field only, expressed by B or H. The definitions (1), (3), (5) and (6) are taken from the literature, whereas (2) and (4) are the equivalents of (5) and (1), respectively.

2.2. Analyzed Machine

In this work, the effect of how the effective permeability is taken into account was examined on the example of a high-speed induction machine with a uniform solid rotor analyzed previously by Garbiec et al. The basic data of the considered motor is presented in Figure 1 and

Table 1. Basic specifications of tested induction motor.

Parameter	Value
Operation frequency (ω)	2000π rad/s
Number of pole pairs (p)	2
RMS phase voltage	50 V
Phase resistance	0.62 Ω
End-winding leakage	98 μH
Rotor conductivity	5.2 MS/m
Machine length (lz)	32 mm
Number of stator slots	24

.

A field-circuit multi-harmonic model with a strong coupling was developed for the analyzed machine. Based on the results of the analysis presented in by Garbiec et al., it was decided to reduce the number of the rotor models to three, i.e., those related to the fundamental harmonic of the magnetic field in the air gap and two most influential higher slot harmonics with ordinal numbers—11 and 13. Neglecting the rotor-end effects the system of equations describing the model takes the form [9]:

$$\left[\begin{array}{ccc}{M}_{11}\left({\mu }_{eff}\right)& M_{12}& M_{13}\\ M_{21}& M_{22}& 0\\ M_{31}& 0& 0\end{array}\right]\left[\begin{array}{c}\underset{\_}{\mathsf{\phi }}\\ \underset{\_}{I}\\ \underset{\_}{\mathsf{\lambda }}\end{array}\right]=\left[\begin{array}{c}0\\ \underset{\_}{E}\\ 0\end{array}\right],$$

(7)

where: $M_{11}$—matrix describing the magnetic and electrical properties of materials, $M_{12}=-M_{21}^T/\left(\mathrm{j}\omega l_z\right)$—matrix describing the distribution and connection method of the stator winding, $M_{13}=M_{31}^T$—matrices describing coupling between the rotor models and the stator model, $M_{22}$—stator winding impedance matrix, $\underset{\_}{\mathsf{\phi }}$—vector of the nodal values of the complex magnetic vector potential, $\underset{\_}{I}$—vector of the amplitudes of the loop currents in the stator winding, $\underset{\_}{\mathsf{\lambda }}$—vector of complex circulations of the magnetic field strength vector, $\underset{\_}{E}$—vector of the complex voltage amplitudes in the loops in the stator winding circuit.

For a comparison, in this paper, calculations were also carried out using a developed complementary time-domain model which uses the DC magnetization curves. The equations for such a model after discretization via Galerkin procedure and the implicit Euler method take the form:

$${\left[\begin{array}{cc}S\left({\mu }_{DC}\right)+G\Delta t^{-1}& -D^{\mathrm{T}}{K}^{\mathrm{T}}\\ l_{z}KD\Delta t^{-1}& K\left(R+L\Delta t^{-1}\right)K^{\mathrm{T}}\end{array}\right]}^n{\left[\begin{array}{c}\phi \\ i\end{array}\right]}^n={\left[\begin{array}{cc}G\Delta t^{-1}& 0\\ l_{z}KD\Delta t^{-1}& KL\Delta t^{-1}{K}^{\mathrm{T}}\end{array}\right]}^{n-1}{\left[\begin{array}{c}\phi \\ i\end{array}\right]}^{n-1}+{\left[\begin{array}{c}0\\ Ke\end{array}\right]}^{n-1},$$

(8)

where: S—reluctivity matrix, G—conductivity matrix, D—matrix describing the winding, K—matrix describing the winding connection method, R—winding resistance matrix, L—winding leakage inductivity matrix, ∆t—integration step, φ—vector of nodal values of the vector magnetic potential, i—vector of instantaneous values of the stator currents, e—vector of the instantaneous supply voltages values. The rotational movement was modeled using the moving band technique. The computing cost for the multi-harmonic model (when choosing an appropriate coupling scheme) is of the order of a per cent compared to that using the time-domain model [9]. This is the main reason for developing the model described by the Equation (7), however, in this paper the authors focus on the accuracy rather than on the computing cost. The nonlinearity in both above formulated computational problems is solved iteratively by updating, respectively ${\mu }_{eff}$ and ${\mu }_{DC}$ from iteration to iteration unless the relative change of solution is below a prescribed value.

Based on the DC magnetization curves for the stator and rotor materials, the effective magnetization curves were determined according to (1)–(6) (Figure 2). For saturated material, each effective magnetization curve is noticeably characterized by higher magnetic flux density values compared to the characteristic for the direct current. At a low level of magnetic saturation, the effective curves ale very close to the DC curve. Furthermore, a relatively small difference is found for different curves determined assuming the sinusoidal variation of the magnetic flux density for the rotor material. Moreover, the curves calculated with the assumption of a sinusoidal variation of the magnetic flux density are closer to the DC curve than the case with the assumption of a sinusoidal variation of the magnetic field strength. Additionally, curves (4)–(6) are closer to each other than (1)–(3). This is due to the completely different shape of the magnetic flux density waveforms (magnetic field strength) in the case when the magnetic field strength (magnetic flux density) is sinusoidal.

3. Calculation Results

The magnetization curves shown in Figure 2 were used to determine three basic operational characteristics of the machine (electromagnetic torque, RMS phase current, power losses in the rotor) vs. rotational speed at two different supply voltages (phase voltages), i.e., 50 V (rated) and 75 V (150% of the rated value) without changing the frequency of power supply. Note that the effective magnetic permeability for the rotor models is determined in each iteration based on the knowledge of the magnetic field strength distribution (expressions (1)–(3)) or magnetic flux density ((4)–(6)) only for the rotor model associated with the fundamental harmonic of the magnetic field distribution. The characteristics calculated were compared with the results of calculations using time-domain model, as described in Section 2. It should be kept in mind that the considered multi-harmonic model, despite the associated mathematical burden, is a simplified approach. This simplification can be observed in e.g., the obtained steady-state distributions of the magnetic flux over the machine cross-section shown in (Figure 3). Based on the results of that comparison, two groups of characteristics were created (Figure 4) showing the total relative percentage calculation error for the three determined quantities (electromagnetic torque, RMS phase current, power loss in the rotor) in relation to the calculation results obtained from the model (8):

$$\epsilon =\frac{{\epsilon }_{T+}{\epsilon }_{I+}{\epsilon }_P }{3}\times 100\%,$$

(9)

where: ${\epsilon }_T=\left|\frac{T_t-T_p}{T_t}\right|$, ${\epsilon }_I=\left|\frac{I_t-I_p}{I_t}\right|$, ${\epsilon }_P=\left|\frac{P_t-P_p}{P_t}\right|$, $T_t,I_t$, $P_t$—respectively, electromagnetic torque, stator RMS phase current, power loss in the rotor calculated with the model formulated in the time domain, $T_p,I_p,P_p$—respectively, electromagnetic torque, stator RMS phase current, and power loss in the rotor calculated using the multi-harmonic model.

Table 2. Averages percent value of the relative error for N = 11 considered points of the curves.

Parameter	μI	μII	μIII	μIV	μV	μVI
$\frac{{{\displaystyle \sum }}_{i=1}^N{\epsilon }_i}{N}$ for 50 V	7.75	12.57	10.54	8.72	8.03	8.64
$\frac{{{\displaystyle \sum }}_{i=1}^N{\epsilon }_i}{N}$ for 75 V	9.75	15.17	12.18	10.46	10.04	10.49

shows the averaged values of that indicator against the entire speed range considered.

Out of the definitions of the effective magnetic permeability considered, noticeably the best results were achieved using μ_I (originally used in calculations). Slightly worse results (except for the locked rotor state) can be obtained with the definitions assuming the sinusoidal variation of the magnetic flux density. On the other hand, an interesting conclusion is that the method how the nonlinearities are taken into account has practically no impact on the results of calculations for near-synchronous rotational speeds. It is caused by a much lower saturation of the rotor surface, crucial from the point of view of the machine properties, and consequently use by the time-harmonic model (7) similar values of effective permeability (see Figure 2). An increased error in the range of nominal operation point (21–27 krpm) is an adverse trend observed. However, note that the waveforms of the global error indicator shown in Figure 4 provide information on the discrepancy in the results of calculations of power losses in the rotor (${\epsilon }_P$). The results of those calculations are distinguished by the greatest discrepancy (for μ_I and N = 11 considered points of the curve $\frac{{{\displaystyle \sum }}_{i=1}^N{\epsilon }_{Pi}}{N}$ = 15.42% at 50 V and 19.32% at 75 V), thus overestimate the global error indicator. In the rotational speed range from 21,000 rpm to 27,000 rpm, the average electromagnetic torque error and the RMS current error is 7.56% and 1.43% at 50 V, 10.93% and 1.67% at 75 V, respectively.

Calculating the effective magnetic permeability based on the distribution of the magnetic field strength or magnetic flux density, determined only for the model associated with the fundamental field harmonic in the machine air gap, can raise doubts. To take into account the impact of the higher harmonics and to assess the error related to this assumption, a new definition of the multidimensional effective magnetic permeability was proposed:

$${\mu }_{VII}\left(H_{1pk},H_{11pk},H_{13pk}\right)=\frac{\sqrt{B_{1pk}^2+B_{11pk}^2+B_{13pk}^2}}{\sqrt{H_{1pk}^2+H_{11pk}^2+H_{13pk}^2}},$$

(10)

where $H_{1pk}$, $H_{11pk}$, $H_{13pk}$ are the magnetic field strength amplitudes related to the fundamental, eleventh and thirteenth harmonics of the field and

$$\begin{gathered} B_{npk}=\frac{2}{\pi }{{\displaystyle \int }}_0^{\pi }{\mu }_{DC}\left(H_{1pk}\sin \alpha +H_{11pk}\sin 11\alpha +H_{13pk}\sin 13\alpha \right)\left(H_{1pk}\sin \alpha +H_{11pk}\sin 11\alpha + \\ +H_{13pk}\sin 13\alpha \right)\sin n\alpha d\alpha \mathrm{for}n=1,11,13 \end{gathered}$$

(11)

For zero values of the higher harmonics, the definition (10) is equivalent to (1). This approach only slightly complicates the model implementation (in each algorithm iteration the distribution of the magnetic field strength in all considered rotor models needs to be obtained). An advantage of this expression, however, is the fact that effective permeability can be expressed before starting the calculations, without knowing the individual amplitudes. For the expression (10), it is possible since $B_{npk}$ does not depend on the phase angles of individual harmonics. In general, this occurs when the frequencies of the higher harmonics are the multiples of the fundamental harmonic. The value of the multidimensional effective permeability is determined only once as a multidimensional table for the assumed variability ranges of individual harmonics of the magnetic field strength. Figure 5 shows the determined distribution of this quantity assuming that $H_{13pk}$ = 0 A/m. The multidimensional permeability determined was used again for calculations and compared with the results obtained for (1). The result of that comparison is shown in

Table 3. Averages percent value of the relative error for N = 11 considered points of the curves.

Parameter	μI	μVII
$\frac{{{\displaystyle \sum }}_{i=1}^N{\epsilon }_i}{N}$ for 50 V	7.75	7.62
$\frac{{{\displaystyle \sum }}_{i=1}^N{\epsilon }_i}{N}$ for 75 V	9.75	9.59

. As can be seen, the application of the multidimensional effective permeability brings only slight improvement to the accuracy of calculation results. This is due to relatively low magnetic field strength values in the areas of the models associated with the slot harmonics, as illustrated in Figure 6.

4. Physical Validation

The analysis demonstrated that, from a practical point of view, it is sufficient to assume the first definition of the effective magnetic permeability, both due to the accuracy of calculations, and ease of implementation. It should be mentioned that, when adopting the sinusoidal variation of the magnetic field strength in determining the effective magnetic permeability, the fixed point method is applied in the multi-harmonic model to solve the nonlinearity problem, which is easy to implement, and it is distinguished by good convergence quality in case of time-harmonic magnetic fields. In the latter case, it is recommended to use the Newton–Raphson method [16].

The model developed, together with the first definition of the effective magnetic permeability, was applied to calculate the RMS current of a real machine (supplied via a quasi-square wave voltage inverter, fundamental frequency of supplying voltage equal to 500 Hz) with three different solid rotors tested by Garbiec and to compare the results with the calculations performed using a model with a weak coupling the same definition of the effective magnetic permeability applied (see Figure 7).

As it can be seen in the figure the idea of using the model with the strong coupling is to simultaneously solve the problem composed of a few submodels (the stator model and rotor models associated with the considered harmonics) coupled by appropriate boundary conditions imposed on the interface boundary in the air gap. This formulation results in one system of equations and determination of the magnetic vector potential distribution for all models via single solution of (7). The interface conditions in the air-gap are formulated with the use of matrix operators of the Fourier transform and the inverse Fourier transform, as shown in Garbiec at. al, enforce the propagation of only the selected harmonics of the air-gap field distribution between the stator model and the considered rotor models. All rotor models act via the same constraint on the stator model. In the considered model with the weak coupling, the monoharmonic field-circuit model is used. The power loss due to slot harmonics is estimated by means of extraction of these harmonics from the calculated distribution of the magnetic vector potential. In this way it is possible to solve the field problem limited to only the area of rotor. In this approach the power factor, torque and stator current can be appropriately corrected, although there is no backward relationship, i.e., the air-gap field remains unaffected by these effects. To this end, it was also necessary to take into account the rotor-end effect [18].

The physical test-stand and result of the comparison of computed and measured data for the machine with different rotors is shown in Figure 8. With the application of the multi-harmonic model with the strong coupling and with the modified sub-model coupling scheme, it was possible to achieve practically identical calculation results, yet in a shorter time by approx. 30%. The errors found result mainly from supplying the physical model via a quasi-square wave voltage inverter (in numerical models, the pure sine-wave voltage supply with an RMS value equal to the RMS value of the measured voltage is assumed).

5. Conclusions

The comparison of various methods for calculating the effective magnetic permeability for a multi-harmonic model of a solid-rotor induction machine carried out in the study proved that, from a practical point of view, the best results are obtained using the basic definition, described by Equation (1). Determining the effective magnetic permeability as the ratio of the magnitudes of the fundamental harmonics of the magnetic flux density and the magnetic field strength calculated assuming the sinusoidal variation of the latter allows for an easy implementation of the algorithm which is important from practical point of view. Furthermore, it was demonstrated that, in the studied case, the most influential magnetic saturation effects come from the fundamental harmonic of the magnetic field which means that the slot harmonics of the air-gap magnetic flux density do not have to be computed in that process. For the same definition of the effective magnetic permeability and an appropriate modification of the coupling scheme, it is possible to obtain practically the same results as in the case of the multi-harmonic model with the weak coupling, but yet in a shorter time. Another step to take will be to conduct detailed research using the developed model to analyze the magnetic saturation effects in the squirrel cage induction machines with both three- and single-phase windings.