Calculation and analysis of losses of magnetic-valve controllable reactor

ZHANG Hui-ying, TIAN Ming-xing, JING Pei, WANG Dong-dong

(1. School of Automation and Electrical Engineering, Lanzhou Jiaotong University, Lanzhou 730070, China；2. Rail Transit Electrical Automation Engineering Laboratory of Gansu Province,Lanzhou Jiaotong University, Lanzhou 730070， China)

Abstract： Magnetic-valve controllable reactor(MCR) has characteristics of DC bias and different types of magnetic flux density in the magnetic circuit and winding current distortion. These characteristics not only lead to loss calculation method of MCR different from that of power transformer, but also make it more difficult to calculate the core loss and wingding loss of MCR accurately. Our study combines core partition method with dynamic inverse J-A model to calculate the core loss of MCR. The winding loss coefficient of MCR is proposed, which takes into account the influence of harmonics and magnetic flux leakage on the winding loss of MCR. The result shows that the proposed core loss calculation method and winding loss coefficient are effective and correct for the loss calculation of MCR.

Key words： magnetic-valve controllable reactor (MCR)； dynamic inverse J-A model； core loss； core partition； winding loss coefficient

0 Introduction

The magnetic-valve controllable reactor (MCR) is a reactive power compensation equipment in power system, which has advantages of low cost, easy maintenance and smooth regulation of reactive power. Losses calculation of MCR is a very valuable research issue[1-2]. The structure and material of MCR are similar to those of transformers, and mainly composed of winding and core. But MCR also has the characteristics different from power transformer as follows: 1) AC and DC co-excitation; 2) only primary winding (working winding), without secondary winding; 3) the core column with winding has small section and magnetic saturation occurs in small section. These differences make the calculation method of core loss and empirical coefficient of winding loss of power transformer not applicable to the loss calculation of MCR.

MCR losses are mainly caused by core and windings, namely core loss and winding loss. Currently, the loss calculation method of power transformer is often used to calculate the loss of MCR[1,3]. Common methods for calculating the core loss of power transformer include: loss table method, Bertotti loss separation method, Steinmetz empirical formula and Finite element method[4-5]. The first three methods are applicable to the calculation of core material loss under sinusoidal excitation; the last method takes a long time because of large amount of computation. The loss empirical coefficient of power transformer winding is determined according to the design and engineering application experience of the transformer because it is related to winding structure, winding current characteristics and leakage magnetic field. Obviously, it is inaccurate to calculate the core loss and winding loss of MCR by using the core loss calculation method and winding loss coefficient of power transformer.

In this paper, the core partition method for calculating the core loss of MCR is proposed. Considering the influence of harmonics and magnetic flux leakage, the winding loss coefficient of MCR is deduced and analyzed, which is used to calculate the winding loss of MCR. At the same time, the result of the new method is compared with that of the method in Ref.[3], and the correctness of the proposed method is verified by the finite element analysis and experiments.

1 Loss calculation of MCR core

MCR is a shunt magnetic saturation reactor, which uses DC bias to adjust the magnetic saturation of the core columns, and changes the value of winding inductance to change the output reactive capacity. The structure diagram of a single-phase MCR is shown in Fig.1. The MCR core is made of ferromagnetic materials (silicon steel sheet), and is of an axisymmetric structure with four identical working windings. Turns of each windingN=N1+N2, tap ratioλ=N2/N；uAis supply voltage;ir1andir2are winding currents. The small section in the core column is called magnetic-valve. The ratio of the magnetic-valve section area to large section area of core column ismv. The rectifier circuit consisting of thyristorsVT1,VT2and diodeVDprovides DC bias current, that is, the control current of MCR. The control current is regulated by changing the trigger angle ofVT1andVT2. From no-load to full load, the corresponding trigger angle ofVT1andVT2varies from 180° to 0°.

Fig.1 Schematic diagram of a single-phase MCR

