Fixed-Point Harmonic-Balanced Method for Nonlinear Eddy Current Problems

A new method to optimally determine the fixed-point reluctivity is presented to ensure the stable and fast convergence of harmonic solutions. Nonlinear system matrix is linearized by using the fixed-point technique, and harmonic solutions can be decoupled by the diagonal reluctivity matrix. The 1-D and 2-D non-linear eddy current problems under DC-biased magnetization are computed by the proposed method. The computational performance of the new algorithm proves the validity and efficiency of the new algorithm. The corresponding decomposed method is proposed to solve the nonlinear differential equation, in which harmonic solutions of magnetic field and exciting current are decoupled in harmonic domain.


Introduction
Non-linear eddy current problems can be solved by the time-stepping method [1] or the harmonic-balanced method [2]. The time-stepping method requires many periods to approach the accurate steady-state solution, while the harmonic-balanced method computes the magnetic field directly in the frequency domain. Compared with the frequency-domain method, the so called brute force method using time-stepping technique spends much more time on the transient process.
A relationship between magnetic field intensity H and the magnetic induction B is represented by introducing the fixed-point reluctivity ν FP [3]. The fixed-point reluctivity ν FP can be regarded as a periodic variable or a constant in harmonic-balanced method. The convergent speed of the solutions depends mainly on the strategies to determine the value of ν FP in non-linear iterations. Different methods to optimally determine ν FP have been presented and investigated in order to ensure the stable and fast convergence of solutions [4][5][6].
When the power transformer works under DC-biased magnetization, the ferromagnetic core will be significantly saturated. Owing to the nonlinearity of magnetic material in electromagnetic devices, there often are high-order harmonics in the exciting current and magnetic field. Therefore, electrical devices such as power transformers and reactors may work abnormally due to the magnetic storm and the transmission of high voltage direct current [7]. In that case, the DC-biased eddy current problem should be solved efficiently and accurately. Furthermore, electromagnetic coupling should always be considered in numerical computation when the strand coils or solid conductors are connected to the voltage source.
In this paper an efficient algorithm to determine the fixed-point reluctivity ν FP is proposed. It is aimed at efficiently computing the non-linear eddy current problem under DC-biased magnetization. The nonlinear magnetic field is computed in harmonic domain. The corresponding decomposed algorithm is presented to solve the nonlinear differential equation sequentially or concurrently, which decreases the memory cost of harmonic-balanced computation of large scale problems.

Fixed-Point Method
Maxwell's equations hold in non-linear eddy current problems,

H J E
where E is the electric field intensity, J 0 is the impressed current density and σ is the conductivity. A relationship between magnetic field intensity H and magnetic flux density B is represented by introducing the fixed-point reluctivity ν FP [8], where M is a magnetization-like quantity which varies nonlinearly with B. Therefore the magnetic field intensity H is split into two parts: the linear part, which is related to the ν FP, and the non-linear part, which varies with the magnetic induction B [3].

Harmonic-Balanced Method
The periodic variables in the electromagnetic field under DC-biased excitation can be approximated by the Fourier-series with a finite number of harmonics [9], where W(t) can be replaced by current density J, vector potential A, magnetic flux density B, magnetic field intensity H, and the magnetization-like quantity M.
Equation (1) can be rewritten in isotropic material by means of harmonic vector in 2-D problems, Each of harmonic vectors B x , B y , M x and M y has a similar expression with (7).

2-D Eddy Current Problems
Non-linear eddy current problems can be solved directly with the prescribed impressed current density. However, the impressed current density is unknown when the solid conductor or strand coil is connected to the voltage source. Therefore, the coupling between the magnetic field and electric circuits should be investigated if the non-linear eddy current problem is solved in the harmonic domain.

Solid Conductor Connected to Voltage Sources
When the solid conductor is fed by the voltage source, the eddy current exists in the solid conductor and the other conducting materials. The non-linear problem can be described by the following equation, where J is the current density. The magnetic vector potential A and the scalar potential V can be linked to the current density J by the equation as follows, The non-linear equation, including the 2-D magnetic and electric fields, can be presented as follows [10], The where Ω ew , Ω ed and Ω en represent the finite element in the whole problem domain, eddy current region and the non-linear region, respectively. N i is the shape function on node i in each finite element, and m is the total number of nodes in one element. D FP and N are the square matrices related to the fixed-point reluctivity and harmonic number [2], respectively. A j is the harmonic vector of the magnetic vector potential on node j and P j is the harmonic vector obtained from M. U is the voltage in harmonic domain per unit length.
By integrating (13) where Ω ec represents the finite element in the conducting region, R is the conductor's resistance per unit length, S cd is the cross-sectional area of the solid conductor.

