磁屏蔽装置结构参数的组合优化方法

2018-11-29吕志峰赵欣张金生王仕成李婷

兵工学报 2018年11期

吕志峰，赵欣，张金生，王仕成，李婷

(火箭军工程大学精确制导与仿真实验室，陕西西安 710025)

0 引言

对于高精度地磁场的测量，为了保证高精度磁测仪器的正常工作，必须消除附近铁磁性物体以及电控设备产生的干扰磁场，这就需要建立近零磁场环境。近零磁场环境在地磁导航、航天工程、国防工业、微粒子精密实验等方面具有重要的应用价值[1-5]。从国内外研究成果来看，大部分近零磁场环境均以被动磁屏蔽为主[6-9]，即采用高磁导率材料(如坡莫合金)构建一个封闭的磁屏蔽装置，从而保证该装置内部空间的磁场接近为“零磁”。

2000年德国国家物理技术研究院建成的零磁场实验室(BMSR-2)由1层射频屏蔽、7层坡莫合金和1层铝共9层材料组成，在0.01 Hz低频磁环境下，其被动磁屏蔽系数已经达到75 000[6]. 2014年，德国慕尼黑大学最新建设的磁屏蔽装置仅采用2层坡莫合金，距屏蔽门60 cm处剩磁小于1 nT，距直径小于6 cm的空洞5 cm处剩磁已经小于3 nT[2]. 从9层材料缩减为2层材料，不论是经济成本还是施工难度，都得到了极大降低。可见，磁屏蔽装置结构参数优化是设计者必须要考虑的问题，磁屏蔽装置结构参数优化研究具有重要的理论意义、工程意义和经济价值。

磁屏蔽效果分析方法是进行磁屏蔽装置结构参数优化的基础和前提。目前对于磁屏蔽效果的分析方法主要有解析法和数值法。解析法计算简单，主要用于工程估算，但是其精度较差[10-12]；数值法主要用于对设计方案可行性的论证及其技术指标的验证[13-14]。对于结构参数的优化，常规思路都是从解析法入手，寻找出解析式的最优解作为其最优结构参数[15-18]，但是解析法精度较低，用纯解析法寻优得到的最优结构参数其可信度较低，且有些形状的磁屏蔽装置很难推导出其解析式，导致无法进行解析寻优。同时，由于数值法计算量大且耗时较多，采用纯数值法进行磁屏蔽装置结构参数优化还存在一定的困难[15]。鉴于此，本文将解析方法精度偏低但计算速度快的特点与数值方法计算过程耗时但精度高的特点相结合，提出一种“解析+数值”组合优化的方法，对磁屏蔽装置结构参数进行优化设计，同时提出将两种方法交互验证的方法，新方法与纯解析优化相比，提高了优化结果的可信度。

1 磁场屏蔽理论

1.1 磁场屏蔽原理

衡量磁场屏蔽效果的参数是磁场屏蔽系数S，其计算公式为

(1)

式中：Bf为屏蔽前磁感应强度值；Bo为屏蔽后屏蔽体中心点处的磁感应强度。

磁场屏蔽可以采用磁路原理来解释。磁路原理类似于电路原理，即将1个高导磁率材料(如硅钢片、坡莫合金等)构成的屏蔽外壳置于干扰磁场中，这时屏蔽外壳与其内部的空气介质组成1个并联的磁路，由于空气介质的相对磁导率接近于1，而屏蔽外壳的相对磁导率能达到几千甚至上万，故空气的磁阻Ra要比屏蔽外壳的磁阻Rm大得多，当外界存在干扰磁场时，绝大部分磁通密度线会沿着磁阻低的屏蔽壳通过，使得屏蔽材料内部的磁感应强度Bi较大，而进入屏蔽体内腔的磁通量很少，从而达到屏蔽磁场的目的。磁场屏蔽原理如图1所示。

与电路原理一样，两个电阻并联中，即使其中1个电阻再大，也会有电流通过，因此，即使屏蔽层的磁导率再高，其内部空间也不可能做到绝对“零磁”。理论上，屏蔽体壳越厚，磁导率越高，其屏蔽效果越好。为了达到最佳的磁屏蔽性能，通常采用多层屏蔽的方法，将残余磁通量逐层屏蔽掉。

