Power Transformer No-Load Loss Prediction with FEM Modeling and Building Factor Optimization

Estimation of power transformer no-load loss is a critical issue in the design of distribution transformers. Any deviation in estimation of the core losses during the design stage can lead to a financial penalty for the transformer manufacturer. In this paper an effective and novel method is proposed to determine all components of the iron core losses applying a combination of the empirical and numerical techniques. In this method at the first stage all computable components of the core losses are calculated, using Finite Element Method (FEM) modeling and analysis of the transformer iron core. This method takes into account magnetic sheets anisotropy, joint losses and stacking holes. Next, a Quadratic Programming (QP) optimization technique is employed to estimate the incomputable components of the core losses. This method provides a chance for improvement of the core loss estimation over the time when more measured data become available. The optimization process handles the singular deviations caused by different manufacturing machineries and labor during the transformer manufacturing and overhaul process. Therefore, application of this method enables different companies to obtain different results for the same designs and materials employed, using their historical data. Effectiveness of this method is verified by inspection of 54 full size distribution transformer measurement data.


Introduction
Estimation of power transformer core losses before manufacturing of the core is a vital issue during the design phase of transformers [1].The objective function of most transformer design programs is Total Owning Cost (TOC) which includes purchase cost, maintenance cost, load losses and core losses cost in the life cycle of transformer [2].Moreover, core losses and efficiency of built transformers should satisfy the minimum requirements according to the existing standards [3,4].On the one hand, if the core loss estimation is higher than the actual value, a higher TOC is the outcome and that specific design may be considered unacceptable in the design stage.On the other hand, lower prediction of core losses can result in financial penalties due to violations from existing standards [5].Therefore, accurate core loss estimation of distribution transformers facilitates manufacturing of realistic transformers with reduced capital cost and increased efficiency [2].Conventional core loss estimation methods can be classified into four main categories [2]:  Empirical methods  Analytical methods

 Artificial intelligence methods  Numerical methods
A method related to each of these categories has its own strength and weakness.Empirical methods work based on experimental measurements and estimation of building factor (BF). BF is defined as the ratio of measured steel core losses to the estimated losses based on the calculation of nominal steel core losses [5].The empirical techniques are recognized for their speed of computation and sophistication due to covering all parts of the core losses.As the building factor depends on several parameters such as the air gap, the overlap areas at joints and the size of stacking holes [5][6][7][8][9], the empirical methods require a large number of measurements [2].Also, due to the continuous evolution of technical characteristics of both the magnetic materials and the core design, the core loss measurements of these distribution transformers should be updated [2].
The core losses can be modeled with different equivalent magnetic circuits by using analytical methods [10].These are based on semi-empirical determination methods of various components of distribution transformer no-load losses (such as hysteresis losses, eddy-current losses, and excess losses) which are functions of power frequency and maximum flux density in the core.Noload losses are simulated by introduction of a resistance to the general equivalent circuit model of the transformer [11].Analytical methods, in addition to their simplicity, are accurate for study of inrush current, ferroresonance, harmonic generation, etc [10,11].However, these methods can not estimate the core losses accurately and usually commercial transformer design programs are employed which utilize numerical or empirical methods [1].
Artificial intelligence methods are often based on neural networks which are used to estimate the core losses as a function of core design parameters [2,12,13].Accuracy of these methods is mostly dependent on the accuracy of training of neural network sets [2].Despite the good performance of neural networks in predicting noload losses of assembled transformers, there are some cases where the estimation error is not acceptable after the completion of transformer construction [2].
Due to the fast rise in speed of computers, numerical methods have become more attractive in the last decades.These methods predict no-load losses by solving Maxwell Equations with numerical techniques such as Finite Element Methods (FEM) or Finite Difference (FD) [1,14,15].The main problems of these methods are high calculation time and inability to model all various components of core losses [5].
In this paper, an effective and novel method is introduced which covers all components of core losses by combination of empirical and numerical methods.All computable components of core losses are calculated by using FEM modeling of transformer core.In this modeling, the following practical conditions and structures are taken into account [1,14,16,17]:  Anisotropy and non-linearity of the magnetic core material;  A voltage boundary condition similar to the standard no-load test of [18];  Core losses in directions other than the rolling direction of the core material;  Joints and the extra core losses resulted from air gaps which caused the distortion of flux distribution;  Stacking holes; As mentioned before, all parts of building factor such as inter-laminar losses and losses due to mechanical stresses cannot be modeled by numerical methods [5].Also, average air gap length related to joints (G) as an important parameter of Additional Localized Joints Losses (ALJL) is often undetermined and differs for each transformer manufacturer.In this paper, a Quadratic Programming (QP) optimization has been introduced to overcome these problems.By using results of empirical methods and experiments of engineers, incomputable parts of BF and G are predicted.This method has the following advantages:  Achieving more accurate estimation;  Simple implementation;  Improvement of core loss prediction continuously over time by performing more experiments on newly built transformers;  Offering a (personalized) core loss estimation program for every distribution transformer factory based on its manufacturing technology; In the rest of this paper, FEM modeling and numerical core loss estimation methods are discussed in Section 2. Combination of numerical and empirical methods is described in Section 3. Section 4 studies QP optimization.In Section 5 the effectiveness of the proposed method is verified against the measurement results on different commercial distribution transformers.Finally summary and conclusions are expressed in Section 6.

