基于磁场积分方程法的磁通门瞬态分析*

侯晓伟，李菊萍，郭俊杰

（宁波中车时代传感技术有限公司，浙江宁波315021）

基于磁场积分方程法的磁通门瞬态分析*

侯晓伟*，李菊萍，郭俊杰

（宁波中车时代传感技术有限公司，浙江宁波315021）

利用静态积分方程法来分析设计磁通门，虽然能够使计算时间大大减少，而且在精度方面也可以达到磁通门的要求，但是静态积分方程法没有准确考虑频率的影响，只能对磁通门进行准静态分析，而对磁通门进行瞬态分析的问题未能得到解决。本文根据磁通门铁芯特点，研究了基于磁场积分方程法的磁通门准静态分析模型，在此基础上，对其数学模型改进，建立了可以对磁通门进行瞬态分析的数学模型。最后以长条形磁通门为例验证了模型，仿真与实验结果符合良好。

磁通门；瞬态分析；积分方程法；准静态分析；铁芯；频率

EEACC:7230doi:10.3969/j.issn.1004-1699.2016.07.014

有限元方法作为一种成熟的方法，用其分析优化磁通门探头结构是一种常用的方法。Baschirotto A等［1-3］采用有限元软件分析设计了一系列的平面磁通门，并详细介绍了有限元软件分析磁通门的步骤以及剖分等关键步骤的处理。从中可以看出有限元方法需要大量的计算资源，计算时间长，例如文献［3］中仅铁芯的剖分单元就为288259，采用P4 2.4 GHz处理器进行一次磁通门瞬态分析所需的时间约为120个小时。如果想知道磁通门线性范围，最佳激励条件等还要进行多次瞬态分析，说明有限元法是非常费时的。尽管积分方程法剖分简单，待求解未知量少，但由于非对称满阵导致占用资源大，计算时间长，所以积分方程法一直被忽略，没有限元应用广。随着快速多极子方法在积分方程法中的应用，非对称满阵问题得到解决［4-6］。采用积分方程法处理电磁问题的研究也越来越多。2010年Vuillermet Y等［7-9］提出采用静态积分方程法分析磁通门，解决有限元法碰到的计算时间长问题。并将该方法与有限元法和传统采用退磁系数考虑铁芯形状的方法进行对比［7-8］，结果显示采用退磁系数的方法计算时间短，但在精度方面远不能与积分方程法和有限元法相比；与有限元法相比，静态积分方程法计算时间大大减少，而且在精度方面也可以达到磁通门的要求。静态积分方程法在采用退磁系数法和有限元法之间提供了一个很好的折中。但是Vuillermet Y在其论文［8］结尾也指出了一个问题，静态积分方程法只对磁通门进行了准静态分析。研究积分方程法瞬态分析磁通门是一个待解决的问题。文献［10］介绍的瞬态积分方程法虽然能够对电磁问题进行瞬态分析，然而这种方法涉及到磁场、磁化强度、涡流密度、磁矢量势、电标量势等未知量，待求解的未知数较多，求解复杂。目前还没有看到采用积分方程法解决磁通门瞬态问题的相关研究。针对此问题，本文根据磁通门铁芯特点［11-12］，研究了基于积分方程法的磁通门准静态分析模型，在此基础上，对其数学模型改进，建立了可以对磁通门进行瞬态分析的数学模型。最后以长条形磁通门为例进行了验证分析。

1　基于磁场积分方程的磁通门分析

考虑如下问题：磁通门铁芯置于外磁场（被测磁场和线圈电流产生的激励磁场）中。已知铁芯的形状参数和铁芯材料的磁化曲线，求解在外磁场中铁芯内部的磁场和磁化强度分布。这里考虑采用积分方程法来解决此问题。磁通门铁芯通常为薄片，考虑到激励磁场和被测磁场沿着铁芯表面，磁通门铁芯的厚度远小于其它尺寸，因此垂直铁芯表面方向的磁化强度可以忽略不计，可近似认为铁芯薄片内部的磁场主要平行于其表面［13］。另外磁通门接收线圈中感应到的是沿铁芯平面的磁场分量产生的磁通，而与之垂直的铁芯横截面内的磁场分量不会对接收线圈输出电压产生影响，因此磁通门铁芯的数学模型可以简化为平面二维模型。

1.1基于磁场积分方程的磁通门准静态分析

根据磁场叠加原理，对于上述问题，磁通门铁芯任意一点的总磁场可以分为三部分：

其中He和Hm为源场，He为激励电流产生的磁场，Hm为被测外磁场，Hd为铁芯材料产生的退磁场，旋度为零。

