周海强，李 勇，陈振茂

（西安交通大学 机械结构强度与振动国家重点实验室，西安 710049）

碳素钢作为一种强磁性物质，具有良好的强度、硬度、塑性和韧性，且价格低廉，已广泛应用于工业生产的各个领域［1］。然而由于腐蚀、磨损、意外载荷等原因可能导致构件表面或内部产生缺陷。目前针对铁磁性材料表面及内部缺陷的无损检测方法有超声检测、漏磁检测［2－3］、电磁超声检测和涡流检测等。由于漏磁检测和电磁超声检测方法都存在外加非均匀强磁场，可能导致材料内部产生非均匀磁场，进而由于磁导率和施加磁场的非线性相关引起磁导率的空间/时间变化。显然，以往采用均匀磁导率不适用于铁磁材料非均匀磁化问题及相关研究。探究非均匀磁化条件下铁磁材料内部磁场分布至关重要。

目前，对电磁超声的研究主要采用试验研究和数值模拟。裴翠祥等［4］建立了非磁性材料电磁超声数值模拟方法，但未考虑铁磁性材料。铁磁材料电磁超声无损检测数值模拟的关键之一是计算材料内部的磁场分布。目前对铁磁性材料内部磁场数值模拟中，通常假设磁导率是常数，不考虑非均匀磁化对材料内部磁场分布的影响。王庆国等［5］通过建立材料的μr－H函数关系计算了铁磁材料内部磁感应强度；高俊吉等［6］采用铁磁区标量磁位和空气区域用矢量H表示的有限元和边界元耦合算法计算了开域静磁场。有鉴于此，文章以麦克斯韦方程为理论基础，基于等效磁极化法开发了铁磁材料内部的磁化强度计算方法和程序。数值模拟中采用有限元与边界元相结合的方法，并采用等效电流模型计算永磁体空间磁场分布。通过对典型问题数值模拟结果的分析及与ANSYS计算结果进行比对，验证了方法的有效性。研究成果为铁磁材料电磁超声无损检测数值模拟研究打下了基础。

1 数值计算方法

1.1 静态磁场的简化计算

永磁体在充磁方向上存在一定的剩余磁化强度，也可认为是磁体内部分布着一定的磁化电流，使得永磁体对外显示磁性。对于一般永磁体，剩磁Br通常为常数，这是由于磁体内各点磁化强度M0相同，磁化电流仅分布于磁体表面，其等效面磁化电流Km＝M0×en（en为磁体表面外单位法向矢量）。当磁体内部磁化强度不均匀时（主要是静磁相互作用的影响），除磁体表面外，磁体内部还存在有体磁流密度Jm＝▽×M0。确定了磁体内磁化电流的分布及大小，就可以根据毕奥—萨伐尔定律确定空间任意点的磁场强度：

1.2 控制方程

通过引入矢量磁位A和库仑规范▽·A＝0，磁准静态问题控制方程可写为：

式中ΩFEM——材料区域；

ΩBEM——空气区域；

μ0——真空中磁导率；

M——铁磁材料磁化强度；

J0——激励电流。

对于强磁性材料，其非线性本构关系可以写为：

根据磁化强度的定义，M－B关系可用下式表达：

在材料区域和空气区域的分界面上，磁场强度切向分量的连续条件为［8］：

式中FEM——材料内部区域；

n1——边界外法向矢量；

BEM——材料外空间；

n2——内边界法向矢量。

根据有限元和边界元混合法思想，分别对式（3）进行有限元离散和对式（4）进行边界元离散，结合边界条件式（6）即可得有限元－边界元控制方程为［9］：

式中 ［P］，［K］，［D］，［G］——系数矩阵；

｛F0｝——和外磁场相关的向量；

｛M｝——虚拟磁化源向量；

［］——和虚拟磁化源相关的系数矩阵；

｛A｝——待求向量位。

在数值模拟算法中，首先假设初始磁化｛M0｝＝0，由式（7）和式（5）计算得｛Mi｝。并判断是否满足收敛条件‖Mi－Mi－1‖V＜ε。如果不满足收敛条件，则用｛Mi｝代替｛Mi－1｝对式（7）进行修正，直到能够满足收敛条件为止。通过B＝▽×A即可求得材料的内部的磁场分布。

2 数值模拟程序验证及计算实例

通过改变永磁体和激励线圈的组合方式，电磁超声可获得不同类型的超声波。而常见的永磁体主要有两种模型：矩形永磁体和U－型永磁体。

2.1 数值计算模型

