基于定点谐波平衡法的铁心磁滞与损耗特性分析

赵小军 崔 灿 李 琳 赵志刚 刘 刚 李慧奇

（1.华北电力大学电力工程学院 保定 071003 2.华北电力大学电气与电子工程学院 北京 102206 3.河北工业大学电磁场与电器可靠性省部共建重点实验室 天津 300130）

1 引言

不同励磁条件下电工材料的磁性能并不相同。与交流励磁相比，当电工钢片或变压器铁心承受直流偏磁时，其磁化特性及损耗特性均会发生明显的改变[1,2]。在实际工程中，电力变压器经常承受直流偏磁或谐波形式的励磁，这会给变压器本身及电网的运行带来一系列危害[3,4]。各种不同励磁条件下电工材料磁特性的准确模拟对于电磁场理论研究和输变电设备的安全运行具有重要的意义。目前有以下相关研究工作:一方面，基于磁化特性、损耗特性的测量和模拟对电工钢片及铁心的磁性能进行研究，如等效磁路长度的测定[5]，铁心中不同区域的有功及无功分离[6]，各种磁滞模型的提出和改进[7-9]，开发不同的测量设备及设计相关的测量方案等[10,11]，这些工作取得的研究成果为面向工程应用的数值仿真奠定了良好的数据基础；另一方面，时步有限元，谐波平衡有限元和时间周期有限元等方法[12-14]被用于非线性磁场的仿真和分析，不同的数值计算方法和多种仿真技术可以为材料性能模拟的深入研究提供准确性及有效性方面的支持。

为了准确计算变压器铁心中的非线性磁场，深入分析不同励磁条件下铁心的磁化特性和损耗特性，需要在场路耦合计算中考虑铁心的磁滞效应[15]。解决该问题的关键在于磁滞模型的选择和如何对磁场中的非线性本构关系进行处理。选择磁滞模型时，应该遵循简单有效、数值计算时易实现的原则。对于含磁滞效应的非线性磁场本构关系的处理，目前主要有两种方法:①引入基于空气磁导率的μ0-B-H-M 关系[16]，②是引入基于定点磁阻率的νFP-B-H-M[17]关系。二者都能解决基于B-H 磁场本构关系中磁阻率ν 的不连续性问题。

本文利用叠片铁心模型分别进行了正弦励磁和直流偏磁励磁下的实验，得到相应的无偏磁磁滞回线和直流偏磁磁滞回线。利用基于损耗函数的磁滞模型对磁滞回线进行拟合和预测，并将该磁滞模型与定点谐波平衡有限元算法相结合，计算无偏磁和有偏磁条件下绕组励磁电流及铁心内的非线性磁场。比较了不同区域内的铁心磁滞特性，分析了不同励磁条件下铁心损耗的分布特征。

2 铁心磁特性实验

图1 叠片铁心模型Fig.1 Laminated core model

图1 所示为叠片铁心模型，铁心上绕有匝数相同的励磁线圈和测量线圈。基于叠片铁心进行空载实验，利用功率分析仪分别测量铁心损耗，励磁线圈的励磁电流i 和测量线圈中的感应电压u。叠片铁心所用硅钢片为30Q140 型取向硅钢片。

在无偏磁条件下，在励磁线圈端口施加不同的交流电压励磁，由下式可以得到相应的无偏磁磁滞回线。

式中，S 和L 分别为铁心截面积和等效磁路长度；φ为铁心磁通；Ncoil为线圈匝数。

本文中等效磁路长度L 取叠片铁心的几何平均磁路长度。

由无偏磁条件下的测量结果，得到不同励磁条件下的磁滞回线，如图2 所示。

图2 无偏磁磁滞回线Fig.2 Hysteresis loops under sinusoidal excitation

将一直流电流源与励磁线圈所在回路相串联，模拟变压器铁心的直流偏磁情况。此时铁心磁通φ包含两部分，即直流磁通φdc和交流磁通φac。φac可按照式（1）进行计算，φdc则可以通过迭代法进行计算[18,19]。