AC and DC hybrid excitation produces various types of magnetic flux densities in magnetic paths of MCR, such as AC flux density, DC flux density, AC and AC superposition flux density, AC and DC superposition flux density, etc. For accurate calculation of core loss, we propose a partition calculation method. The rules and steps of the core partition are as follows: Firstly, according to the type of magnetic flux density, the core is divided into several parts. Secondly, parts of the same type magnetic flux density are further subdivided into different parts based on the size of magnetic flux density. According to the above core partition rules, the core of MCR is divided into four partsS1-S4, as shown in Fig.1.

J-A model[6]is a hysteresis model based on the domain wall theory of ferromagnetic materials, and can well describe the hysteresis characteristics of ferromagnetic materials during magnetization. In order to describe the effect of eddy current, Baghel et al. improved the J-A model according to the energy conservation principle, and obtained the dynamic inverse J-A model as[7]

with the following complementary relationships

He=H+αM,

B=μ0(M+H),

Substituting the complementary relationships into Eq.(1), the detailed expression of the dynamic inverse J-A model is given by

whereManis anhysteretic magnetization,Mis total magnetization,μ0is the permeability of free space,δis a directional parameter,His magnetic field,Bis magnetic flux density, andδMis introduced to avoid nonphysical negative susceptibilities. According to the measured hysteresis loop of iron core material, it is very convenient to identify the parameters of dynamic inverse J-A model by using genetic algorithm or simulated annealing algorithm[8]. To simplify description, Eq.(2) and complementary relationships are replaced byfD(B,H)=0.

The dynamic inverse J-A model takesBas input and accords with the fact that MCR uses voltage source as excitation source. Therefore, it is a more suitable hysteresis model to describe the magnetization characteristics of ferromagnetic materials. The area ofB-Hhysteresis loop calculated by the dynamic inverse J-A model includes the hysteresis loss, eddy current loss and abnormal loss of core material, whereas the area ofB-Hhysteresis loop calculated by Preisach model or polynomial model only includes hysteresis loss[9].

In an alternating magnetic field, the sum of the hysteresis loss, eddy-current loss and abnormal loss in the unit mass core material is called the iron loss per unit mass. The difference of flux density types throughS1-S4leads to the different iron loss per unit mass, so the iron loss per unit mass ofS1-S4is calculated separately.

The loss per unit mass of core material is calculated by[7]

(3)

wherefis the frequency of an alternating magnetic,ρmis the density of core material. ∮cHdBrepresents the energy loss of unit volume in a cycle.

Hence, the core loss of MCR can be calculated by

PFe=m1pfe1+m2pfe2+m3pfe3+m4pfe4,

(4)

wherem1-m4are the quality ofS1-S4, respectively, obtained by design parameters;pfe1-pfe4are the loss per unit mass ofS1-S4, respectively.

In core of MCR, AC magnetic densityBac=UA/(2ωNSA)[10], whereUAis the amplitude ofuA,SAis the area of large section of core column, andω=2πf.

The relationship between thyristor trigger angleθand DC biasBdcis[11]

XarccosX)]-1,

(5)

The core loss of MCR is calculated by the proposed core partition method. Firstly,pfe1-pfe4are calculated, respectively.

1) The magnetic density inS1is equal toB1=mv(Bac+Bdc), theB-Hhysteresis loop is calculated by substitutingB1into Eq.(1), and then the loss per unit massPfe1ofS1is obtained by Eq.(3).

2) Settingke=ka=0, theB-Hhysteresis loop is calculated by substitutingBdcinto Eq.(1), and then the loss per unit massPfe2ofS2is obtained by Eq.(3).

3) The magnetic density inS3is equal toB2=Bac+Bdc, theB-Hhysteresis loop is calculated by substitutingB2into Eq.(1), and then the loss per unit massPfe3ofS3is obtained by Eq.(3).

4) TheB-Hhysteresis loop is calculated by substitutingBacinto Eq.(1), and then the loss per unit masspfe4ofS4is obtained by Eq.(3).

Secondly, substitutingm1-m4andpfe1-pfe4into Eq.(4),PFeis obtained.