FEM Modeling and Loss Calculation of Core
The problem main equation is magnetostatic Maxwell equation written as follows [1].
where B and r are magnetic flux density vector and reluctivity tensor respectively.Common model of the reluctivity tensor is [1]: where r xx , r yy and r zz are highly non-linear functions of B [1].Although the non-diagonal elements of reluctivity tensor are non-zero due to the crystalline structure of the core material [19], the upper assumption is good enough to model the anisotropy degree of the material to reduce the complexity of problem and computation time [1].In this paper, Gaussian curve fitting with 3 terms is used for modeling the non-linear behavior of relative permeability elements as in Equation (3).
where µ 0 and µ r are vacuum and relative permeability respectively.a i , b i and c i are parameters which are determined with curve fitting algorithm.Gaussian method has a powerful ability for modeling the extremely non-linear behavior of material permeability [14].For instance, Figure 1 shows the measured data and the fitted curve for M5 CRGO Steel.Solving three-dimensional (3-D) models causes several problems such as very long simulation time and needs great amounts of computer resources (RAM, hard disk) [15].Especially exact 3-D modeling of the core with its laminations becomes impossible.Considering these problems, engineers and transformer designers try to replace the 3-D core model with some two-dimensional (2-D) simulations [1,[14][15][16].In this paper 2-D modeling of the core is used which is justified in commercial programs.So r zz in Equation ( 2) is not of interest in this method [1].
For 2-D modeling of core losses, two different models are proposed in the published literature as follows [1,14,16,17,20,21]:  XY-plane model  Z-plane model These two models are discussed in the rest of this section.

XY-Plane Model
XY-plane Model of the transformer core includes yokes and limbs of the core.In this model the upper half of the core is modeled and Equation ( 1) is solved with FEM over the structure.The imposed voltages on the transformer coils have a sinusoidal form, under the standard test condition.Flux boundary conditions on the centre of limbs are used to guarantee the sinusoidal waveform for the applied voltage [1].Using Faraday's equation, the total flux in the cross section of the limbs is a given value at each instant of time as follows [1]: where V, λ, n and S are sinusoidal imposed voltage, flux linkage, normal vector of limb cross section, and area of limb cross section, respectively.The presence of core stacking holes can increase the core losses up to 3% [5].These holes distort the flux wave shape, causing higher harmonic content and higher losses [17].
In this simulation the holes are modeled to obtain more accurate results, while [17] has guaranteed that the results of modeling the stacking holes in XY-model are in agreement with the measured data.
FEM program solves 2-D model of XY-plane for every pocket and returns the magnetic flux density distribution in the XY-plane at any instant of time.In Figure 2 the magnetic flux density distribution at two different instants of time is shown.
To calculate the magnetic flux density variation with respect to time in any element of transformer core, similar to [1], 12 equal space instants of time are assumed.The XY-plane model is solved for these 12 instants with their specific flux boundary conditions.For instance, Figure 3 shows magnetic flux density waveforms in the rolling direction and the cross rolling direction at point (P) in Figure 2.

Z-Plane Modeling
Due to presence of air gaps at the joints, magnetic flux travels between the overlapping laminations of yokes and limbs.When these laminations are saturated, flux crosses the air gaps [1].This flux distortion causes additional localized losses.This additional loss cannot be evaluated from XY-plane models.Therefore Z-plane model is utilized to solve this problem taking into account air gap at the joints [1,16,20,21].In    The variation of magnetic flux density in the SSL and MSL joints along the line L, in Figure 4, is plotted and shown in Figure 5.As it is shown in Figure 5, the flux distortion takes place only in a distance of a few centimeters around the air gaps.Therefore the application of z-plane model is limited to this area and it is independent of the overall structure of the transformer core [16].
For several magnetic flux densities, 2-D structure of Figure 4 is solved and the extra localized joint loss is calculated.By integrating losses over the sheets and dividing by thickness of sheets, ALJL den (W/m 2 ) is obtained [14,17].By multiplying ALJL den and surface of joints, total additional localized joint losses can be obtained [14].

Core Loss Model for Non-Sinusoidal Waveforms
As shown in Figure 3, the instantaneous magnetic flux density waveforms in the elements of transformer core demonstrate a non-sinusoidal variation.Therefore, in this modeling the equations related to the sinusoidal waveforms cannot be used.In recent years, accurate equations are developed to model the core loss when a non-sinusoidal flux density waveform exists in the magnetic materials.In this paper Equation ( 5) is used to calculate the losses in each element of the transformer core as proposed by [22,23]: (5 where B m , a, σ, ρ and T represent maximum flux density, the sheet thickness, its conductivity, its mass density and the time period respectively.K h , K exc and n a,b,c are the constant coefficients of the magnetic material.These coefficients can be obtained using the nominal loss of the magnetic sheets in specific frequencies and specific magnetic flux densities.Using Equation ( 5), spatial distribution of core loss of XY-plane is shown in Figure 6.Stacking holes are ignored in Figure 6 for a better visualization.

Combination of Empirical and Numerical Methods
Empirical methods predict no-load loss based on building factor estimation. No-load loss and consequently BF consist of various parts.Reference [5] has categorized different parts of BF and effect of each part by numerous experimental researches as indicated in Table 1.factor [24].General equation which is proposed by empirical methods can be written as follows:

As shown in
where CL is predicted core losses.G, OL and N are air gap length, overlap distance and number of steps at joints respectively.SL is nominal loss of magnetic sheets.V yol , BF yol , V cl , BF cl , V Lj , BF Lj , V Tj and BF Tj are volumes and building factors of yokes and outer limbs, center limb, L-joints and T-joints respectively.Reference [5] has considered volume of joints as overlap volume plus 15 millimeters of adjacent volume of core.Some parts of different BFs shown in Table 1 can be calculated by numerical method introduced in Section 2. These computable parts are indicated by '+' sign in Table 1.By elimination of these parts and inserting of cal-culated losses through Equation ( 6), this equation can be rewritten as follows: where W is the total weight of core.XYL and ZL are losses computed form XY-plane model and Z-plane model.As shown in Table 1, by elimination of computable parts from Equation ( 6), other parts are common between yokes and limbs.Therefore, BF yl is introduced as building factor of all limbs and yokes.Signs of '%' indicate ratio of the part volume to the overall core volume.
As mentioned before, SL•W is nominal loss of transformer core (NL) that can be predicted by Epstein coefficient.Also, in common designs of power transformer core such as mitred 45˚, the ratio of %V Lj to %V Tj is constant.Therefore Equation ( 7) can be rewritten finally by these simplifications as follows: where %V yl and BF yl are percent volume and building factor of yokes and limbs.%V j and BF j are percent volume and the combination of building factors of T-joints and L-joints.
Variables of Equation ( 8) are G, BF yl and BF j .These variables should be set by experimental data obtained from manufactured cores.Other parameters of Equation ( 8) can be either calculated in Section 2 or obtained using geometrical data of core design.

Building Factor Optimization
As shown in Section 3, after calculation of XY-plane losses and Z-plane losses with numerical methods, an optimization for selection of building factor (BF) and air gap length (G) is needed.First of all, relation between ZL and G should be determined.Figure 7 shows variation of ZL/NL versus G for various transformers for two example cores and a measured data from [25].It can be seen in [25] that there is a linear relationship between ZL/NL and G in the desired range where G varies from 1 to 3 millimeters.This is the result of an experimental measurement.Simulation results indicate the same thing (Figure 7).Therefore, ZL/NL can be written as follows: where a and b are constant variables obtained by curve fitting.As illustrated in Figure 7, these parameters differ for various cores.Using Equations ( 9) and ( 8) can be rewritten as follows: Assume the method has been provided by n number of transformer core design and measured no-load loss of them.The objective function of optimization is defined as summation of square of relative errors: where M i is the measured core loss of ith transformer.According to Equation (10), the objective function can be written as follows: where: By some mathematical works, Equation ( 12) can be written as a standard Quadratic Programming equation as follows [26]: . where Error is average of relative errors and var(Error) indicates variance of them.
By solving QP [26], the best possible solution for building factors is obtained.Also, a good estimation for air gap average length which is inevitable in manufacturing of core is attained.These parameters minimize average and standard deviation of prediction error as much as possible.As shown in Figure 8, by using XYL as no-load loss prediction, probability distribution of estimation has a bias error and undesired standard deviation.By introducing losses due to joints as ZL in prediction, average error reduces slightly but standard deviation does not change considerably.Finally by considering corrective building factors the core loss estimation becomes reasonable and average error and standard deviation reduce to an acceptable range.Also these optimized parameters can be used for core loss estimation of upcoming designs and can be improved every time required by new measured data.Therefore, the factory can use its personalized data to improve the no-load loss prediction without changing complex numerical program, only by a simple optimization using measured data.

Verification with Literature
Different parts of the proposed method can be verified against the published works in the related literature.For example, Figure 3 and Figure 6 demonstrate a very good similarity with the results in [1].Also Figure 4 and Figure 5 are verified with the results of [16] when CRGO steel is used.

Verification with Constructed Cores
Reference [27] reported complete data of a constructed core and measured core loss.This core is made of CRGO (grade M5) and by mitred 45˚ joints.The overlap length on joints is 5 mm, the core has one pocket by width of 70 mm, height and length of core are 300 mm and 350 mm.Core losses measured by a well-designed circuit and existing building factor over different flux densities is reported to have tolerance of 4%.As shown in Figure 9, estimation has a significant error by using only XY-plane losses for loss prediction.As mentioned in [27], air gap length (G) is tried to be less than 1 mm. Figure 9 indicates that even with assuming 1 mm air gap, estimated loss is less than the measured data.However, by introducing corrective building factor equal to 0.09, all pre-  dicted losses will come in the expected range.Measured data has an error of 4% which is recognized by [27] as measuring circuit error.This case can show that numerical methods separately cannot model transformer core loss, when real air gap length is used.Although by using greater air gap this weakness can be solved for numerical methods, but transformer design programs will use incorrect estimation about importance of air gap joints.

Verification with Full-Size Transformers
The authors had opportunity to use 54 full size transformers of one company for verification of the proposed method.These transformers are of 12 different designs.
For every type of design four or five manufactured transformers have been tested and the mean of their core losses has been nominated for its design.All transformers are made of cold rolled grain oriented (M5T30) material.The cores were three phase three limb stacked core with structures of mitred 45˚.Multi-Step lap joints have been designed for all cases.From these 12 different types, eight types have been used for training of method and computing of building factors and air gap length.Four types have been used for testing the proposed method named as the test group here.Optimization shows that the best value of building factors for joints is 37.01% and its best value for yokes and limbs is 0.21%.Also, the mean air gap length for the studied transformers of that company is 3 millimeters.By these parameters, the predicted core losses for the test group are shown in Table 2.As it can be seen in Table 2, this hybrid combination of numerical and empirical methods has an excellent ability to estimate core losses of transformers.

Summery and Conclusions
In this paper a hybrid method which is the combination of numerical and empirical methods has been used for accurate estimation of core losses of three phase stacked In the proposed method, first all of the computable parts of building factor are calculated by FEM models, and then other components of BF are optimized by quadratic programming.If more measured data become available, the core loss estimation is possible to be improved over time by using this method.Meanwhile, different companies can obtain different results for the same designs and materials employed, based on their historical data.
In the future work, the authors will use magnetizing current for estimation of the effective air gap length and will improve the method by introducing an accurate measured air gap length.

Figure 1 .
Figure 1.Gaussian curve fitting for modeling non-linear behavior of M5 CRGO Steel.

Figure 4 ,
the magnetic flux density behavior for Single Step Lap (SSL) and Multi-Step Lap (MSL) structures are shown using CRGO steel while operating under a magnetic flux density of 1.7 Tesla.

Figure 2 .
Figure 2. Magnetic flux density distribution at 2 instants of time.

Figure 3 .
Figure 3. Magnetic flux density waveform at point P indicated in Figure 2.

Figure 4 .
Figure 4. Magnetic flux distribution at area of joints due to presence of air gaps while operating under magnetic flux density of 1.7 Tesla in CRGO steel, a) Single Step Lap, b) Multi-Step Lap.

Figure 6 .
Figure 6.Spatial distribution of core loss of XY-plane model in upper half of 3-phase 3-limb core.

Figure 7 .
Figure 7. Percentage increases of power losses as a function of air gap length.

17 )
The last term in Equation (17) is a constant value and can be eliminated from the objective function.Finally, the QP problem is written as follows: the variance definition, the objective function can be assumed equal to Equation (19):

Figure 8 .
Figure 8. Probability distributions of estimated core losses by and without using building factors optimization.

Figure 9 .
Figure 9. Predicted core losses by and without using building factors optimization at different flux densities.

Table 1 ,
core volume includes different sections (yokes and outer limbs, center limb, corner joints, T-joints) and every section has special building