由直流偏磁条件下的测量结果，得到不同励磁条件下的磁滞回线，如图3 所示。

图3 直流偏磁磁滞回线（Idc=0.426A）Fig.3 Hysteresis loops under DC-biased excitation

3 磁滞模型与磁滞回线的拟合

3.1 基于损耗函数的磁滞模型

在正弦励磁条件下，如图4 所示，可以将励磁电流i 分为两部分，一部分im与磁通φ 同相位，一部分ih与感应电动势e 同相位[20]。于是，在i-φ 关系中可以引入以下函数

式中，i1对应于图4 中im；i2对应于图4 中ih。

图4 无偏磁条件下的励磁电流Fig.4 Magnetizing current under sinusoidal excitation

基于以上分析，对于如图5 所示的i-φ 磁滞回线，曲线“aOc”为中间磁化曲线（即i-φ 磁滞回线的中点轨迹），其所代表的励磁特性与式（4）中的i1相对应，i-φ 磁滞回线外围上任意一点与中间磁化曲线“aOc”的间距“ef”反映了铁心的磁滞效应且与损耗相关，与式（4）中的i2相对应，称之为损耗函数[21,22]。损耗函数的变化趋势与中间磁化曲线相反，当磁化曲线上升至最高点时，损耗函数值为零；当磁化曲线下降至零点时，损耗函数值则为最大，与“Ob”相对应。根据损耗函数的特点，可作如下定义

图5 基于损耗函数的磁滞模型Fig.5 Hysteresis model based on consuming function

式中，φm为磁通的幅值；D 是与损耗函数相关的系数。

损耗函数i2可进一步写为以下形式

式中，Iob为损耗系数，与磁通幅值φm及频率f 均相关，可以通过实验中的测量结果计算得到。

由式（2）、式（3）可知，对于叠片铁心，可以提出基于损耗函数的B-H 磁滞模型为

式中，Bm为无偏磁条件下交流磁通密度的幅值；H为铁心中总的磁场强度。

同理，在直流偏磁条件下，基于损耗函数的B-H磁滞模型如下

式中，Bdc为直流磁通密度，Bdc=φdc/S。

由于各中间磁化曲线均通过磁滞回线的顶点，其轨迹为基本磁化曲线，在有限元计算中可以分别选择无偏磁和直流偏磁下的基本磁化曲线作为式（8）中的H1(B) 和H1d(B)[19]。

3.2 磁滞模型的验证

通过实验，可以对测量结果进行处理，进而得到不同励磁条件下磁滞回线所对应的损耗系数Hob，结果见表1～表3。

表1 无偏磁条件下的损耗系数（Idc=0A）Tab.1 Consuming coefficients under sinusoidal excitation（Idc=0A）

表2 直流偏磁条件下的损耗系数（Idc=0.426A）Tab.2 Consuming coefficients under DC-biased excitation（Idc=0.426A）

表3 直流偏磁条件下的损耗系数（Idc=0.847A）Tab.3 Consuming coefficients under DC-biased excitation（Idc=0.847A）

由此可利用基于损耗函数的磁滞模型对测量得到磁滞回线进行仿真，无偏磁磁滞回线的结果如图6 所示，直流偏磁磁滞回线的结果如图7 所示。通过比较可以看出，基于损耗函数的磁滞模型能够较好的模拟无偏磁和直流偏磁条件下叠片铁心的磁滞效应。

图6 无偏磁磁滞回线的测量与仿真结果Fig.6 Simulated and measured hysteresis loop under sinusoidal excitation

图7 直流偏磁磁滞回线的测量与仿真结果Fig.7 Simulated and measured hysteresis loop under DC-biased excitation

4 计算结果及其分析

4.1 谐波平衡法

在无偏磁条件下，变压器端口承受稳态电压励磁，励磁电流中只含有奇次谐波，因此励磁电流i可以表达成如下形式

式中，I 为励磁电流密度i 的谐波矢量表达式。

