基于极限梯度提升的完美匹配单层智能算法实现航空瞬变电磁问题高效吸收*

冯乃星 王欢 朱子贤 董纯志 李宏杨 张玉贤† 杨利霞 黄志祥‡

1) (安徽大学电子信息工程学院,合肥 230601)

2) (安徽大学,智能计算与信号处理教育部重点实验室,合肥 230601)

对于航空瞬变电磁的低频探地问题,除了精度和效率需要考虑,深地探测问题的复杂度也不容忽视,特别是对于低频复杂问题存在异常体与背景间的多尺度效应.为了模拟开域问题,有限厚度区域的完全匹配层被用于截断计算域,然而这也无形中增大了整个模型,造成计算复杂度增加.鉴于此,提出了一种新的基于极限梯度提升(extreme gradient boosting,XGB)的完美匹配单层模型,并将该模型集成到时域有限差分求解器中,以进一步提高时域有限差分仿真的性能.所提出的基于XGB 的完美匹配单层模型通过特征注意力集成学习方法可以获得更高的精度,同时占用更少的内存、消耗更少的时间.此外,由于该模型依托于传统机器学习模型,因此它在模型训练的稳定性和轻量级方面具有显著的优势.最后,通过对航空瞬变电磁应用进行三维数值模拟,验证了该方法的有效性和稳定性.该模型不仅在精度、效率和问题复杂性方面具有优势,而且还可以成功地集成到时域有限差分求解器中,解决低频航空瞬变电磁问题.

1 引言

基于极低频激发源的航空瞬变电磁(airborne transient electromagnetics,ATEMs)问 题[1–3]可以用时域有限差分(finite-difference time-domain,FDTD)方法[4–11]来解决,这是应用最广泛的先进数值方法之一.目前,显式和隐式的FDTD 方法分别被用于处理ATEM 问题,而这两种方法本身面临着不同的挑战.对于显式方法,主要存在两大挑战: 1) 低频地下探测问题需要非常高的空间采样密度,这样就需要大量的时间步长;2) 准确、高效地求解电磁波与开放区域内的极损介质相互作用需要通过吸收边界条件(absorbing boundary condition,ABC[12–15])来截断物理区域.对于隐式方法,受色散误差、分割误差、近似误差和非对称效应[16–19]的影响,会衍生出严重的计算精度误差问题;此外,稀疏矩阵的求逆计算也是不可避免的.无论是显式方法,还是隐式方法,都有其各自要面对的难题,而本文就是针对前两种方法面对的挑战,通过引入新的方法或技术,使得无论在精度、效率还是仿真问题复杂度上都能够取得好的性能.

近年来,随着人工智能(artificial intelligence,AI)的蓬勃发展,基于机器学习(machine learning,ML)的实现[20]被广泛应用于许多不同的工程和科学问题上.ML 能够从具有相同模式的训练数据中提取潜在的映射规律,从而能够预测并产生一个新的输出.因其自适应、非线性建模和数据驱动等特性,ML 可以移植到多种应用中,如基于超材料的逆设计[21,22]、超表面成像仪[23]、可编程光子电路[24]、深度学习皮肤传感器[25]、一种可基于石墨烯在吸收和透射间快速静电调制的频选反射器[26],基于层状石墨烯光子结构的优化单向角度不敏感超宽带吸收器[27]等.

最近,基于ML 的完美匹配层(perfectly matched layer,PML)分裂场方案的双曲正切基函数(hyperbolic tangent basis function,HTBF)[26]模型被用于取代传统的有限厚度的ABC 区域,规避了根据当前和之前步骤中的局部和邻近场来获取局部场的整个模拟过程中的圆形计算处理.这充分证明,ML 技术可以被应用到数值方法中,以避免在ABC 区域内的重复过程.

随后,为了进一步有效地削弱强隐失波,并在更大程度上节省计算机资源,开发了基于深度可微森林(deep differentiable forest,DDF)的完美匹配单层(perfectly matched monolayer,PMM)模型和基于非分裂场的复频率偏移(complex-frequencyshifted,CFS)CFS-PML 方案用于二维电磁问题模拟[27].通过仿真结果验证,基于DDF 的PMM 不仅可以取代传统的PML,而且具有较高计算精度并可以成功集成到FDTD 的计算过程中.