2 Loss calculation of MCR winding

In Ref.[12], the equivalent magnetization model of MCR is proposed based on the single value magnetization curve. With the same idea, based on the dynamic inverse J-A model, the equivalent magnetization model of MCR is given by

(6)

whereBe,nis equivalent magnetic density;He,nis equivalent magnetic field;Bv,nis magnetic density in a magnetic valve;n=1 indicates column 1, andn=2 indicates column 2. It is the basis for calculating and analyzing MCR working current and winding current.

SubstitutingBe2andBe1into Eq.(6) to calculateHe1andHe2, respectively, according to Ampere loop law,He1andHe2are substituted into Eq.(7) to calculate winding currentir1andir2, respectively. Due to DC bias and magnetic saturation, the winding current contains a large number of harmonics. The frequency and amplitude of harmonics contained in winding current are obtained by fast Fourier transformation (FFT) method as

(7)

The winding loss of MCR including winding resistance loss and winding eddy current loss are calculated as follows.

2.1 Winding resistance loss coefficient

The loss of winding resistance is the basic loss of winding loss, and the loss of winding resistance will be increased by harmonics[13]. In order to consider the influence of harmonics on winding resistance loss under non-sinusoidal current, the harmonic loss correction coefficient of winding resistance loss is defined as

(8)

wherePDis the winding resistance loss at non-sinusoidal current with peak valueIP, andPDsis the winding resistance loss at sinusoidal current with peak valueIP.

According to the principle of harmonic balance, the non-sinusoidal winding current is decomposed into the sum ofkharmonic currents as

(9)

Hence, considering the influence of harmonics, the winding resistance loss under non-sinusoidal current can be expressed as

(10)

The resistance loss of a single winding under sinusoidal current is

(11)

Substituting Eqs.(10) and (11) into Eq.(8), we get the correction coefficient of winding resistance loss as

(12)

It can be seen from Eq.(12) that when the ratio ofI1toIPis constant, the value ofA1has nothing to do withIP.

2.2 Winding eddy loss coefficient

The winding eddy loss is related to the magnetic flux leakage intensity in the space at which the winding is located. There are two kinds of magnetic flux leakage that affect the winding eddy loss，axial magnetic flux leakage and radial magnetic flux leakage along the core limb. Because of unbalanced winding ampere-turn and magnetic saturation of magnetic value, the axial and radial magnetic leakages of MCR are much larger than that of transformer. Next, we study the winding eddy loss coefficient of MCR in detail.

In practice, the winding eddy loss coefficient is defined as the ratio of winding eddy loss to winding resistance loss, namely

(13)

(14)

The axial eddy current loss coefficient and radial eddy current loss coefficients areA2andA3, respectively;PAECis the wind eddy loss caused by axial magnetic flux leakage; andPRECis the winding eddy loss caused by radial magnetic flux leakage.

For a conductor in a uniform leakage magnetic field, its unit volume eddy loss can be calculated by[14]

(15)

whereBδis leakage flux density,fis winding current frequency,zis the thickness of the conductor perpendicular toBδandρis the resistivity of the conductor.

The four windings in Fig.1 have the same structural parameters and distribution of magnetic flux leakage. Therefore, we only analyze the eddy current loss of one winding. The longitudinal profile of winding and the distribution of axial magnetic flux leakage along the radial direction are shown in Fig.2.

The axial leakage flux density is triangular distributed along the radial direction of the winding, as shown in Fig.2. Therefore, the eddy losses of the unit volume windings caused by axial magnetic leakage along the radial of the winding are not equal. The axial leakage flux density at pointEis

(16)

Fig.2 Distribution of axial magnetic flux leakage along radial direction

Substitutingz=aandBδ=Byinto Eq.(15)，the unit volume eddy loss caused by axial magnetic flux leakage is given by

(17)

It can be seen from Eq.(17) thatpaceincreases with the increase ofy.

The eddy loss caused by axial leakage flux in one winding of Fig.1 is calculated by

(18)

whereVis the volume of one winding.