将磁通门铁芯剖分为n个单元，分析时的假设条件为：①每个单元内的磁化强度是均匀的；②磁荷只在单元之间的交界面以及铁芯与空气的边界出现。可得：

其中hi为第i个单元的总磁场强度，为激励电流在第i个单元产生的磁场，为被测外磁场在第i个单元产生的磁场，为铁芯材料在第i个单元产生的退磁场，mi为第i个单元的磁化强度，下标x和y分别表示在铁芯平面x轴和y轴方向的磁场分量。

1.1.1激励磁场计算

根据毕奥-萨伐定律，可以求得激励电流在场点r处产生的磁场强度为：

其中nJ为激励电流方向，J为激励电流密度，r为磁场的坐标，r′为电流的坐标，dV为激励电流体积微元。

1.1.2退磁场计算

由铁芯材料产生的退磁场所满足的方程为：

退磁场旋度为零，存在一个标量满足：Hd=-∇ϕ。而标量所满足的方程为：

铁芯材料所有剖分单元在第i个单元产生的磁场为：

其中ϕj为第j个单元对应的磁标势，r为源点和场点之间的距离，mj为第j个单元的磁化强度，sj为第j个单元的外表面，其中dsj=dsjn，n为第j个单元沿外表面法线单位矢量。将式（3）和式（6）代入式（2）可得：

退磁矩阵wij的计算详见文献［14］。联立式（7）与铁芯材料的磁化曲线，建立磁通门准静态分析模型，即可求解被测磁场和激励电流产生的激励磁场中铁芯内部的磁场和磁化强度分布，从而求得磁通门准静态输出信号。

1.2基于磁场积分方程的磁通门瞬态分析

1.1节所研究的方法由于没有考虑频率对电磁场计算结果的影响，仅限于对电磁场进行准静态分析，无法对磁通门进行瞬态分析。本节在1.1.1节的磁通门准静态分析方法基础上，研究引入频率的影响，建立磁通门瞬态分析模型，对磁通门输出信号进行瞬态分析。

在低频情况下可以忽略位移电流［10］，电流密度i可分解为涡流ie和给定的激励电流密度i0，即i=ie+io。根据磁场叠加原理，考虑涡流的影响后，磁通门铁芯内磁场强度Hs（与前面采用H表示铁芯内磁场不同，这里为方便起见，用Hs表示铁芯内磁场）可以分解成：

H0是源磁场，包括给定的电流产生的磁场强度和被测外磁场，Heddy是涡流产生的寄生磁场，Hd是铁芯材料产生的退磁场。

将涡流产生的磁场强度移项到左端，并记H=H0+Hd，可得 Hs=H+Heddy，铁芯内磁场方程组为：

其中铁芯材料数学模型M（H）为：

对于一些特殊形状的铁芯，如薄片状，圆柱状等，可以用铁芯形状参数和磁感应强度时变率表示涡流产生的寄生磁场［15-17］。

以分析长条形磁通门为例，忽略横向磁场影响，采用文献［18］介绍的铁芯材料的数学模型。在外磁场H中，铁芯材料内的磁场Hs与磁化强度M的关系由以下方程组描述：

以上研究了采用静态磁场积分方程怎样分析瞬态问题，瞬态磁场积分方程组与静态磁场积分方程组形式相同，频率的影响可以归结为调用考虑频率的铁芯材料数学模型。

2　磁通门瞬态分析模型验证分析

以长条形磁通门（图1所示）为例对1.2节提出的磁通门瞬态分析模型进行验证分析。激励线圈和接收线圈都采用均匀绕线的方式。激励线圈匝数为330，接收线圈匝数为1 100。铁芯材料是长度、宽度和厚度分别为30 mm、1 mm和0.1 mm的坡莫合金。铁芯材料模型由式（11）描述，其模型参数为 a=1.069 4×10-6，c=0.3，k=1.714 3，MS=5.98×105，Ru=0.399 7，H0=0.000 1，H1=1.1263×10-20，a0=0.082 7，a1=0.020 1，N=0.001 4。图2为正弦电压激励的长条形磁通门等效电路。其中铁芯材料由式（11）所示的模型描述。AC为交流激励电流源，R是值为13.3 Ω的激励线圈电阻，Hx为沿着铁芯表面x轴方向的被测磁场，Vout为磁通门的输出信号。

图1　长条铁芯磁通门

图2　长条铁芯磁通门等效电路

