Dynamic hysteresis modeling of silicon steel sheet considering excess eddy-current loss

Abstract

In this paper, we investigate modeling methods of magnetic properties considering the hysteretic properties and the excess eddy-current loss. First, a simple method in which an excess magnetic field given by the Bertotti’s model is added to a surface magnetic field as a post-processing of the 1-D finite-element analysis is compared with a 1-D finite-element method directly considering an excess magnetic field. Additionally, homogenization methods of magnetic property of an electrical steel sheet with the Bertotti’s model are discussed from the standpoints of computational accuracy and costs.

Keywords

Magnetic hysteresis excess eddy-current loss homogenization

1. Introduction

As electrical machines become more efficient and smaller, demands for improvement of the computational accuracy and costs of magnetic field analyses are increasing. Considering the complex properties of an electrical steel sheet including the hysteretic properties is necessary to improve the accuracy of a magnetic field analysis. Reference [1] proposes a fast and accurate post-processing method to estimate the iron loss combining the play model [2], which is one of the hysteresis models, with a 1-D finite-element method (FEM) [3].

To further improve the accuracy of the magnetic field analysis, it is essential to consider the excess eddy-current loss [4] in addition to the hysteretic properties and the skin effect. The Bertotti’s model [5] is widely used to represent the iron loss of an electrical steel sheet including the excess eddy-current loss [6]. The accuracy of the method in [1] is expected to be improved by combining it with the Bertotti’s model.

In this paper, we investigate two types of combinations of the Bertotti’s model with the 1-D FEM in [1]: one is a simple method in which the excess magnetic field is added to the surface magnetic field as a post-processing of the 1-D FEM, the other is the direct incorporation of the excess magnetic field in the governing equation of the 1-D FEM. Furthermore, we apply the former simple method to the Cauer circuit model [7–9] and a homogenization method [10] based on Legendre polynomials [11,12] with the Bertotti’s model and examine their accuracy.

2. Modeling method of magnetic property

2.1. Play model [2]

The scalar play model represents the hysteretic properties as follows: $\begin{eqnarray}\displaystyle p_{{\zeta}_{n}}(b)=b-\frac{(b-p_{{\zeta}_{n}}^{\ast })}{\max (|b-p_{{\zeta}_{n}}^{\ast }|,{\zeta}_{n})}{\zeta}_{n}, & & \displaystyle\end{eqnarray}$ (1) where b is the magnetic flux density, p _ζn is a play hysteron operator with a width of ζ_n, and $p_{{\zeta}n}^{\ast }$ is the operator at the previous time step. By using p _ζn, the magnetic field strength h is given by $\begin{eqnarray}\displaystyle h_{\text{DC}}=\mathop{\sum }_{n=0}^{N_{p}-1}f_{{\zeta}_{n}}(p_{{\zeta}_{n}}(b)), & & \displaystyle\end{eqnarray}$ (2) where N _P is the number of play hysterons, and f _ζn is a shape function for the operator p _ζn. The play model is mathematically equivalent to the Preisach model [13]. Therefore the same identification method for the Preisach model can be applied [14].

2.2. Bertotti’s model [5]

The Bertotti’s model represents the magnetic field strength h as follows: $\begin{eqnarray}\displaystyle h=h_{\text{st}}+h_{\text{ec}}+h_{\text{ex}}=h_{\text{DC}}(b)+\frac{{\sigma}d^{2}}{12}\frac{\text{d}b}{\text{d}t}+C_{\text{brt}}\left(\left|\frac{\text{d}b}{\text{d}t}\right|\right)^{1/2}\cdot \text{sign}\left(\frac{\text{d}b}{\text{d}t}\right), & & \displaystyle\end{eqnarray}$ (3) where σ is the conductivity, d is the thickness of an electrical steel sheet, and C _brt is the coefficient for the excess field. The three terms in the right-hand side correspond to the static hysteretic property h _st represented by the scalar play model h _DC, the classical eddy-current magnetic field h _ec, and the excess magnetic field h _ex, respectively. Because the classical and excess eddy-current losses corresponding to (3) are analytically calculated under sinusoidal flux excitation, the excess loss is separated from iron losses measured at two frequencies and therefore C _brt can be determined [15].