继上述二维EM 问题,将本算法应用于三维微带传输线问题的求解可进一步验证ML 方法的通用性、准确性与高效性.因此,选择可以利用可区分的注意机制来访问外部内存库中数据/程序的神经图灵机(neural turing machine,NTM)模型合并到基于DDT[28]的FDTD 求解器中来提高计算效率.所提出的基于DDT 的NTM 模型对三维微带传输线可获得约95%的高精度.

然而,在提出的基于DDT 的NTM 模型[28]中,只考虑了单极CFS-PML 方案,因此它在吸收低频传播波方面的效果远不如基于单极拉伸坐标的PML(stretched-coordinate based PML,SCPML)方案,特别是对于三维低频ATEM 应用.为了克服上述有/无CFS 方案的局限性,提出了一些多极PML[29–31]可以明显衰弱强低频传播和隐失波.此外,双极PML(two-pole PML,TP-PML)被验证是一种更好的选择[32],与多极PMLs(multipole PML,MP-PML)(n>2)相比,它不仅需要更少的CPU 时间和内存,而且保持了几乎相同的吸收精度.

本文主要关注进一步提高低频ATEM 问题的算法精度和降低问题复杂度,因此考虑了基于Z变换技术的TP-PML 方案对开放区域问题进行建模.为了实现更高的预测精度和进一步降低问题的复杂度,本文提出了一种高效、准确的基于极限梯度提升(extreme gradient boosting,XGB)的PMM 模型来取代传统的有限厚度PML.所提出的基于XGB 的PMM 模型通过特征注意的集成学习方法,可以获得更高的精度,同时实现更少的内存和CPU 时间消耗.此外,基于传统ML 模型的特点,该模型在模型训练的稳定性和模型轻量级方面具有显著的优势.这种方法在这里被称为基于XGB 的PMM 模型.在这项工作中,该方案不仅对三维低频ATEM 问题具有较高的精度和效率,而且可以集成到FDTD 求解器中,并与一般差分方法相兼容性.

全文的研究流程框架安排如下,第1 节为全文的引言,介绍关于深地探测低频问题建模求解过程中,显/隐式数值算法在求解过程中所面临的挑战,并通过引入当前机器学习技术克服相关的困难,从而引出本文的工作要点;其次,第2 节是对基于机器学习的完全匹配层吸收边界条件模型算法的提出、训练、预测及优化进行阐述;第3 节则是通过相关的深地探测问题去验证所提出算法模型的有效性及准确性,用于测试所提出模型算法的稳定性和精度是否达标;最后,对全文工作的总结.

2 基于XGB 的完美匹配单层模型

自从由Berenger[33]提出了完全匹配层吸收边界条件以来,各种完全匹配层实现算法已经用于截断时域有限差分网格.拉伸坐标完全匹配层[34]是将麦克斯韦方程组中的两个旋度方程映射到复数拉伸坐标空间,其优点是在完全匹配层的角和棱上的实现非常简单.与最初Berenger 的完全匹配层一样,拉伸坐标完全匹配层对隐失波的吸收无效[34].随后,复频率偏移完全匹配层[35]的出现备受关注,这是由于复频率偏移完全匹配层能够更有效地衰减低频隐失波以及减少电磁相互作用过程中后期的反射[36],其实现是通过简单地将复平面的极点从实轴移动到负虚半平面上实现的.在文献[36]中,基于拉伸坐标完全匹配层公式和卷积定理提出的卷积完全匹配层能够有效地实现复频率偏移完全匹配层.然而,复频率偏移完全匹配层对低频行波的吸收效果较差[37–39].因此,为了可以同时吸收行波和隐失波,双极PML 被提出并证实其具备同时吸收二者的性能.在本文中,所采用的就是双极PML 吸收边界条件进行区域截断并收集数据集,用于后期神经网络模型的训练.

PML ABC 通常被用来截断计算区域,从而得到一个开放域.一般来说,使用的是有限厚度层,比如8 层、10 层或更多层,这将导致额外的资源消耗,扩大了问题的复杂性.PML 最优参数(如α,κ和σ等)的选择主要凭借经验以获得满意的吸收精度,因此在试错时会耗费大量时间成本.为了克服这些问题,提出了基于XGB 的PMM 模型,从而取代传统的PML 方案进行计算区域的截断.

在第1 节中已经讨论过,基于XGB 的PMM模型在模型训练稳定性和模型轻量级方面具有明显的优势,其训练过程,如图1 所示.

图1 XGB 的训练过程Fig.1.Training process of XGB.

