缆索作为斜拉桥、悬索桥等桥梁的关键部件之一，对其进行检测具有重要意义[1]．目前磁性检测和磁致伸缩导波等电磁声检测方法已在缆索检测中有着广泛的应用[2-3]，如何确定合适的缆索内部磁化强度是上述方法应用中遇到的难题．为此，须要建立缆索磁化模型，准确获取缆索内部磁场空间分布，为磁化器结构优化和传感机理的定量研究提供理论基础．在磁性检测和导波检测中，由于缆索结构的复杂性，其精细化磁化器设计一般只能采用有限元法．但闭合式磁化器结构优化的设计变量多，且变化范围广，因此有必要先确定合理的变量初值及其变化范围以缩短计算时间．目前一般先采用磁路法初步确定缆索磁化器结构[4-5]，但磁路法无法获取缆索内部的磁场分布，造成多层次钢丝束中断丝和腐蚀等缺陷难以定位的问题；同时，磁路法无法给出磁化器和缆索各个参数之间的定量化关系，故当缆索规格变化时，磁化器须要重新设计．为此，本研究提出缆索闭合式磁化的解析模型，利用解析法[6-7]建立磁场分布与模型结构和磁参数的显式表达式，更加精确地描述空间恒定磁场分布，并且当缆索规格变化时，可以直接获取磁化器相关参数，降低磁化器设计工作量，同时具有计算时间短、不存在数值计算不稳定问题的优点．为此本研究建立基于等效面磁荷的缆索闭合式磁化解析模型，应用截断区域特征函数展开法(TREE)[8]及传递矩阵[9]求解磁标位方程，实现对缆索内部磁场的计算；同时用有限元仿真验证解析模型的正确性，并对误差来源进行分析．1 缆索闭合磁化模型解析计算1.1　缆索闭合磁化模型建模图1为长直缆索外套有环形磁化器的漏磁检测示意图，其中缆索内部是由多根钢丝组合而成的钢丝束，其截面为不规则多边形，缆索外层由聚乙烯(Polyethylene，PE)包裹．10.13245/j.hust.230609.F001图1长直缆索漏磁检测示意图建模计算实际缆索中静态磁场分布极其困难，考虑到钢丝间的气隙比较小，在建模中将缆索钢丝束简化处理为实心的长直钢棒，即认为钢丝间不存在间隙，且钢丝束截面近似为圆形．环形磁化器穿过缆索，且与长直钢棒同轴，故简化后的模型如图2所示，具有轴对称性，其磁场分布与角度无关，可以在二维柱坐标系下进行求解．由于永磁体的相对磁导率近似为1，因此可以不考虑永磁体沿其轴向的边界面条件．而衔铁相对磁导率很高，须要考虑其上下两端的边界面条件．但在漏磁检测中，磁化器极靴中间区域的磁场对检测更为重要，即[-z2，z2]之间的磁场分布，从而忽略衔铁上下端面的边界问题，将图2所示的二维模型进一步简化为一维模型．10.13245/j.hust.230609.F002图2二维柱坐标下缆索闭合式磁化模型在对模型进行简化后，将模型由内向外划分为图2所示的5个求解区域．其中区域1为无限长的钢棒，半径为r1，其金属截面积等于缆索钢丝束金属截面积[10]；区域2为PE层，外径r2，在计算中相对磁导率设为1；区域3为永磁体层，外径为r3，在[-z2，-z1]和[z1，z2]区域内为永磁体，整层的磁导率可视为1；区域4为衔铁层，外径为r4，由于不考虑其轴向的边界条件，可视为无限长直金属层；区域5为空气层，相对磁导率为1．各区域的磁导率分别记为μi，i为区域序号，i=1，2，…，5．1.2　恒定磁场解析计算在整个区域内使用磁标位对磁场进行计算，磁标位与电位类似，表示磁场在某位置的势，为与电位区分，使用φm表示(m表示磁场)．永磁体区域采用等效面磁荷法进行处理，认为磁荷仅分布在环形永磁体的内外表面，即r=r2和r=r3的[-z2，-z1]∪[z1，z2]区域．此外，其他区域均无电流与磁荷存在，假设各层材料为均匀线性材料，对于区域1～5，磁标位均满足拉普拉斯方程∇2φm=0．(1)在柱坐标系下使用分离变量法可以得到磁标位的通解φm(i)=∫0∞AiI0(αr)+BiK0(αr)∙Cicos(αz)+Disin(αz)dα, (2)式中：α为任意实数；r和z为柱坐标系下的位置坐标；I0(αr)和K0(αr)分别为第一类和第二类虚宗量贝塞尔函数．根据截断区域特征函数展开法原理，在z=±h处截断，并在截断边界上施加狄利克雷边界条件φm(r，±h)=0，从而积分解可以转化为下式的级数形式φm(i)(ρ,z)=∑n=1NsAiI0(αnr)+BiK0(αnr)∙Cicos(αnz)+Disin(αnz), (3)式中：NS为级数累加项数；α=nπ/h，n=1，2，…．由于K0(αr)在r=0处发散，I0(αr)在r趋向于无穷大时也是发散的，因此在区域1有B1=0，在区域5有A5=0；同时，由于磁场Bz关于z=0对称，磁场Br关于z=0反对称，因此有Ci=0．那么各个区域的磁标位表达式如下：φm(1)(r,z)=∑n=1NsA1I0(αnr)sin(αnz)；(4)φm(2)(r,z)=∑n=1NsA2I0(αnr)+B2K0(αnr)sin(αnz)；(5)φm(3)(r,z)=∑n=1NsA3I0(αnr)+B3K0(αnr)sin(αnz)；(6)φm(4)(r,z)=∑n=1NsA4I0(αnr)+B4K0(αnr)sin(αnz)；(7)φm(5)(r,z)=∑n=1NsB5K0(αnr)sin(αnz)．(8)接下来根据各层磁场的交界面条件求解待定系数Ai和Bi，先考虑除区域3永磁体以外的边界条件，根据磁场在交界面的法向边界条件和切向的连续性，有：μi∂φm(i)(r,z)∂rr=ri=μi+1∂φm(i+1)(r,z)∂rr=ri;φm(i)(r,z)r=ri=φm(i+1)(r,z)r=ri. (9)将磁标位表达式代入交界面条件，联立式(9)可得从第i层向第i+1层推导的系数传递矩阵和第i+1层向第i层推导的系数传递矩阵分别为Ai+1Bi+1=ΓAi→Ai+1ΓBi→Ai+1ΓAi→Bi+1ΓBi→Bi+1AiBi；(10)AiBi=ΓAi+1→AiΓBi+1→AiΓAi+1→BiΓBi+1→BiAi+1Bi+1，(11)式中：ΓAi+1→Ai=[μiI0(αnri)K1(αnri)+μi+1I1(αnri)K0(αnri)]/Si;ΓBi+1→Ai=(μi-μi+1)K0(αnri)K1(αnri)/Si；ΓAi+1→Bi=(μi-μi+1)I0(αnri)I1(αnri)/Si；ΓBi+1→Bi=[μi+1I0(αnri)K1(αnri)+μiI1(αnri)K0(αnri)]/Si;Si=μi[I0(αnri)K1(αnri)+I1(αnri)K0(αnri)]；ΓAi→Ai+1=[μi+1I0(αnri)K1(αnri)+μiI1(αnri)K0(αnri)]/Ti;ΓBi→Ai+1=(μi+1-μi)K0(αnri)K1(αnri)/Ti；ΓAi→Bi+1=(μi+1-μi)I0(αnri)I1(αnri)/Ti；ΓBi→Bi+1=[μiI0(αnri)K1(αnri)+μi+1I1(αnri)K0(αnri)]/Ti;Ti=μi+1[I0(αnri)K1(αnri)+I1(αnri)K0(αnri)]．1.3　等效面磁荷法下面采用等效磁荷法处理区域3永磁体内外表面的交界面条件，设永磁体内外表面的面磁荷分别为σm2和σm3．但是对于模型中的环形径向磁化永磁体，除磁极端面上存在面磁荷外，其内部的体磁荷密度不能视为0．为简化模型便于分析，以正负磁荷总量相等为基本原则，重新考虑两个内外磁极端面上的面磁荷密度[11]，有σm2/σm3=r3/r2．故不妨令σm2=kσ0(r2+r3)/r2(z∈[-z2,-z1]⋃[z1,z2]);σm3=-kσ0(r2+r3)/r3(z∈[-z2,-z1]⋃[z1,z2]), (12)式中：k为0～1之间的常数，称为磁荷分配系数，与永磁体的结构尺寸相关；σ0=μ0Hc，其中Hc为永磁体的矫顽力．据此可以写出区域3的两侧法向交界面条件μi∂φm(i)(r,z)∂rr=ri=μi+1∂φm(i+1)(r,z)∂rr=ri+σmi，(13)进而可以得到从第2层和第4层向第3层的系数传递矩阵分别为A3B3=ΓA2→A3ΓB2→A3ΓA2→B3ΓB2→B3A2B2+P2Q2；(14)A3B3=ΓA4→A3ΓB4→A3ΓA4→B3ΓB4→B3A4B4+M3N3，(15)式中：M3=2K0(αnr3)αnhT3∫z1z2σm3sin(αnz)dz；N3=-2I0(αnr3)αnhT3∫z1z2σm3sin(αnz)dz；P2=2K0(αnr2)αnhS2∫z1z2σm2sin(αnz)dz；Q2=2I0(αnr2)αnhS2∫z1z2σm2sin(αnz)dz．利用系数传递矩阵分别从第1层和第5层向第3层推导，联立式(10)、(11)、()和(15)可以得到A1=(N3-Q2)ΓB5→A3+(P2-M3)ΓB5→B3ΓB5→A3ΓA1→B3-ΓB5→B3ΓA1→A3；B5=(N3-Q2)ΓA1→A3+(P2-M3)ΓA1→B3(ΓB5→A3ΓA1→B3-ΓB5→B3ΓA1→A3)；A2=A1ΓA1→A2；B2=A1ΓA1→B2；A3=A1ΓA1→A3+P2；B3=A1ΓA1→B3+Q2；A4=B5ΓB5→A4；B4=B5ΓB5→B4，式中：ΓA1→A3=ΓA1→A2ΓA2→A3+ΓA1→B2ΓB2→A3；ΓA1→B3=ΓA1→A2ΓA2→B3+ΓA1→B2ΓB2→B3；ΓB5→A3=ΓB5→A4ΓA4→A3+ΓB5→B4ΓB4→A3；ΓB5→B3=ΓB5→A4ΓA4→B3+ΓB5→B4ΓB4→B3．将待定系数代入可以计算出各层磁场为：Bz(i)(r,z)=-μi∑n=0NsαnAiI0(αnr)+BiK0(αnr)cos(αnz); (16)Br(i)(r,z)=μi∑n=0NsαnBiK1(αnr)-AiI1(αnr)sin(αnz). (17)根据上述解析表达式，通过数学软件可以求解缆索内外部任意一点(r，z)处的磁感应强度，进一步地，可以计算出磁场的空间分布．2 有限元仿真验证2.1　模型验证由于实验直接测量缆索内部磁场比较困难，同时有限元方法广泛应用于解析模型验证[12]，因此采用有限元仿真方法对解析结果进行验证．解析计算中，使用Matlab R2021a编程实现式(16)和(17)中磁场的计算；有限元计算中，在Ansys Electronics 2017中建立与图2相同的二位轴对称有限元模型进行仿真求解，使用静态磁场求解方案．有限元仿真中使用如图3所示的磁化器结构，其中两块永磁体完全相同，材料设置为N52，矫顽力为955 kA/m，方向相反．磁化器的尺寸上，解析计算和有限元仿真相同，永磁体宽度Wm和厚度Tm分别为30 mm和20 mm，磁极间距Gm为50 mm，衔铁厚度Ta为20 mm．10.13245/j.hust.230609.F003图3磁化器结构示意图被测缆索钢丝束外径均为109 mm，根据实验经验相对磁导率设为300，PE护套厚度为5 mm．同时，解析计算中，设置截断区域长度h=100z2，累加项数Ns=1 000，磁荷分配系数为0.55．首先验证式(16)的轴向磁感应强度Bz，分别提取仿真结果中r=0，20，40 mm三条线上的轴向磁感应强度数据与相应的解析计算结果作对比，结果如图4(a)所示，图中：实线表示解析计算；虚线表示有限元仿真；三种颜色分别对应径向不同位置的计算结果．同理，提取仿真结果中z=0，10，20 mm三条线上的轴向磁感应强度数据，与解析计算结果比较，结果如图4(b)所示．10.13245/j.hust.230609.F004图4轴向磁感应强度的验证由图4可知解析求得的轴向磁感应强度分量Bz随r坐标和z坐标变化的趋势与仿真结果基本一致．进一步分析可知：轴向磁感应强度Bz在磁化器中间区域[-55，55] mm的计算误差较小，但沿缆索轴向远离磁化器方向计算误差越来越大，甚至磁场方向发生反转．根据1.1节中的简化，该误差主要来源于未考虑衔铁上下表面的边界面条件．但由于磁化器中间区域的磁场在检测中更为重要，因此只考虑[-z2，z2]区域磁感应强度的情况下，解析计算得到磁场分布趋势正确，存在一定程度的误差，后面对误差进行进一步讨论．同时，提取仿真结果中r=10，30，50 mm三条线上的径向磁感应强度数据，与根据式(17)计算得到的结果进行对比，结果如图6(a)所示．提取仿真结果中z=10，20 mm和r=25 mm三条线上的径向磁感应强度数据，与解析计算结果比较，结果如图5(b)所示．10.13245/j.hust.230609.F005图5径向磁感应强度的验证10.13245/j.hust.230609.F006图6轴向磁感应强度Bz的误差由图5(a)可以看出：径向磁感应强度Br在永磁体中心z=(z1+z2)/2处的计算绝对误差最大，由永磁体中心向两侧误差逐渐减小．图5(b)中，径向磁感应强度Br在两磁极中间区域[-55，55] mm的范围内随r坐标变化的趋势与仿真结果也是基本一致的．径向磁感应强度Br和轴向磁感应强度Bz的误差来源是不同的，后续将逐一进行分析．2.2　模型误差与讨论下面对模型的计算误差进行更为详细的讨论．前面已经分析过由于未考虑衔铁上下端面边界的简化处理，轴向磁感应强度Bz仅在磁化器区域[-55，55] mm的计算误差较小．那么在该区域内，Bz的相对计算误差e的变化情况如图6(a)所示．由图6(a)可以看出：磁感应强度Bz的解析计算结果相对计算误差沿r坐标和z坐标增大而增大，也就是说缆索中心磁场的计算精度更高，而沿缆索轴向和径向向外的计算误差则逐渐增大．但若将计算区域进一步限制在两磁极之间区域[-25，25] mm时，如图6(b)所示，缆索内部磁场的计算误差均可控制在3%以内．同样地，在磁化器中间区域[-55，55] mm之间，径向磁感应强度Br的相对计算误差e如图7所示．由图7可以看出：径向磁感应强度Br沿轴向和径向分布的解析计算结果均偏大，缆索内部的计算误差均不大于12%．缆索内部外层次的相对计算误差在永磁体区域较大，其主要源为等效均匀面磁荷的假设，实际情况下面磁荷不能视为均匀分布，因此该区域的计算误差较大．10.13245/j.hust.230609.F007图7径向磁感应强度的误差3 结语本研究建立了缆索磁化闭合式磁化解析模型，首先根据实际情况对模型进行简化，采用截断区域特征函数展开法和传递矩阵求解磁标位方程，并采用等效面磁荷法对永磁体进行等效，得到缆索内部空间磁场分布表达式，然后使用有限元仿真对解析计算结果进行验证，结果表明：解析求得的磁感应强度分量空间变化的趋势与仿真结果基本一致；轴向磁感应强度的计算误差沿缆索轴向远离磁化器方向计算误差增大，磁极间隙区域内缆索各层次磁场的计算误差可控制在3%以内；而径向磁感应强度在永磁体中心处的计算绝对误差最大，但各层次磁场的计算误差均不大于12%，故该解析模型可用于求解缆索内部磁场分布，为缆索磁化强度选择提供了一种新方法．须要指出的是：本研究采用金属截面积相等原则将缆索等效为钢棒，忽略了缆索几何形状与其钢丝束的影响，但其可作为初步模型，为磁化器精细化设计提供初值．