同理可知在直流偏磁条件下，i 的谐波矢量表达式为

在稳态励磁下，电磁场中的各物理量均具有周期性，因此各变量均可表达为式（10）、式（11）所示的形式。

4.2 基于B-H-M 的定点谐波平衡方程

基于巴拿赫不动点定理，当不考虑各向异性时，可在二维非线性磁场中引入如下关系[17]

式中，νFP为定点磁阻率，谐波矢量Bx、By、Hx、Hy、Mx、My的表达式同式（10）、式（11）。

由式（12）可以看出，定点磁阻率νFP将磁场强度H 分为了线性和非线性两部分，线性部分与定点磁阻率相关，非线性部分与类磁化强度矢量M相关。

基于定点磁阻率的二维非线性磁场方程为

利用伽辽金法可以得到关于式（13）的定点谐波平衡有限元方程

式中，Q 为系数矩阵；Uk为线圈k 的端口电压谐波矢量；Zk为阻抗矩阵；Ik为流过线圈k 的励磁电流密度的谐波矢量；Ck为场路耦合矩阵；P 是与类磁化强度矢量M 相关的谐波矢量[23]。

4.3 νFP的取值与非线性场的求解

依据式（14）可对励磁电流和磁矢量位的谐波矢量同时进行求解，在求解过程中，谐波解的收敛性取决于定点磁阻率νFP的取值。在非线性静态场的计算中，定点磁阻率可按照下式确定[24]

式中，νdmax和νdmin分别为磁化曲线H=F(B)上微分磁阻率的最大值和最小值。

在时域有限元计算中，Dlala 等也给出了时步迭代中定点磁阻率的确定方案[25]。

在定点谐波平衡有限元计算中，初始迭代过程中采用较大的定点磁阻率（ν1=3ν0～ν0）可以保证谐波解的稳定收敛，之后采用较小的定点磁阻率（ν2=ν0/10～ν0/40）则能够实现谐波解的快速收敛。

图8 所示为考虑磁滞效应时的计算流程图，其中H=F(B) 与不同励磁条件下的磁滞回线相对应。忽略叠片铁心的各向异性，Bx、By、Hx、Hy满足以下关系

图8 计算流程Fig.8 Flow chart of the computational process

通过式（16），可以实现迭代计算中由B 求解H的过程。

5 计算结果及分析

5.1 励磁电流

对无偏磁和有偏磁条件下的励磁电流和磁场进行计算，励磁电流的计算结果如图9～图11 所示。从比较结果可以看出，计算结果与测量结果吻合较好。无偏磁条件下，由于磁滞效应的影响，励磁电流的前半周与后半周沿横轴（时间）不对称，但是沿纵轴（电流）呈反对称的特点。在直流偏磁条件下，由于励磁电流中的各次谐波分量迅速增大，励磁电流波形呈尖顶波形状，此时磁滞效应对于电流波形的影响不再明显，且电流波形也不再具有对称性。由图9 可以看出，在交流励磁较小，磁滞效应的影响相对较大时，励磁电流的计算结果与测量结果仍然存在一定的误差。这是因为在无偏磁条件下，交流励磁较小时叠片铁心的磁滞回线呈椭圆形，而本文提出的基于损耗函数的磁滞模型不能准确模拟椭圆形磁滞回线，因此数值计算结果的误差相对较大。要得到更加准确的结果，需要对该磁滞模型作进一步的改进和完善。

图9 无偏磁条件下的励磁电流Fig.9 Magnetizing current under sinusoidal excitation

图10 直流偏磁条件下的励磁电流（Idc=0.426A）Fig.10 Magnetizing current under DC-biased excitation（Idc=0.426A）

图11 直流偏磁条件下的励磁电流（Idc=0.847A）Fig.11 Magnetizing current under DC-biased excitation（Idc=0.847A）

5.2 铁心内不同区域磁滞特性