现在,对相关的表达式和符号作详细推导,以便能够清晰、正确地描述模型的细节.本文有N个样本,其数据集表示为X={x1,x2,···,xN},对应的目标值定义为Y={y1,y2,···,yN}.根据不同的增强模式,XGB 模型可以学习K个决策树,这些决策树的装配表示为{T1,T2,···,TN},从而训练一个新的决策树并将其添加到模型中以减少迭代训练过程中的损失值,如下所示:

对于在算法流程图中构建单个决策树的构建过程,只关注以下表达式:

YA,YB是通过将数据预测值以大小排序并随机取值将样品分为小于随机值和大于随机值的两部分:

(2)式的具体作用是指导决策树单个节点对数据集的划分取值,它在(2)式划分叶节点的具体方法是多次取不同数据预测值的中间值将数据集划分为不同的两组并用c1,c2作为两组数据的预测值,可以根据实际情况选择损失函数Loss,预测值c1和c2可以通过求解(3)式得到最优的结果.通过多次将数据一分为二,二分为四将构成类似树状结构,而c1和c2也将被计算的越来越精确,因为数据的规模越来越小,达到最小的一两个数据为一组时,几乎没有误差.

上述介绍决策树构建的简单流程,接下来将确定XGB 的目标函数来规定算法整体公式: 已知训练数据集T={(x1,y1),(x2,y2),···,(xn,yn)},损失函数,正则化项Ω(fk),则构造整体目标函数为

其中i是第i个样本,k是第k棵树.是第i个样本的预测值:

将上式代入目标函数并记一阶导数为g,二阶导数为h,得到最终结果为

正则化项由两个部分构成,叶子节点数量T与叶子节点权重向量ω,公式如下:

以此方式约束叶子节点的数量和权重范围.

在单棵树的分裂过程中,定义如下两个参数:

而针对叶子节点划分方式的最优计算则求解如下函数的最优解:

其中λ和γ为超参数,用于约束模型的复杂度和深度.

增强的迭代方法表明,在下一次迭代启动之前需要前一次迭代的残差,残差将作为下一个决策树的训练集用于训练下一个决策树,并直接与基础模型相加来直接降低模型的偏差,因此,XGB 模型选择负梯度近似残差,如下式所示:

(4)式显示了残差的计算方法,在每个单个决策树构建完成后,首先根据(4)式将单个决策树整合到整个模型中,然后用(5)式计算新的残差以获得新的残差数据M={m1,m2,···,mN}.

这些残差数据与样本数据集X={x1,x2,···,xN}一一对应.除了第1 个决策树是基于X和Y之间的对应关系所构造之外,后续的决策树都是基于X和M之间的对应关系构造的.

在传统的开域问题的FDTD 模拟中,有限的PML 层(如8 或10 层)通常用于终止物理域,这导致模拟过程需要大量的CPU 时间和内存来处理吸收边界.为了解决这两个问题,本文将重点介绍用基于XGB 的PML 模型来开发和实现FDTD域的终止.

本节构建了一个基于多个决策树和注意机制的XGB 模型,如表1 所列.一般来说,XGB 是一个集成学习方法,通过迭代来训练决策树和注意样本与之前的模型预测误差,在整个仿真期间,根据当前及前一时间步内的局部和邻近场所取得的局部场的循环计算过程可以避免.在前面的介绍中,已经详细说明了XGB 的具体逻辑,其中损失函数的选择和一些关键的超参数是至关重要的.结果表明,当单个决策树的最大深度为10,最小样本的分割数为8,单个采样比为90%并采用均方误差作为损失函数时,模型的性能更加理想.

表1 基于XGB 的ABC 模型算法Table 1.XGB-based ABC algorithms.

为了进一步验证所提出的基于XGB 的ABC模型,进行了数值计算.如图2 所示,三维模型的几何形状分别由空气、岩石和矿体组成.该矿体是一个深度为536 m,宽度为536 m,高度为536 m的立方体,且其中心位于(2010 m,2010 m,670 m),其电导率为4 S/m,岩石的电导率为0.005 S/m.空间用的均匀网格在三个方向均为67 m 离散化.仿真是在4020 m × 4020 m × 1876 m 的区域下完成的.

图2 三维ATEM 问题的几何模型Fig.2.Geometry of 3D ATEM Problem.