1.2 磁场屏蔽的解析分析方法

磁屏蔽的解析方法主要有磁路法或磁标位法。对于单层无限长磁性材料构成的圆柱形屏蔽罩，其横向磁屏蔽系数[19]为

(2)

式中：μr为材料相对磁导率；t为屏蔽层厚度；R为屏蔽层半径。

对于n层磁屏蔽层，其横向磁屏蔽系数为

(3)

对于单层有限长磁性材料构成的圆柱形屏蔽罩，其纵向磁屏蔽系数[19]为

(4)

式中：K为中间变量，

(5)

a=L/R，L为磁屏蔽罩长度，α和β分别可以通过测量圆柱端和侧壁的标准磁通量分布确定，通常取α=0.85±0.03，β=1.83±0.06.

对于n层磁屏蔽层，其纵向磁屏蔽系数为

(6)

1.3 磁场屏蔽的数值分析方法

分析磁场屏蔽问题的经典理论是求解麦克斯韦方程，但是由于磁场屏蔽问题中涉及到边界条件，使得该问题的分析复杂且计算繁琐，因此在20世纪90年代之前很少有人使用数值方法求解磁屏蔽问题。但随着计算机硬件的迅猛发展以及有限元法(FEM)理论的不断完善，数值分析方法求解磁场屏蔽问题得到了越来越广泛的应用。

FEM是以变分原理为基础逐步发展起来的，由于电磁场可以以拉普拉斯方程和泊松方程描述，而这类方程与泛函极值问题有着紧密联系，故可以采用FEM来分析磁场屏蔽问题。目前，世界上已有多款商用电磁场有限元分析软件，其基本思想都是基于麦克斯韦微分方程，采用有限元离散形式，将模型中的电磁场计算转变为庞大的矩阵求解。本文中数值计算采用的是Ansoft Maxwell软件，它不仅具备常规电磁场有限元分析软件的特点，还具有参数化建模的功能[20]，该功能对于屏蔽层结构参数的优化具有重要作用。

2 磁屏蔽装置结构参数的组合优化方法

解析方法虽然计算简单，但是由于在推导过程中做了一些近似处理，其计算精度不高，况且并不是所有形状的磁屏蔽装置都可以通过磁路法或磁标位法推导得到解析解，因此解析方法用于结构参数的优化存在一定局限性；有限元数值计算方法在计算过程中严格按照麦克斯韦微分方程进行求解，理论上计算精度高，从国内外研究成果看，其计算结果与实际系统吻合程度也较好[13-14]，但是由于求解过程中涉及到庞大的矩阵求解，计算耗时费力，故数值法一般多用于对设计方案的可行性和技术指标进行论证和验证。

鉴于解析法和数值法的优点和不足，本文提出一种“解析+数值”组合优化方法：对于给定的初始结构条件和待优化结构参数，首先确定装置的磁屏蔽性能计算是否存在解析表达式，若难以推导其解析表达式，则利用人工神经网络能够无限逼近非线性函数的优点，建立其磁屏蔽性能计算模型，得到其“黑箱解析式”；然后利用智能优化算法(如粒子群优化(PSO)算法、遗传算法(GA)、模拟退火(SA)算法及这些算法的改进形式等)强大的寻优能力，在短时间内对解析式或“黑箱解析式”进行解析优化，找到解析最优解，达到“粗”寻优的目的；基于解析法得到的优化结果在Ansoft Maxwell软件中进行数值建模，利用软件的参数化建模功能，将待优化参数设定为参数变量，以较小的步长在解析优化结果附近进行数值扫描计算，即可得到1条或1族曲线，通过曲线分析得到数值最优解，从而实现“精”寻优；最后，通过双重交互验证的方法对优化结果的可信度进行评估，将解析计算得到的“粗”优解与数值计算得到的“精”优解进行比较，判断两种解是否一致，在两种解接近的情况下再将“精”优解作为最终的结构参数代入解析法中，将解析法求得的屏蔽系数与数值法求得的屏蔽系数进行比较：如果二者偏差在可接受范围内，则说明两种方法具有一致性，通过组合优化得到的“精”优解可以作为最优解；如果二者偏差过大，则说明优化过程存在问题，需要重新优化。造成“粗”优解与“精”优解偏差过大的原因可能是解析寻优或数值寻优过程出现了问题，因此重新优化可以从如下两方面进行：1)解析寻优是通过智能优化算法实现的，由于智能优化算法寻优可能会落入局部最优(如PSO算法)，可以采用智能优化算法的改进形式或增加种群规模，跳出局部最优，使最终的解析寻优解更加接近实际最优解；2)数值寻优采用的是FEM，由于FEM的计算精度与模型的网格剖分大小有关，剖分越细，计算精度越高，可以减小网格的剖分，重新进行数值寻优。组合优化流程如图2所示。