由计算得到的磁场，可以得到铁心内不同区域的磁滞特性。如图12 所示，为用于计算的叠片铁心几何模型，其中阴影部分代表线圈所在区域。根据铁心结构及磁场在铁心中的分布特点，可以将铁心分为两个区域，柱-轭区和接缝区。在两个区域中分别选择A 点和B 点，基于计算结果可以得到各点的磁场及对应的磁滞特性。

图12 叠片铁心的几何计算模型Fig.12 Geometric model of the laminated core for computation

由图13～图18 可以看出，在无偏磁条件下，叠片铁心内各点的磁滞回线是对称的，而直流偏磁条件下，叠片铁心内各点的磁滞回线呈现出明显不对称的特征。这是导致无偏磁和有偏磁条件下叠片铁心励磁特性和损耗特性存在差异的主要原因。由图14 和图15 的比较还可以看出，A 点处不同方向的磁滞特性存在差异，主要原因是在叠片铁心接缝区内的同一点处，磁通密度在不同方向上的分量并不相等。图17 和图18 的比较结果也反映了直流偏磁条件下接缝区内同一点处不同方向上磁滞特性的差异。

图13 无偏磁条件下的B 点y 方向磁滞回线Fig.13 Hysteresis loop on point B under sinusoidal excitation

图14 无偏磁条件下的A 点y 方向磁滞回线（Hy-By）Fig.14 Hysteresis loop on point A under sinusoidal excitation（Hy-By）

图15 无偏磁条件下的A 点x 方向磁滞回线（Hx-Bx）Fig.15 Hysteresis loop on point A under sinusoidal excitation（Hx-Bx）

图16 直流偏磁条件下的B 点y 方向磁滞回线（Idc=0.426A）Fig.16 Hysteresis loop on point B under DC-biased excitation（Idc=0.426A）

图17 直流偏磁条件下的A 点x 方向磁滞回线（Hx-Bx,Idc=0.426A）Fig.17 Hysteresis loop on point A under DC-biased excitation（Hx-Bx,Idc=0.426A）

图18 直流偏磁条件下的A 点y 方向磁滞回线（Hy-By,Idc=0.426A）Fig.18 Hysteresis loop on point A under DC-biased excitation（Hy-By,Idc=0.426A）

5.3 铁心的损耗特性

由计算得到每个单元的B 和H，可以按照下式计算铁心的损耗[26,27]

式中，Pi为第i 个单元中的比损耗；ρ 为铁心叠片的密度；T 为周期；Ns和ds分别为铁心中的叠片数量和单片厚度；Ne为单元总数；Si为单元面积；Pt为铁心损耗。

表4 中给出了无偏磁条件下铁心损耗的计算结果与测量结果，其中Pc为铁心比损耗的计算值，Pm为铁心比损耗的测量值，Bacm为铁心中的交流磁通密度的幅值，可由式（2）计算得到。由比较可以看出，计算结果与测量结果相吻合。

表4 无偏磁条件下的铁心损耗（Idc=0）Tab.4 The iron loss under sinusoidal excitation

由图12 中A、B 两点的磁场计算结果可以看出，铁心中的磁场分布并不均匀，尤其是接缝区和柱-轭区的磁场分布明显不同，由此可知，铁心损耗在不同区域的分布也将不同。分别计算铁心接缝区和柱-轭区的比损耗，观察不同励磁条件和励磁形式对二者的影响，结果见表5 和表6，其中Pj为接缝区比损耗，PL为柱-轭区比损耗。

表5 无偏磁下铁心不同区域比损耗（Idc=0）Tab.5 Specific loss in different areas in the laminated core under sinusoidal excitation

表6 偏磁下铁心不同区域比损耗（Idc=0.426A）Tab.6 Specific loss in different areas in the laminated core under dc-biased excitation

由计算结果可以看出，当铁心中交流磁通密度相同时，直流偏磁条件下的铁心损耗大于无偏磁条件下铁心损耗。同时，铁心内的损耗并不均匀，主要体现在接缝区与柱-轭区的差异。在相同励磁条件下，柱-轭区比损耗大于接缝区的比损耗。

6 结论