2.3. Cauer circuit model [7–9]

The magnetic field strength h in an electrical steel sheet can be represented as follows [10]: $\begin{eqnarray}\displaystyle \frac{\partial ^{2}h}{\partial z^{2}}={\sigma}\frac{\partial b}{\partial t}\quad \left(-\frac{d}{2}\leq z\leq \frac{d}{2}\right), & & \displaystyle\end{eqnarray}$ (4) where the z denotes the thickness direction.

In a linear theory, if the spatial distribution of the magnetic flux density is symmetric with respect to the center of an electrical steel sheet in the z direction, the magnetic flux density can be expanded by using Legendre polynomials P _2n(z) (n = 0,1, …) [7,11]: $\begin{eqnarray}\displaystyle b=b_{\text{a}}P_{0}+b_{2}P_{2}+b_{4}P_{4}+\cdots ={\mu}[h_{0}P_{0}+h_{2}P_{2}+h_{4}P_{4}+\cdots ∼]={\mu}h, & & \displaystyle\end{eqnarray}$ (5) where b _a is the average magnetic flux density and 𝜇 is the permeability.

Simultaneous equations obtained by substituting (5) into (4) can be represented by a ladder circuit with infinite stages shown in Fig. 1 [7], where h _s is the surface magnetic field strength, R = 4∕(σd ²), and L = 𝜇. The circuit is called the standard Cauer circuit. The magnetic field strengths h _2n and db _a∕dt correspond to the currents and the input voltage, respectively. The current passing through 3R corresponds to the eddy current generated by the temporal change of the main magnetic flux 𝜙₀. 𝜙₂ corresponds to the secondary flux generated by the eddy current. The current passing through 7R corresponds to the eddy current induced by the temporal change of 𝜙₂. In practical use, we need to truncate the standard Cauer circuit in finite stages. In this paper, we use the standard Cauer circuit with two stages and resistive termination as shown in Fig. 2, in which Legendre polynomials over the 6th orders are ignored. The relationship between Cauer circuit representation of eddy-current fields and its Legendre Expansion for a steel sheet is discussed in detail in [7].

Fig. 1.

Standard Cauer circuit.

Fig. 2.

Standard Cauer circuit (two stages, resistive termination).

Actually, the electrical steel sheet has a strong nonlinearity in its magnetic property. When the static magnetic property is represented as h _DC, the property of the first inductor L can be expressed as $\begin{eqnarray}\displaystyle h_{0}=h_{\text{DC}}({\Phi}_{0}). & & \displaystyle\end{eqnarray}$ (6) However, it is not straightforward to represent the property of the second inductor L∕5 because of the hysteretic properties. Therefore, the property of L∕5 is represented by using the following approximations [7,8]: $\begin{eqnarray}\displaystyle h_{2}\text{} & =\text{} & \displaystyle h_{\text{DC}}(5{\Phi}_{2}),\end{eqnarray}$ (7) $\begin{eqnarray}\displaystyle h_{2}\text{} & =\text{} & \displaystyle 5{\Phi}_{2}/{\mu}_{\text{d}},\end{eqnarray}$ (8) where 𝜇_d = [dh _app∕db]⁻¹ is the differential permeability and h _app is the approximated non-hysteretic curve of h _DC. As for h _app, we adopted the reversible component of h _DC given by the play model.

2.4. Method proposed by Gyselinck et al.

The method proposed by Gyselinck et al. [11,12] represents magnetic properties of an electrical steel sheet by using Legendre polynomials expansion of magnetic flux density like the Cauer circuit model. Therefore they are mathematically equivalent to each other in the linear theory. The equivalence of the truncated standard Cauer circuit to the method in [11] and [12] is proved in [7]. By substituting (5) into (4), the following equations are obtained: $\begin{eqnarray}\displaystyle h=h_{\text{s}}-{\sigma}d^{2}\left[{\beta}_{2}\frac{\partial b_{\text{a}}}{\partial t}+{\beta}_{4}\frac{\partial b_{2}}{\partial t}+\cdots \right], & & \displaystyle\end{eqnarray}$ (9) where b _a is the average magnetic flux density, d²𝛽_2n(z)∕dz ² = −P _2n−2(z)∕d ² (𝛽_2n(±d∕2) = 0), and b _2n are determined so that the following equation can be satisfied: $\begin{eqnarray}\displaystyle \frac{1}{d}\int _{-d/2}^{d/2}(h-h_{\text{DC}})P_{2n}∼\text{d}z=0. & & \displaystyle\end{eqnarray}$ (10)