本文采用的样本量共为295200 组,每次采样包括15 个元素,分别是ex_it1,ex_it2,ex_it3,ex_it4,ey_it1,ey_it2,ey_it3,ey_it4,hz_it1,hz_it2,hz_it3,hz_it4,ex_ot,ey_ot,hz_ot.这些数据被用于学习一个新模型,用于预测电场与磁场的数值.每组所耗时332.06 s,共耗时近27.23 h.文中使用了295200 个样本作为训练数据,训练集的样本量为20 万,验证集和测试集的样本量为47600 个.该模型大约需要81.3 s 来训练60 次迭代,而预测295200 测试集样本只需要0.298 s.虽然模型在训练过程中只能进行连续训练,但在推理阶段也可以并行执行,以显著提高速度.该模型的内存使用量为4.08 MB.

在初始模型的设计中,我们的目的是设计一系列的实验来确定模型的三个重要的超参数,包括树的最大深度,叶片的最小样本分裂数,以及单个样本的采样率.如图3 所示,实验结果表明,更深的树状结构可以更好地学习数据之间的潜在相关性,因此,我们选择的最大深度为8.为了确保模型具有更强的泛化能力,在这里选择了90%的单样本采样率.

图3 训练过程的精度 (a)树深度;(b)叶子的最小样本分割数;(c)单个样本的抽样比Fig.3.Accuracies of training process: (a) Tree depth;(b) minimum number of sample splits for the leaves;(c) sampling ratio for a single sample.

用上述的超参数对所提出模型进行训练,然后应用该模型来预测样本的Ex,Ey和Hz.通过预测值与真实值之间的相对误差e来衡量精度A,如下式所示:

从图4 可以看到,以(6)式计算的指标作为参考,所提模型的Ex,Ey,Hz的精度分别为97.6%,98.7%,98%.XGB 模型长期以来广泛应用于各个领域,但在数值分布稳定均匀的数值模拟领域仍占有一席之地.在未来,还可以通过XGB 的特定调试来实现更高的精度和效率.

图4 训练过程中Ex,Ey,Hz 的相对误差Fig.4.Relative errors of Ex,Ey,and Hz during the training process.

3 数值实例及讨论

为了验证所提出的基于XGB 模型的准确率和效率,考虑了两种不同源的更复杂的ATEM 问题,分别为案例A——接地线源;案例B——机载磁偶极子源.图5 中的“案例A 和案例B”代表本节中的两个不同的例子.案例A 中,接地线源Jx被置于岩石上方,案例B 中,磁点偶极子源Mx位于空气中.

图5 三维ATEM 问题的几何模型Fig.5.Geometry of 3D ATEM problems.

从图5 可以看到,信号源是由重复频率为25 Hz的双极方波脉冲所激发产生.图6 展示了一个周期(即0.04 s)的脉冲波形.

图6 双极方波脉冲,持续时间为0.04 s,上升时间为0.0025 s,下降时间为0.0125 sFig.6.Bipolar square wave pulse with duration of 0.04 s,rise time 0.0025 s,and fall time 0.0025 s.

案例A地面上的线源和目标

案例A 中进行模拟的水平电流线源Jx位于地表上,这是典型的大地航空瞬变电磁,如图5 中的“案例A”所示.线源Jx从(1943 m,2010 m,1474 m)到(2077 m,2010 m,1474 m),其几何外观分别由空气、岩石和矿石组成,如图5 所示.每个矿山都是深度为536 m、宽度为536 m、高度为536 m 的立方体,其电导率σore为4 S/m,岩石的电导率为σrock为0.005 S/m.整个空间用 ∆x=∆y=∆z=67 m的均匀网格离散化.物理域的深度为4020 m、宽度为4020 m、高度为1876 m.

如图5 所示,为了使结构复杂,在岩石层中嵌入了3 个矿体,以验证所提出的基于XGB 的模型的鲁棒性和有效性.接下来,用五种方法对案例A进行分析,这些方法包括具有10 层PML 截断的传统FDTD 方法,具有5 层CFS-PML 截断的隐式FDTD 方法,内嵌复频率因子的深度可微森林完全匹配单层模型,内嵌复频率因子的深度决策树神经图灵机模型,及内嵌双极因子的梯度增强决策树完全匹配单层模型.

如图7(a)所示,在窗口的水平轴上设置时间长度为0.16 s 以观察接收机所在位置的二次磁场(即矿体目标的总场Hz减去无矿石的主场Hz0)的稳定性.从图7(a)和表2 可以看出,所提出的基于XGB 的PMM 模型不仅比基于DDF 的PMM模型和基于DDT 的NTM 模型有着更高的精度,而且计算效率比传统的隐式FDTD 方法更快.此外,在基于XGB 的PMM 模型中,只考虑了单层PML,从而降低了问题的复杂性.