3 磁屏蔽装置结构参数优化算例分析

磁屏蔽装置的应用极其广泛，对于不同的工程应用背景，磁屏蔽装置结构可能差别很大(如圆柱形、球形或矩形，单层或多层等)，限于篇幅，这里仅对有解析式和无解析式两种情况，各通过一个仿真优化算例来说明本文方法的有效性，其他形状的磁屏蔽装置优化均可按照图2所示的优化过程进行。

3.1 解析式已知情况下的结构参数优化

3.1.1 结构参数初始化及待优化参数的确定

由1.2节可知，圆柱形磁屏蔽装置的磁屏蔽性能计算存在解析式，且这种形状的磁屏蔽装置应用很广泛。下面对圆柱形磁屏蔽装置的结构参数进行优化。

3.1.2 智能优化算法快速“粗”寻优

目前，常用的智能优化算法有PSO算法、GA、SA算法及这些算法的改进形式等，本文采用自适应权重PSO算法进行解析寻优，它是一种基于群体的随机优化技术，具有群体智能、迭代格式简单、内在并行性及可快速收敛到最优解附近等优点，尤其适用于多变量优化。算法的基本原理不再赘述，详情可参考文献[21].

将Smin=min(St,Sa)作为目标函数，以R1、R2和L2作为待优化变量，采用自适应权重PSO算法进行寻优。为了避免内层圆柱形装置尺寸过长或者过扁，引入限定条件15 mm≤R1≤35 mm，经求解得到最优结构参数为：R1=29.76 mm，L1=35.92 mm，R2=44.24 mm，L2=99.71 mm，此时装置的磁屏蔽系数达到最大值，为S=10 833.

3.1.3 数值法扫描“精”寻优

考虑到数值法计算复杂且需要消耗大量的内存，本文采用图形工作站进行处理。根据3.1.2节中的优化结果，在Ansoft Maxwell软件中建立内层半径R1=29.76 mm、长度L1=35.92 mm和外层半径R2=44.24 mm、长度L2=99.71 mm的两层空腔圆柱形磁屏蔽装置，屏蔽层厚度均为t=1 mm，相对磁导率设置为μr=10 000. 同时在圆柱形屏蔽装置的横向和纵向分别设置两块平行放置的永磁体材料NdFe30，采用平行充磁的方式产生稳定的横向磁场和纵向磁场，为了更接近实际情况，将两个方向的磁场设置为地磁场值附近，均为52 278.3 nT，双层圆柱形磁屏蔽装置三维数值仿真模型如图3所示。

为了精细寻优，需要对R1、R2和L2进行参数化设置。本文以3.1.2节中的“粗”优解为搜索中心，对其±3.00 mm范围(即R1为26.76～32.76 mm、R2为41.24～47.24 mm、L2为96.71～102.71 mm)进行扫描计算，扫描步长设定为0.5 mm. 将计算结果导出，在MATLAB中绘图。由于每个参数的取值有13种情况，3个参数对应的计算有133=2 197个计算结果，限于篇幅，这里仅对部分计算结果(R1分别为26.76 mm、28.76 mm、30.76 mm、32.76 mm)进行展示，如图4所示。

对数值法扫描计算得到的结果进行统计后可知：当R1=28.76 mm、L1=38.48 mm、R2=45.24 mm、L2=98.21 mm时，Smin=min(St，Sa)达到最大值，为S=9 246. 至此，采用数值法得到了结构参数的“精”优解。

3.1.4 交互验证

将解析法的“粗”优解与数值法的“精”优解进行比较，结果如表1所示。