3. 1-D finite-element method considering excess eddy-current loss

3.1. 1-D finite-element methods with Bertotti’s model

Two types of 1-D FEMs incorporating the Bertotti’s model are discussed. One is the direct incorporation of the excess magnetic field h _ex given by the Bertotti’s model in the governing equation of 1-D FEM [3] (Method A) as follows: $\begin{eqnarray}\displaystyle \frac{\partial (h_{\text{DC}}+h_{\text{ex}})}{\partial z}={\sigma}\frac{\partial A}{\partial t}, & & \displaystyle\end{eqnarray}$ (11) where A is the magnetic vector potential, h _dc is given by (2) and h _ex is defined in (3). The discretized form of (11) is represented as follows by using the Galerkin method: $\begin{eqnarray}\displaystyle G=\int \left[(h_{\text{DC}}+h_{\text{ex}})\frac{\partial N_{i}}{\partial z}+{\sigma}N_{i}\frac{\partial A}{\partial t}\right]\text{d}z=0, & & \displaystyle\end{eqnarray}$ (12) where G and N _i are the residual and the interpolation function, respectively. Method A can consider the skin effect regarding the excess magnetic field. Since h _ex is not differentiable everywhere, the fixed-point method is used as a nonlinear iteration method for Method A.

The other is a very simple method in which h _ex is added to the surface magnetic field h _FEM as a post-processing of the normal 1-D FEM (Method B). Method B has an advantage in computational costs because the Newton–Raphson method can be used as a nonlinear iteration method.

Table 1

Analysis condition

Parameter	Value
Thickness [mm]	0.5
Conductivity [S/m]	2.564 ×10⁶
Density [kg/m³]	7,700
C _brt	0.41
Frequency [Hz]	50 and 2 k

Table 2

Calculation time (B _m = 1. 5 T)

	Calculation time [s]
Method	50 Hz	2 kHz
A	83.1	78.2
B	2.22	2.92

Fig. 3.

Calculated C _brt.

3.2. Accuracy verification of 1-D finiteelement methods considering excess eddy-current loss

We carried out hysteretic magnetic field analyses of a ring specimen (IEC: M470-50A5) with Methods A and B. Table 1 shows the analysis conditions. The excess eddy-current loss coefficient C _brt is calculated from iron losses measured at 20 and 50 Hz. Although C _brt depends on the amplitude of magnetic flux density waveform B _m, we assumed that C _brt is constant and determined by averaging C _brt in the range of 0.8 T ≤ B _m ≤ 1.8 T in which the variation of C _brt is small as shown in Fig. 3. Because the modeling of excess loss based on constant C _brt has a room for improvement, we will investigate the better treatment of C _brt in future. The DC hysteresis loops required for identifying the play mode are estimated by subtracting the second and third terms in (3) from the measured result at 50 Hz. 40 symmetric loops (0.05 T ≤ B _m ≤ 2.00 T at intervals of 0.05 T) are used for the identification.

Figures 4 and 5 show the analysis and measurement results of magnetic properties of the specimen. Table 2 shows calculation time of each methods (B _m = 1. 5 T). Because Method A adopts the fixed-point method as a nonlinear iteration method, it needs much more calculation time than Method B. At 50 Hz, both methods give sufficiently good accuracy. At 2 kHz, however, both methods overestimate the magnetic field strength and the iron loss at high magnetic flux density although the condition that B _m is over 1.0 T at 2 kHz is not realistic in a practical electric machine. At low magnetic flux density less than about 0.5 T, both methods give accurate results. Therefore, they have sufficient accuracy for practical use. In addition, the difference between both methods is very small. From the standpoint of computational costs, Method B is a useful approach to consider the excess eddy-current loss.