The axial eddy loss coefficient of MCR winding is obtained by substituting Eqs.(11) and (18) into Eq.(13) as

(19)

Note thatA2is independent of winding current and related to the winding current frequency and winding parameters such as resistance, turns and dimensions.

Due to magnetic saturation, the magnetic valve segment generates a large number of magnetic flux leakage. Because of the air gap edge effect, the distribution of radial magnetic flux leakage along the axial direction is as shown in Fig.3.

Fig.3 Distribution of radial magnetic flux leakage along axial direction

Note that whenxchanges in the range (-le/2,le/2), the corresponding radial magnetic flux leakageBxis symmetric about the coordinate origin. Therefore, we only analysis the winding eddy loss generated by radial magnetic flux leakage in the range (0,le/2). Supposeing that there is a pointGin the winding, the distance from pointGto the horizontal center line of magnetic-valve isx, to the vertical center line of the magnetic-valve isr, the radial flux density at pointGisBx, the magnetic flux density of the pointFof the magnetic-valve surface isB0, according to the flux continuity theorem, there is

2πr0B0dx=2πrBxdx,

that is

(20)

wherer0is the radius of the magnetic-valve.

Through circuit analysis the expression ofB0can be obtained by substituting Eq.(20) into Eq.(21) as

(21)

Then we get

(22)

Next, the expression ofBxcan be obtained by substituting Eq.(22) into Eq.(21) as

(23)

The values of parameterμ,N,h,r1andr0are constants. When max(x)=le/2, min(r)=r1-aand winding current is equal to the peak valueIp, the maximum valueBxmis got as

(24)

Whenxis in the range (-le/2,le/2), the radial leakage flux densityBxat any point in the winding is described by

(25)

Substitutingz=le/2 andBδ=Bxinto Eq.(15)，the unit volume eddy loss caused by radial magnetic flux leakage is given by

(26)

The eddy loss caused by radial leakage flux in one winding of MCR is calculated by

(27)

The radial eddy loss coefficient of MCR winding is obtained by substituting Eqs.(11) and (27) into Eq.(14) as

(28)

Note thatA3is independent of winding current and related to the winding current frequency, magnetic valve parameters and winding parameters such as resistance, turns and dimensions.

Hence, considering the influence of harmonics and magnetic flux leakage, the winding loss coefficient of MCR is

λW=A1+A2+A3.

(29)

Based on the above analysis, we get the loss sum of MCR four windings of MCR as

PCu=4(1+λW)PDs.

(30)

The losses of MCR is equal to the sum of the core loss and the winding loss namely

PLoss=PCu+PFe=

m1pfe1+m2pfe2+m3pfe3+m4pfe4.

(31)

The essential steps of MCR losses calculation are shown in Fig.4. By changingθ(the thyristor trigger angle), the losses of MCR under different working conditions can be obtained.

Fig.4 Flow chart of MCR losses calculation

3 Model validation and discussion

An example of MCR is used to verify the proposed core loss calculation method and the winding loss coefficient. In this example, MCR parameters are as follows:UN=220 V，SN=2.86 kVA，fN=50 Hz， and other parameters are detailed in Table 1.

Table 1 Parameter values of MCR prototype

The hysteresis loop data of core material are measured by using Epstein frame method[15]. Based on the hysteresis loop data from experimental measurement, the simulated annealing algorithm is used to identify the parameter values of the dynamic inverse J-A model. The parameter values of the dynamic inverse J-A model are as follows:Ms=1 566 538 A/m,α=0.7×10-5，a=5.5 A/m,k=35 A/m,c=0.06,ke=0.016 andka=0.45. Using the dynamic inverse J-A model, theB-Hhysteresis loops of core materials under different excitation conditions are obtained, as shown in Fig.5.

Fig.5 B-H hysteresis loops of core material under different excitations

It can be seen from Fig.5 that the shapes of hysteresis loops change from symmetry to asymmetry with the increase of bias when AC magnetic density remains unchanged, and the areas of the hysteresis loops also change, which indicates that the magnetization characteristic and unit volume loss of the core material change under different excitation conditions.