表2 以0.04 s 为周期的计算时间Table 2.Computational time for 0.04 s as a period.

图7 对于案例A,利用五种不同方法所取得的二次磁场 (a)时间轴上二次场数值解与预测解对比;(b)相对反射误差计算Fig.7.Secondary Hz achieved by five different methods for case A: (a) Comparison of secondary Hz field between numerical methods and machine learning methods;(b) computation of relative reflection errors.

为了保证吸收性能可以满足工程问题的要求,进行了相对反射误差与时间关系的研究.从图7(b)可以看到,具有10 层PML 截断的传统FDTD 方法,具有5 层CFS-PML 截断的隐式FDTD 方法,内嵌复频率因子的深度可微森林完全匹配单层模型,内嵌复频率因子的深度决策树神经图灵机模型,及内嵌双极因子的梯度增强决策树完全匹配单层模的最大相对误差分别为–47.49 dB,–46.02 dB,–35.68 dB,–36.59 dB 和–41.59 dB.最后,结果表明,与基于DDF 的PMM 模型和基于DDT 的NTM模型相比,所提出的基于XGB 的PMM 模型达到了工程应用的要求(低于–40 dB).

从案例A 可以得出如下结论: 1) 所提出模型可以实现良好吸收精度类似于传统FDTD 方法和隐式FDTD 方法;2) 所提出模型仅采用单层PML,可以大大降低问题的复杂度;3) 所提出模型与传统的FDTD 方法和隐式FDTD 方法在二次磁场Hz上可以得到良好的一致性;4) 所提出模型的计算时间比其他方法少,尤其是传统的FDTD方法;5) 所提出模型能获得更好的数值精度和吸收精度.

案例B机载磁偶极子源和目标

为了进一步验证所提出模型的鲁棒性,对案例B 也进行了研究.如图5 所示,案例B 的配置与案例A 的配置相同.磁点偶极子源Mx位于表面上方,由一个重复频率为25 Hz 的双极方波脉冲激发.

由图8(a)可知,所提出的基于XGB 的PMM模型不仅与传统的FDTD 方法和隐式FDTD 方法有很好的一致性,而且在数值精度方面优于其他机器学习方法,这与案例A 相似.此外,通过对相对反射误差与时间关系的研究来观察这些方法的吸收精度.如图8(b)所示,所提出的基于XGB 的模型的最大相对误差与传统的FDTD 方法和隐式FDTD 方法相似,这些方法可以满足电磁工程问题中对吸收精度要求的最低标准.

图8 对于案例B,利用五种不同方法所取得的二次磁场 (a)时间轴上二次场数值解与预测解对比;(b)相对反射误差计算Fig.8.Secondary Hz achieved by five different methods for case B: (a) Comparison of secondary Hz field between numerical methods and machine learning methods;(b) computation of relative reflection errors.

为了保证吸收性能可以满足工程问题的要求,进行了相对反射误差与时间关系的研究.从图8(b)可以看到,具有10 层PML 截断的传统FDTD 方法,具有5 层CFS-PML 截断的隐式FDTD 方法,内嵌复频率因子的深度可微森林完全匹配单层模型,内嵌复频率因子的深度决策树神经图灵机模型,及内嵌双极因子的梯度增强决策树完全匹配单层模的最大相对误差分别为–47.26 dB,–47.02 dB,–35.29 dB,–36.98 dB 和–41.98 dB.最终结果表明,与基于DDF 的PMM 模型和基于DDT 的NTM模型相比,所提出的基于XGB 的PMM 模型符合工程应用的要求(低于–40 dB).

4 结论

本文提出了一种基于XGB 模型可替代PMM方法来解决极低频ATEM 地下探测问题.XGB 技术可以获得较高的精度和效率,最终数值结果表明,所提出的PMM 模型不仅能在数值和吸收的准确率上与传统的FDTD 方法和隐式FDTD 方法达到较好的一致性,而且在很大程度上降低了问题的复杂度,从而节省了更多的计算资源.此外,由于将基于Z变换方法的TP-PML 方案纳入基于XGB的PMM 模型中,使得系统的求解精度得到了进一步的提高.鉴于此,可以看到,所提出的模型可以成功地集成到FDTD 求解器中来解决ATEM 问题的建模.