从表1中可以发现：对于待优化参数R1、R2和L2，通过解析法和数值法求得3个参数的偏差分别为1.00 mm、1.00 mm、1.50 mm，均未超过±3.00 mm搜索范围的一半，可以看做两种解具有一致性；对于因变量参数L1，两种方法求得的偏差为2.56 mm，由于该参数随着R1的改变而改变，可能会出现较大的偏差，但是2.56 mm的偏差相比于解析解35.92 mm和数值解38.48 mm，波动较小，故两种解也可以看做是一致的。这一结论表明“精”寻优结果是正确的。

将“精”优解代入解析(2)式和(4)式中，求得总体磁屏蔽系数Smin=min(St，Sa)=10 128，3.1.3节中数值法求得的总体磁屏蔽系数为9 246，二者相对误差为8.71%. 需要说明的是，由于解析法中横向磁屏蔽系数计算公式(即(2)式)的适用范围为无限长圆柱形磁屏蔽罩，纵向磁屏蔽系数计算公式(即(3)式和(4)式)中的α和β使用的是经验值，其真值需要测量圆柱端和侧壁的标准磁通量分布才能准确确定，故从理论层面而言，解析法求得的屏蔽系数与数值法求得的屏蔽系数必然存在较大偏差，但是从工程层面而言，二者相对误差为8.71%，是工程上可以接受的误差范围，可以看作是一致的，从而验证了“精”寻优结果的正确性。

通过上述解析法与数值法的交互验证，有力地说明了“精”寻优结果的正确性，并最终可以确定其最优结构参数为R1=28.76 mm，L1=38.48 mm，R2=45.24 mm，L2=98.21 mm. 同时，通过这个仿真案例可以发现：由于解析法计算精度较低，传统的单纯用解析法进行寻优得到的优化结果，其可信度必然也较低，而本文提出将计算精度较高的数值法引入磁屏蔽装置结构参数的优化，并通过交互验证，极大地提高了寻优结果的可信度。

3.2 解析式未知情况下的结构参数优化

3.2.1 结构参数初始化及待优化参数确定

目前，矩形磁屏蔽装置的应用非常广泛，但是其磁屏蔽性能的解析表达式很难推导，下面对解析式未知的矩形磁屏蔽装置结构参数进行优化。

首先设定其初始条件，设定装置为单层结构，其长、宽、高分别为X、Y、Z，内部体积满足V=X×Y×Z=100 cm3. 设定屏蔽材料的相对磁导率μr=10 000，屏蔽层厚度t=1 mm.

对于单层矩形磁屏蔽装置，其结构参数有3个，即X、Y和Z，由于V=X×Y×Z=100 cm3这一限定条件，即X、Y和Z中只有两个变量是独立的，设定X和Y为自变量，Z=V/(X×Y)为因变量，故待优化参数为X和Y.

3.2.2 人工神经网络建立“黑箱解析式”

由于矩形磁屏蔽装置的结构特点，导致磁路法很难用于推导其解析解，而人工神经网络通过调整权重和阈值，可以“学习”或发现变量之间的关系，尤其适用于复杂的非线性映射关系，因此，可以采用人工神经网络建立矩形装置的磁屏蔽性能计算模型，将该网络模块作为其“黑箱解析式”。目前，常用的人工神经网络模型主要有BP神经网络、径向基函数(RBF)神经网络、反馈型神经网络等，本文采用RBF神经网络进行建模，它具有逼近精度高、网络规模小和学习速度快等优点。模型原理不再赘述，详情可参考文献[22].

采用RBF神经网络建模，必须要有一定的输入数据和输出数据用于网络训练。针对磁屏蔽性能计算的RBF神经网络模型，其输入为装置的结构参数，输出为磁屏蔽系数。理论上，不同结构尺寸的装置对应不同的磁屏蔽系数，理想情况下，期望在实际中构造大量不同结构尺寸的装置，再进行实际测量，以得到其磁屏蔽系数，从而获取足够的训练数据，但这必然需要投入大量人力、物力和财力，显然是不合理的。鉴于数值法的计算结果与实际系统吻合程度较高，因此可以采用数值法获取训练数据，以大大降低经济成本和时间成本。