After obtaining the hysteresis loops of the core material under different flux densities, the core losses of MCR under different working conditions are calculated according to the calculation steps given in Fig.4. The core losses of MCR under no load (θ=180°), half load (θ=90°) and full load (θ=0°) was calculated by using the proposed method and the loss table method, respectively. The loss table method calculates the MCR core loss by looking up the unit mass loss value of the core material with different sinusoidal flux densities provided by the material manufacturer. Without the measuring device of MCR core loss, we built 3D finite element model of the MCR and simulated it to verify the correctness of the proposed calculation method. The core losses of MCR at differentθvalues are calculated by using different method, as shown in Table 2.

Table 2 Parameter values of MCR prototype

The results in Table 2 show that under the same load condition, the results of the proposed method are close to those of the finite element method, which illustrates that the proposed method for calculating core loss of MCR is correct and reliable. In addition, we can see that the results of the loss table method are larger than those of the other two methods.

Using Eqs.(6) and (7), we calculate the winding current of MCR, and then use FFT method to decompose the winding current to get harmonic components. The variation of the amplitude of each harmonic with the trigger angle is shown in Fig.6.

Fig.6 Harmonic amplitudes of winding current

The amplitude of the 9-order harmonic is about 1% of the fundamental amplitude. The amplitude of harmonic greater than 9-order harmonic is smaller than that of the 9-order harmonic. Therefore, we calculate the winding loss without the influence of harmonic greater than 9-order harmonic. Substituting the harmonic components into Eq.(12), the loss coefficientA1is obtained. The MCR parameter values are substituted into Eqs.(19) and (28), respectively, and the loss coefficientsA2andA3can be obtained.

The winding loss empirical coefficient of power transformer is used to calculate the winding loss of MCR[5]. In this paper, the winding loss of MCR is calculated by using the proposed winding loss coefficients and the empirical coefficient of power transformer, respectively, as shown in Fig.7.

Fig.7 Winding losses of MCR calculated by two methods

At the same trigger angle, the winding loss calculated by the proposed coefficient is larger than that calculated by the empirical coefficient. Atθ=0°, the winding loss by the proposed coefficient is about 15% larger than that of the empirical coefficient.

The core losses calculated by the loss table method and the winding loss calculated by the empirical coefficient are substituted into Eq.(31), respectively, and then we obtain the loss of MCR calculated by the method used in Ref.[5]. By using the same method, the loss of MCR calculated by the proposed method can be obtained. In Fig.8, the loss of MCR calculated by the proposed method and the method used in Ref.[5] are given, and the measurement results of MCR losses are also given. The result is normalized, and the normalized base value is the measurement value of MCR losses at rated operation.

Fig.8 Comparison between calculated value and measured value of MCR losses

It is clear from Fig.8 that the method proposed in this paper is effective and more accurate for calculating the losses of MCR. The proposed method does not take into account the additional losses generated by control loops, winding leads and structural parts, so the calculation result is smaller than the measured value. As the trigger angle decreases, the difference between the calculated and measured values becomes larger and larger, which indicates that the larger the output capacity of the reactor, the larger the proportion of additional losses without being taken into account in the MCR losses.

4 Conclusion

This paper presents a method for calculating the core loss of MCR based on the core partition and the dynamic inverse J-A model. The new method takes into account difference types of magnetic flux density in MCR magnetic circuit and the influence of DC bias. The rules of core partition are given in detail, so it is very convenient to apply this method to the core loss calculation of other complex structural core with mixed excitation.

Through analysis and deduction, the expression of winding loss coefficient of MCR is obtained. The new coefficient takes into account the effects of harmonics and magnetic flux leakage on the MCR winding loss, and varies with harmonic components.

The accuracy and effectiveness of the proposed method are verified by comparing the calculated result of MCR losses with the measured result. The result also shows that the core loss calculation method and winding loss empirical coefficient of power transformer are not suitable for the loss calculation of MCR.