Strand Coil Connected to Voltage Sources
The strand coil consists of fine wires where the eddy current is generally too small to be considered for computation. The supplied voltage U and the exciting current I in the coil can be linked by Kirchhoff's Law and Faraday's Law [11], where N coil is the turn number of the strand coil. The frequency-domain system equation considering electromagnetic coupling can be obtained according to the harmonic-balanced theory, where I u is a unit matrix of the same size with D FP and N. Consequently, the harmonic solutions of the magnetic field and magnetizing current can be computed simultaneously by solving (21) when the solid conductor and strand coil are both connected to the voltage source,

S + T G A K P C Z I U
where G is related to the spatial distribution of the magnetizing current when the strand coil is fed by the voltage source, while K appears on the right side of the equation when the solid conductor is connected to the voltage source.

Determination of Fixed-Point Reluctivity
The fixed-point reluctivity ν FP can be regarded as a periodic variable when it is determined in each time step. Consequently, the harmonic coefficients of ν FP can be used to calculate D FP in the harmonic-balanced method [2]. All elements in the square matrix D FP are non-zero, which indicates harmonic solutions are coupled with each other. In that case the memory demand will increase significantly in the large-scale computation, although fast convergence is achieved. In fact the ν FP can be a constant in the harmonic domain, and is determined as follows, where B max represents the maximum value of the magnetic induction in one period. The mean (e mean ) and maximum (e max ) variation of the reluctivity defined by ν=H/B can be observed to check the convergence of the harmonic solutions. In this paper the stopping criterions are set to e mean = 0.1% and e max = 1% in the numerical computation of the one-dimensional and two-dimensional eddy current problems.

Laminated Steel Sheet
As shown in Fig.1, A thin electrical steel sheet carrying eddy current is modelled and computed to observe variation of the magnetic induction when the lamination operates under different types of magnetization.
The 30Q140 oriented steel sheet of 0.3 mm thickness is first tested under sinusoidal flux in 50 Hz, and then the DC flux is provided. The conductivity of the sheet is σ = 2.22×10 6 S/m.  The spatial distribution of magnetic induction in the thin sheet is well clarified in Fig. 2. The magnetic flux for the boundary condition is B av,ac = 0.994T and B av,dc = 0.7036T. Fig.  3 compares the waveforms of magnetic induction under sinusoidal (indicated by "ac") and dc-biased (indicated by "dc") magnetizations. Notice that the eddy current in the sheet leads to the non-sinusoidal waveform of the magnetic induction. Furthermore, the distribution of the magnetic induction varies with the depth of the sheet in the x-direction.

Copper Conductor Surrounded by Ferromagnetic Screen
The two-dimensional problem consists of a solid copper conductor and an iron screen with an air gap. As shown in Fig.  4, the iron screen surrounds the conductor. The eddy current exists in both the copper conductor and iron screen. The conductivities of the copper and iron are σ= 5.7×10 7 S/m and σ = 1.0×10 6 S/m, respectively. The copper conductor is connected to a voltage source of 50 Hz. The B-H curve is detailed in [3]. 892 second-order elements with 2781 nodes are used in the numerical computation. Computational costs of the proposed method and the traditional method [2] are compared in Table I. M c and T c represent the memory demand and computational time, respectively. N h is the truncated harmonic number. Compared to the traditional method, the proposed method significantly reduces memory requirements with a slight increase in computational time due to a few more non-linear iterations. The calculated magnetizing current is compared with that obtained by using time-stepping method, and the good congruency proves the validity of the proposed method.

Decomposed Algorithm
When the fixed-point reluctivity is computed according to (24), the nonlinear equation in (23) can be linearized by ν FP which is space-dependent and time-independent. Therefore harmonic solutions can be decoupled and calculated separately. Equation (11) can be decomposed as follows, where k is the harmonic number. The decoupled equation system can be solved sequentially and concurrently, updating harmonic solutions by Gauss-Seidel and Jacobi iterative method respectively.

Conclusion
The convergent speed of harmonic solutions highly depends on the determination of optimal fixed point reluctivity in the fixed-point harmonic-balanced method. Differential reluctivity can be used to guarantee the stable and fast convergence of harmonic solutions. Due to the linearized system equation, harmonic solutions can be decoupled and computed in parallel, which can improve the computational efficiency and reduce the memory cost greatly. The proposed algorithm is more efficient than the traditional harmonic-balanced method.