在Ansoft Maxwell软件中建立X=Y=Z=46.41 mm的单层空腔矩形磁屏蔽装置，屏蔽层厚度t=1 mm，相对磁导率设置为μr=10 000，使用永磁体材料NdFe30产生52 278.3 nT的背景磁场，单层矩形磁屏蔽装置三维数值仿真模型如图5所示。

为了尽可能多地获取训练数据，对长X和宽Y进行参数化设置，两个参数的取值范围均为10～100 mm，变化步长设为5 mm，则可以得到192=361组数据，计算结果如图6所示。

将361组数据分为两部分：随机抽取其中90%用于RBF神经网络训练，剩下10%用于验证网络精度。验证结果如图7所示。

从图7(a)可以看出，采用RBF神经网络建立的网络模块求解结果与数值法求解结果基本一致。图7(b)是二者之间相对误差的统计结果。由图7(b)可见，二者最大相对误差仅为3.09%，表明RBF神经网络建立的网络模块精度很高，可以用作磁屏蔽性能计算的“黑箱解析式”。

3.2.3 智能优化算法快速“粗”寻优

以矩形磁屏蔽装置的磁屏蔽系数S为目标函数，以长X和宽Y为待优化变量，采用自适应权重PSO算法，对3.2.2节中建立的磁屏蔽性能计算“黑箱解析式”进行寻优。为了避免矩形装置尺寸过长或者过扁，引入限定条件20 mm≤X≤70 mm、20 mm≤Y≤70 mm、20 mm≤Z≤70 mm，经求解，得到最优结构参数：X=32.44 mm，Y=45.23 mm，Z=69.98 mm. 此时装置的磁屏蔽系数达到最大值，为S=202.96.

3.2.4 数值法扫描“精”寻优

数值计算仍采用图5所示的三维仿真模型。为了精细寻优，需要对长X和宽Y进行参数化设置。以3.2.3节中的“粗”优解为搜索中心，对其±3.00 mm范围(即X为28.94～34.94 mm、Y为41.73～47.73 mm)进行扫描计算，扫描步长设定为0.50 mm. 将计算结果导出，在MATLAB中绘图。数值扫描计算结果如图8所示。

对数值法扫描计算得到的结果进行统计后可知：当X=32.44 mm、Y=45.23 mm、Z=68.15 mm时，矩形磁屏蔽装置的磁屏蔽系数达到最大值，为S=208.61. 至此，采用数值法得到了结构参数的“精”优解。

3.2.5 交互验证

将3.2.3节中解析法的“粗”优解与3.2.4节中数值法的“精”优解进行比较，比较结果如表2所示。

表2 解析法与数值法求解结果对比

从表2中可以发现：对于待优化参数X和Y，通过解析法和数值法求得的两个参数其偏差均为0.50 mm，相比于±3.00 mm搜索范围要小得多，可见两种解具有一致性；对于因变量参数Z，两种方法求得的偏差为1.83 mm，由于该参数随着X和Y的改变而改变，可能会出现较大偏差，但是1.83 mm的偏差相比于解析解69.98 mm和数值解68.15 mm，波动较小，故两种解也可以看做是一致的。这一结果表明“精”寻优结果是正确的。

将“精”优解代入3.2.2节建立的“黑箱解析式”中，求得装置的磁屏蔽系数S=200.18，3.2.4节中数值法求得的磁屏蔽系数为208.61，二者相对误差为4.04%，是工程上可以接受的误差范围，二者可以看做是一致的，从而也验证了“精”寻优结果的正确性。

通过解析法与数值法的交互验证，有力地说明了“精”寻优结果的正确性，并最终可以确定其最优结构参数为X=32.44 mm，Y=45.23 mm，Z=68.15 mm.

通过与3.1.4节的交互验证结果进行对比可以发现：对采用人工神经网络建立的“黑箱解析式”进行“粗”寻优得到的结果与数值法“精”寻优得到的结果，其相对误差和对应的磁屏蔽系数相对误差都要比3.1.4节中的相对误差低很多，这是因为人工神经网络能够高精度地拟合复杂的非线性映射关系，且网络训练数据来源于精度较高的数值计算结果，因此本文建立的“黑箱解析式”精度要优于传统经过近似处理推导得到的解析式精度，使得“粗”寻优的结果更加接近真实的最优解，大大降低了后续“精”寻优陷入局部最优的可能。