（1）基于叠片铁心的实验和计算表明，在无偏磁条件下和直流偏磁条件下，变压器铁心的磁滞特性和损耗特性明显不同，磁滞特性表现为偏磁后的磁滞回线的不规则与不对称性，损耗特性表现为偏磁后铁心损耗的增大。

（2）基于损耗函数的磁滞模型简单有效，能够模拟铁心在无偏磁和有偏磁条件下的磁滞效应。将其与定点技术和谐波平衡法相结合时，在数值计算中易实现，适用于频域有限元计算。

（3）在无偏磁和有偏磁条件下，铁心接缝区与柱-轭区磁场分布并不均匀，不同位置处的磁滞特性不同，并由此导致铁心损耗分布的不均匀性。在铁心接缝区内，同一位置处不同方向上的磁滞特性存在差异。

[1]赵志刚,刘福贵,张俊杰,等.直流偏磁条件下变压器励磁电流的实验与分析[J].电工技术学报,2010,25(4):71-76.Zhao Zhigang,Liu Fugui,Zhang Junjie,et al.Measurement and analysis of magnetizing current in DC-biased transformers[J].Transactions of China Electrotechnical Society,2010,25(4):71-76.

[2]郭满生,梅桂华,刘东升,等.直流偏磁条件下电力变压器铁心B-H 曲线及非对称励磁电流[J].电工技术学报,2009,24(5):46-52.Guo Mansheng,Mei Guihua,Liu Dongsheng,et al.B-H curve based on core and asymmetric magnetizing current in DC-Biased transformers[J].Transactions of China Electrotechnical Society,2009,24(5):46-52.

[3]段旭.直流偏磁下变压器振动机理与振动信号分析研究[D].重庆:重庆大学,2013.

[4]李泓志.直流接地极接地性能及入地电流对变压器影响的研究[D].北京:华北电力大学,2010.

[5]Cheng Z,Takahashi N,Forghani B,et al.Modeling of magnetic properties of GO electrical steel based on Epstein combination and loss data weighted processing[J].IEEE Transactions on Magnetics,2014,50(1):6300209.

[6]杜永,程志光,颜威利,等.电力变压器全斜接缝叠片铁心工作条件下的磁性能模拟[J].电工技术学报,2010,25(3):14-19.Du Yong,Cheng Zhiguang,Yan Weili,et al.Working magnetic property modelling of the power transformer laminated core with mitred joints[J].Transactions of China Electrotechnical Society,2010,25(3):14-19.

[7]曹林,何金良,张波,等.直流偏磁状态下电力变压器动态磁滞损耗模型及验证[J].中国电机工程学报,2008,28(24):141-146.Cao Lin,He Jinliang,Zhang Bo,et al.Dynamic hysteresis loss model of power transformer under DC current biasing and its verification[J].Proceedings of the CSEE,2008,28(24):141-146.

[8]张艳丽,刘洋,谢德馨,等.耦合改进型矢量磁滞模型的变压器磁场分析及实验研究[J].中国电机工程学报,2010,30(21):109-113.Zhang Yanli,Liu Yang,Xie Dexin,et al.Finite element analysis of magnetic field in transformer core coupled with improved vector hysteresis model and its experimental verification[J].Proceedings of the CSEE,2010,30(21):109-113.

[9]刘洋,张艳丽,谢德馨,等.考虑硅钢片矢量磁特性的复数 E&S 模型[J].中国电机工程学报,2012,32(3):144-149.Liu Yang,Zhang Yanli,Xie Dexin,et al.A complex E&S model considering 2D vector magnetic properties of silicon steel sheet[J].Proceedings of the CSEE,2012,32(3):144-149.

[10]Miyagi D,Yoshida T,Nakano M,et al.Development of measuring equipment of dc-biased magnetic properties using open-type single-sheet tester[J].IEEE Transactions on Magnetics,2006,42(10):2846-2848.

[11]Ishikawa S,Yanase S,Okazaki Y.AC magnetic properties of electrical steel sheet under twodimensional dc-biased magnetization[J].IEEE Transactions on Magnetics,2012,48(4):1413-1416.