采用正弦电压信号对磁通门进行激励，比较了考虑三种不同频率条件下（1 kHz、3 kHz和5 kHz）的基于磁场积分方程的磁通门瞬态分析模型的仿真结果和实验结果。图3、图7和图11分别表示在三种不同频率和激励信号幅值条件下，激励线圈电流的仿真结果和实验结果比较。图4、图8和图12分别表示在三种不同频率和激励信号幅值条件下，接收线圈输出电压的仿真结果和实验结果比较。图5、图9和图13分别表示在三种不同激励频率条件下，磁通门铁芯的仿真磁滞回线和实验磁滞回线的比较。图6、图10和图14分别表示被测磁场为30 μT时及三种不同的频率条件下，随着激励电压增大，输出二次谐波电压幅值与激励线圈电流幅值的关系。

图3　激励线圈电流的仿真结果和实验结果比较

图4　接收线圈两端输出电压的仿真结果和实验结果比较

图5　激励频率为1 kHz时的磁通门铁芯的磁滞回线的仿真和实验比较

图6　激励频率为1 kHz时二次谐波幅值与激励线圈电流幅值的关系

图7　激励线圈电流的仿真结果和实验结果比较

图8　接收线圈两端输出电压的仿真结果和实验结果比较

图9　激励频率为3 kHz时的磁通门铁芯的磁滞回线的仿真和实验比较

图10　激励频率为3 kHz时二次谐波幅值与激励线圈电流幅值的关系

图11　激励线圈电流的仿真结果和实验结果比较

图12　接收线圈两端输出电压的仿真结果和实验结果比较

图13　激励频率为5 kHz时的磁通门铁芯的磁滞回线的仿真和实验比较

图14　激励频率为5 kHz时二次谐波幅值与激励线圈电流幅值的关系

由图3、图4、图7、图8、图11和图12可以看出，因为实际中的激励电流需要一定的响应时间，在开始的一段时间内，激励电流和输出信号的实验值滞后于仿真模型中的激励电流和输出信号。在后续的时间条件下，仿真中的激励电流和输出信号与实验中的激励电流和输出信号吻合良好。从图5、图9和图13可以看出，仿真磁滞回线与实验得到的磁滞回线几乎一致。从图6可以看出，在小于95 mA的激励电流条件下，由仿真模型得到磁通门二次谐波幅值略低于实验测得的磁通门二次谐波幅值；在大于95 mA的激励电流条件下，由仿真模型得到的磁通门二次谐波幅值高于实验测得的磁通门二次谐波幅值。而图10和图14中随着激励电流的增大，二次谐波的仿真值略低于实验值。表1给出了在被测磁场为30 μT时及三种不同的频率条件下，二次谐波最大值的仿真结果与实验结果的对比。

表1　二次谐波最大值的仿真结果与实验结果对比

从仿真结果与实验结果的对比可以看出，随着激励电压和激励频率的改变，仿真结果准确的反应了线圈电流、输出电压等参数的变化。表明磁通门瞬态分析模型准确地综合考虑了探头结构、铁芯材料、激励电压、激励频率等因素的影响，准确地考虑了随着频率增大而增强的磁感应强度滞后现象。如果采用文献［8-9］的方法，铁芯材料由磁化曲线描述，无法考虑频率引起的磁感应强度滞后现象，只能对磁通门进行准静态分析。

3　小结

采用静态积分方程法对磁通门准静态分析不能准确考虑频率的影响，瞬态积分方程法虽然能够对磁通门进行瞬态分析，然而这种方法中的未知量磁矢势为静态积分方程法未知量磁标势的3倍，使求解麻烦。针对这个问题，本文根据磁通门铁芯特点，分析了解决静磁问题的磁场积分方程，在此基础上考虑了涡流对积分方程的影响，对解决静磁问题的静态积分方程改进，建立了可以对磁通门进行瞬态分析的数学模型。最后以长条形磁通门为例进行模型验证，仿真了激励电压、激励频率改变时的接收线圈输出电压、输出二次谐波电压幅值与激励线圈电流幅值的关系和磁滞回线，这些仿真结果与实验结果符合良好。改进了的模型可以准确地反应频率的影响，并对磁通门进行瞬态分析。

［1］ Baschirotto A，Dallago E，Malcovati P，et al.An Integrated Micro-Fluxgate Magnetic Sensor with Front-End Circuitry［J］.IEEE Transactions on Instrumentation and Measurement，2009，58（9）：3269-3275.

［2］ Baschirotto A，Dallago E，Ferri M，et al.A 2D Micro-Fluxgate Earth Magnetic Field Measurement Systems with Fully Automated Acquisition Setup［J］.Measurement，2010，43（1）：46-53.

［3］ Baschirotto A，Dallago E，Malcovati P，et al.Fluxgate Magnetic Sensor in PCB Technology［C］//IEEE Instrumentantion on and Measurement Technology Conference Proceedings，2004（2）：808-812.