Fig. 4.

Analysis and measurement results at 50 Hz.

Fig. 5.

Analysis and measurement results at 2 kHz.

4. Homogenization methods considering excess eddy-current loss

4.1. Homogenization methods with Bertotti’s model

The computational accuracy of the 1-D FEM decreases as the frequency and magnetic flux density increase. Moreover, the 1-D FEM requires a fine mesh to consider the skin effect accurately. Consequently, the 1-D FEM needs huge computational costs. To overcome this problem, the Cauer circuit model and the method by Gyselinck et al. have been investigated as an alternative of the 1-D FEM. In this paper, we introduce Bertotti’s model into these methods. In the previous section, we compared two approaches for applying Bertotti’s model to the 1-D FEM: one is incorporating h _ex in the formulation directly and the other is just adding h _ex as a post-processing. As a result, both methods give almost the same accuracy, although the computational cost of the latter approach is much lower. Therefore, we adopt the latter method (adding h _ex to surface magnetic field as a post-processing) for Cauer circuit model and the method by Gyselinck.

4.2. Accuracy verification of homogenization methods considering excess eddy-current loss

The analysis conditions are the same as the previous chapter. With regard to the Cauer circuit model, we used the circuit shown in Fig. 2. The representation methods of the property of L∕5 in (7) and (8) are investigated (Method I and Method II, respectively). As for the method by Gyselinck et al., we used up to the 4th order of Legendre polynomials (Method III).

Figures 6 and 7 show the analysis and measurement results of magnetic properties of the specimen under sinusoidal excitation. Table 3 shows calculation time of each methods (B _m = 1. 5 T). The analysis results of the 1-D FEM (Method B in the previous chapter) are also shown for comparison. At 50 Hz, all methods give the accurate results. At 2 kHz, however, each method has characteristic shapes of the loops and overestimates the magnetic fields. At high magnetic flux density, all methods give the iron loss more accurately than the 1-D FEM. In addition, at low magnetic flux density, Method I gives the most accurate results. Method I is based on the simple circuit and requires lower computational costs than Method III and the 1-D FEM. Therefore, by using Method I as a post-processing to estimate the iron loss, we can reduce the computational cost of the analysis while maintaining the accuracy.

Fig. 6.

Analysis and measurement results at 50 Hz.

Fig. 7.

Analysis and measurement results at 2 kHz.

Fig. 8.

Analysis and measurement results under non-sinusoidal excitation.

Fig. 9.

Analysis and measurement results under PWM excitation.

Table 3

Calculation time under sinusoidal excitation (B _m = 1. 5 T)

	Calculation time [s]
Method	50 Hz	2 kHz
I	0.24	0.21
II	0.18	0.18
III	0.56	0.54
1-D FEM	2.22	2.92

Table 4

Iron loss and calculation time under non-sinusoidal excitation

Method	Iron loss [W/kg]	ϵ [%]	Calc. time [s]
Measured	10.93	-	-
I	10.49	−4.04	0.22
II	10.87	−0.60	0.20
III	10.76	−1.57	0.53
1-D FEM	10.75	−1.68	2.78

Table 5

Iron loss and calculation time under PWM excitation

Method	Iron loss [W/kg]	ϵ [%]	Calc. time [s]
Measured	4.17	-	-
I	4.27	2.38	1.70
II	4.31	3.30	1.21
III	4.20	0.62	3.73
1-D FEM	4.20	0.65	31.5

Figure 8 shows the analysis and measurement results of magnetic properties of the specimen under harmonic excitation. The fundamental frequency is 100 Hz, and the order and amplitude of the harmonic are the 13th and 0.1 T, respectively. Table 4 shows the values and relative errors ϵ against the measurement results of the iron loss per period given by each method and calculation time. Methods II and III give more accurate results than 1-D FEM. Therefore, if the fundamental frequency is sufficiently low, the Cauer circuit model is effective as a post-processing to estimate the iron loss even if higher harmonics are superimposed to the excitation waveform.