图1为电磁超声无损检测中常见的两种永磁体模型。模型1中永磁体的尺寸为50mm×25mm×12.5mm，上下面为磁极。模型2中永磁体N极（S极）的尺寸为10mm×20mm×12.5mm，两极中心间距d0＝40mm，图中试件为碳素钢，其尺寸为200mm×200mm×15mm，该材料的B－H曲线如图2所示。

2.2 程序验证及计算结果

2.2.1 程序验证

通过有限元分析软件ANSYS对模型1进行了计算，得到了材料上表面的磁场强度分布。图3给出了y＝0截面距表面h＝0.5mm处ANSYS计算结果和程序计算结果的比较。横坐标为距永磁体中心的距离；纵坐标为磁感应强度。由图可见该程序的结果与ANSYS计算结果一致性较好，误差在允许范围，说明了方法的有效性。

2.2.2 等效磁化电流模型试验验证

以长，宽，高分别为a＝50mm，b＝25mm，h＝12.5mm的矩形永磁体为例，充磁方向沿厚度方向，磁体上下面为磁极。首先通过高斯计测得该矩形永磁体上表面中心处的磁感应强度Bz＝0.331 4T，其对应的等效面电流大小｜KM｜＝9.35×105A/m。根据式（1）计算所得距永磁体上表面中心不同距离处的磁场与实测量值的比较见图4。基于等效电流模型的计算结果与试验测试结果具有较好的一致性。

图4 距永磁体上表面中心h的磁场强度变化曲线

2.2.3 铁磁体内部磁场的分布

通过将永磁体等效为等效电流模型，可求得试件中磁场的分布。在条形永磁体作用下，试件上表面的磁场B的分布如图5（a）所示，试件上表面的相对磁导率分布如图5（b）所示。在U型永磁体作用下，试件上表面的磁场B的分布如图6（a），试件上表面的相对磁导率分布如图6（b）所示。由图5（a）和6（a）可得，在永磁体的正下方，铁磁材料内部的磁场强度较大，且磁体的边缘处磁场强度达到最大，中间磁场强度较小。同时，随着离磁体距离的增加，磁场强度衰减较快。由图5（b）和6（b）可得，在材料没达到磁饱和状态时，磁场较大处的相对磁导率也较大，且远高于平均值。所以在计算磁性材料内部磁场时，将磁导率设为常数将导致较大的误差。

图5 矩形永磁体磁感应强度和相对磁导率的分布

图6 U－型永磁体磁感应强度和相对磁导率的分布

3 结论

基于麦克斯韦方程组和等效磁极化方法，建立了铁磁体内部磁场的数值模拟的有限元与边界元法混合法程序。利用开发程序计算了不同永磁体作用下铁磁材料内部磁感应强度和相对磁导率的分布情况。在外加非均匀磁场作用下，铁磁材料内部的磁场分布也表现出较大的非均匀性。在永磁体正下方，磁场强度较大。在铁磁材料没有达到磁饱和的情况下，磁场较大处，其相对磁导率也较大，且远远高于其平均值。所以，在计算铁磁材料内部磁场时，将磁导率设为常数将会带来较大的误差。

［1］田民波.磁性材料［M］.北京：清华大学出版社，2001.

［2］Yong Li，John Wilson，Gui Yun Tian.Experiment and simulation study of 3Dmagnetic field sensing for magnetic flux leakage defect characterisation［J］.NDT＆ E International，2007，40（2）：179－184.

［3］Yong Li，Gui Yun Tian，Steve Ward.Numerical simulation on Magnetic Flux Leakage Evaluation at high speed［J］.NDT ＆ E International，2006，39（5）：367－373.

［4］裴翠祥，陈振茂.电磁超声的数值模拟方法［J］.无损检测，2008，30（9）：603－607.

［5］王庆国，张灵振，郝爱民，等.非均匀磁化条件下铁磁材料内部的磁场分［J］.测量技术学报，1998，12（3）：429－433.

［6］高俊吉，刘大明.边界棱边元FEM—BEM法在开域静磁场计算中的应用［J］.哈尔滨工程大学学报，2009，30（11）：1311－1316.

［7］Hantila I F，Preda G，Vasillu M.Polarization method for state fields［J］.IEEE Transactions on Magnetics，2000，36（4）：672－675.

［8］Fetzer J，Kurtz S，Lehner G.Comparison of analytical and numerical integration techniques for the boundary integrals in the BEM－FEM coupling considering TEAM workshop problem no.13［J］.IEEE Transactions on Magnetics，1997，33（2）：1227－1230.

［9］Mihalache O，Preda G，Uchimoto T，et al.Crack reconstruction in ferromagneic materials using nonlinear FEM－BEM scheme and neural networks［M］.in Electromagnetic Nonde－structive Evaluation（V），J.Pavo etl.（Eds.）：IOS Press，2001：67.