［4］ Mayergoyz I，Andrei P，Dimian M.Nonlinear Magnetostatic Calculations Based on Fast Multipole Method［J］.IEEE Transactions on Magnetics，2003，39（3）：1103-1106.

［5］ Balasubramanian S，Lalgudi S，Shanker B.Fast-Integral-Equation Scheme for Computing Magnetostatic Fields in Nonlinear Media［J］.IEEE Transactions on Magnetics，2002，38（5）：3426-3432.

［6］ Hafla W，Buchau A，Groh F，et al.Efficient Integral Equation Method for the Solution of 3-D Magnetostatic Problems［J］.IEEE Transactions on Magnetics，2005，41（5）：1408-1411.

［7］ Vuillermet Y.About the Use of Numerical Integral Methods to Simulate a Fluxgate Magnetometer：The Ring-Core Example［J］. Sensors&Actuators A：Physical，2010，163（1）：48-53.

［8］ Vuillermet Y，Audoin M，Cuchet R.et al.Application of a Non-linear Method of Moments to Predict Microfluxgates Output［J］.Sensors&Actuators A：Physical，2010，158（2）：212-216.

［9］ Vuillermet Y，Audoin M.Application of a Non-Linear Numerical Integral Method to Predict Micro Fluxgate Output［C］//14th Biennial IEEE Conference on Electromagnetic Field Computation（CEFC），2010：1-1.

［10］Koizumi M，Onizawa M.Computational Method of Three Dimensional Eddy Current by Using Volume Integral Equation Method［J］.IEEE Transactions on Magnetics，1991，27（5）：4077-4080.

［11］杨尚林，刘诗斌，郭博，等.铁芯参数对磁通门输入输出特性影响分析［J］.传感技术学报，2014，27（1）：40-47.

［12］郭博，刘诗斌，杨尚林，等.微型正交激励磁通门结构设计［J］.传感技术学报，2014，27（7）：910-915.

［13］Canova A，Repetto M.Integral Solution of Nonlinear Magnetostatic Field Problems［J］.IEEE Transactions on Magnetics，2001，37（3）：1070-1077.

［14］周国华，肖昌汉，刘胜道，等.基于六面体单元表面磁场积分法求解三维静磁场［J］.电工技术学报，2009（3）：1-7.

［15］Mandache L，Topan D.Managing Eddy Current Losses and Ferromagn/Etic Material Nonlinearities in Distorting Regimes［C］// Electric Machines and Drives Conference，2009.IEMDC'09. IEEE International，2009：1449-1454.

［16］Jiles D.Frequency Dependence of Hysteresis Curves in‘Non-Conducting'Magnetic Materials［J］.IEEE Transactions on Magnetics，1993，29（6）：3490-3492.

［17］Chikazumi S，Charap S H.Physics of Magnetism［M］.RE Krieger Publishing Company Malabar，FL：1978.

［18］Mandache L，Topan D.Managing Eddy Current Losses and Ferromagnetic Material Nonlinearities in Distorting Regimes［C］//IEEE International Electric Machines and Drives Conference，2009：1449-1454.

侯晓伟（1985-），河南洛阳人，工程师。研究方向为微电子与固体电子学，主要致力于微型磁通门传感器、专用集成电路 ASIC等方面的研究，xwhou2009@ 126.com；

李菊萍（1982-），山西运城人，工程师。研究方向为微电子与固体电子学，主要从事微型磁通门传感器等方面的仿真研究，li2008100096@mail.nwpu.edu.cn。

Transient Analysis of Fluxgate Based on Magnetic Field Integral Equation Method*

HOU Xiaowei*，LI Juping，GUO Junjie
（Ningbo CRRC Times Transducer Technique Co，Ltd，Ningbo Zhejiang 315021，China）

The static integral equation method can be used to simulate and design the fluxgate，which reduces the computation time and shows a satisfactory accuracy.However，the effect of frequency on fluxgate was not considered accurately based on static integral equation method，so it applies only to the quasi static analysis and can't be used to carry out the transient analysis of fluxgate.In order to solve this problem，the paper firstly studied the quasi static analysis of fluxgate with integral equation method according to the characteristics of fluxgate core.Based on the above work，the mathematical model was modified for fluxgate transient analysis.Finally a rod core fluxgate as an example was analyzed to validate the modified model.The simulation results agree well with experimental results.

fluxgate；transient analysis；integral equation method；quasi static analysis；core；frequency

TP212.1

A

1004-1699（2016）07-1026-06

项目来源：国家自然科学基金项目（F040508）；高等学校博士点专项基金项目（20126102110031）

2015-12-23修改日期：2016-02-27