[12]Yao Y,Koh C,Ni G.3-D nonlinear transient eddy current calculation of online power transformer under DC bias[J].IEEE Transactions on Magnetics,2005,41(5):1840-1843.

[13]Zhao X,Li L,Lu J,et al.Characteristics analysis of the square laminated core under dc-biased magnetization by the fixed-point harmonic-balanced FEM[J].IEEE Transactions on Magnetics,2012,48(2):747-750.

[14]Takahashi Y,Iwashita T,Nakashima H,et al.Parallel time-periodic finite-element method for steady-state analysis of rotating machines[J].IEEE Transactions on Magnetics,2012,48(2):1019-1022.

[15]刘福贵,杨庆新,颜威利,等.一种考虑磁滞特性的磁场数值计算方法[J].华北电力大学学报,2005,32(6):59-61.Liu Fugui,Yang Qingxin,Yan Weili,et al.A numerical method of magnetic field computation considering magnetic hysteresis property[J].Journal of North China Electric Power University,2005,32(6):59-61.

[16]Matsuo T,Shimasaki M.Time-periodic finite-element method for hysteretic eddy-current analysis[J].IEEE Transactions on Magnetics,2012,38(2):549-552.

[17]Mathekga M,Mcmahon R.Application of the fixed point method for solution in time stepping finite element analysis using the inverse vector Jiles-Atherton model[J].IEEE Transactions on Magnetics,2011,47(10):3048-3051.

[18]王祥珩,徐伯雄.变压器的偏磁问题[J].变压器,1992,29(8):11-14.Wang Xiangheng,Xu Boxiong.The problems on biased transformer[J].Transformer,1992,29(8):11-14.

[19]赵小军,李琳,程志光,等.基于直流偏磁实验的叠片铁心磁化特性分析[J].电工技术学报,2011,26(1):7-13.Zhao Xiaojun,Li Lin,Cheng Zhiguang,et al.Analysis of magnetizing characteristic of laminated core based on the dc-biasing experiment[J].Transactions of China Electrotechnical Society,2011,26(1):7-13.

[20]Kulkarni S V,Khaparde S A.Transformer engineering:design and practice[M].New York:Marcel Dekker,2004.

[21]Lin C E,Wei J B,Huang C L,et al.A new method for representation of hysteresis loops[J].IEEE Transactions on Power Delivery,1989,4(1):413-420.

[22]Faiz J,Sharifian M.Hysteresis loop modeling techniques and hysteresis loss estimation of soft magnetic materials[J].COMPEL,2001,20(4):981-1001.

[23]赵小军,李琳,程志光,等.定点谐波平衡有限元法与叠片铁心直流偏磁磁化特性研究[J].中国电机工程学报,2011,31(9):126-132.Zhao Xiaojun,Li Lin,Cheng Zhiguang,et al.Fixed-point harmonic-balanced finite element methodand dc-biasing magnetizing characteristics of laminated core[J].Proceedings of the CSEE,2011,31(9):126-132.

[24]Hantila F,Preda G,Vasiliu M.Polarization method for static field[J].IEEE Transactions on Magnetics,2000,36(4):672-675.

[25]Dlala E,Arkkio A.Analysis of the convergence of the fixed-point method for solving time-stepping nonlinear field problems[J].IEEE Transactions on Magnetics,2008,44(4):473-478.

[26]张艳丽,李玉梅,刘洋,等.考虑不同磁特性模型的感应电机铁心损耗分析[J].中国电机工程学报,2013,33(27):120-127.Zhang Yanli,Li Yumei,Liu Yang,et al.Analysis of core losses in induction motors considering different magetic property models[J].Proceedings of the CSEE,2013,33(27):120-127.

[27]Enokizono M,Shimoji H,Horibe T.Loss evaluation of induction motor by using magnetic hysteresis E&S2model[J].IEEE Transactions on Magnetics,2002,38(5):2379-2382.