Figure 9 shows the analysis and measurement results of magnetic properties of the specimen under PWM excitation. Table 5 shows the values and relative errors ϵ against the measurement results of the iron loss per period given by each method and calculation time. Method III and the 1-D FEM give the accurate iron losses, but their loops differ slightly from the measured loop. Methods I and II give the closest shape to the measured loop, but overestimates minor loops. Under PWM excitation, the method by Gyselinck et al. is more appropriate than the Cauer circuit model.

5. Conclusion

This paper investigates the modeling method of dynamic hysteresis incorporating the Bertotti’s model. First, we investigated two types of the 1-D FEM incorporating the Bertotti’s model. The results revealed that the accuracy of the 1-D FEM can be improved sufficiently by just adding the excess magnetic field to the surface magnetic field as a post-processing. Moreover, we applied the Bertotti’s model to the homogenization methods: the Cauer circut model and the method proposed by Gyselinck et al. It is revealed that the Cauer circuit model gives the most accurate iron loss under sinusoidal and harmonic excitations. The method proposed by Gyselinck et al. gives the most accurate iron loss under PWM excitation. Hence, it is expected that the computational costs can be reduced by using these methods as the post-processing to estimate the iron loss.

References

Takeda

Takahashi

Fujiwara

Ahagon

and Matsuo

, Iron loss estimation method for rotating machines taking account of hysteretic property, IEEE Trans. Magn. 51(3) (2015), 7300504.

Matsuo

and Shimasaki

, An identification method of play model with input-dependent shape function, IEEE Trans. Magn. 41(10) (2005), 3112–3114.

Yamazaki

and Fukushima

, Iron loss model for rotating machines using direct eddy current analysis in electrical steel sheets, IEEE Trans. Energy Conversion 25(3) (2010), 633-341.

Zirka

S.E.

Moroz

Y.I.

Marketos

and Moses

A.J.

, Dynamic hysteresis modeling, Physica B 343(1–4) (2003), 90–95.

Bertotti

, Hysteresis in Magnetism, Academic Press, 1998.

Yamazaki

and Fukushima

, Iron-loss modeling for rotating machines: Comparison between Bertotti’s three-term expression and 3-D eddy-current analysis, IEEE Trans. Magn. 42(8) (2010), 3121–3124.

Shindo

Miyazaki

and Matsuo

, Cauer circuit representation of the homogenized eddy-current field based on the Legendre expansion for a magnetic sheet, IEEE Trans. Magn. 52(3) (2016), 6300504.

Miyazaki

Mifune

Matsuo

Shindo

Takahashi

and Fujiwara

, Equivalent circuit modeling of dynamic hysteretic property of silicon steel under pulse width modulation excitation, J. Appl. Phys. 177 (2015), 17D110.

Hui

S.Y.R.

Zhu

J.G.

and Ramsden

V.S.

, A generalized dynamic circuit model of magnetic cores for low- and high-frequency applications —Part II: Circuit model formulation and implementation, IEEE Trans. Power Electronics 11(2) (1996), 251–259.

10.

Krah

J.H.

Bergqvist

and Engdahl

, Homogenization algorithms for calculation of field quantities in laminated magnetic materials, IEEE Trans. Magn. 38(5) (2002), 2385–2387.

11.

Gyselinck

and Dular

, A time-domain homogenization technique for laminated iron cores in 3-D finite-element models, IEEE Trans. Magn. 40(2) (2004), 856–859.

12.

Gyselinck

Sabariego

R.V.

and Dular

, A nonlinear time-domain homogenization technique for laminated iron cores in three-dimensional finite-element models, IEEE Trans. Magn. 42(4) (2006), 763–766.

13.

Mayergoyz

, Mathematical models of hysteresis, IEEE Trans. Magn. 22(5) (1986), 603–608.

14.

Bobbio

Miano

Serpico

and Visone

, Models of magnetic hysteresis based on play and stop hysterons, IEEE Trans. Magn. 33(6) (1997), 4417–4426.

15.

Suehiro

Mifune

Matsuo

Kitao

and Nakano

, Ladder circuit modeling of dynamic hysteretic property representing excess eddy-current loss, IEEE Trans. Magn. 54(3) (2018), 7300704.