摘要
为了评估潜艇水下腐蚀相关稳恒电场分布特性,从潜艇水下腐蚀相关稳恒电场产生机理出发,基于等效点电流源建立潜艇水下腐蚀相关稳恒电场正演模型,并利用Tikhonov正则化根据已知电场数据求解等效点电流源电流强度,对潜艇周围海水空间的腐蚀相关稳恒电场进行推算。将某型潜艇的COMSOL软件仿真结果作为模拟试验数据,对所提方法有效性进行验证。结果表明:由水深38 m电场值向水深42.5 m和水深33.5 m进行反演时,相对均方根误差、最大值相对误差、峰峰值相对误差均不超过0.06;由较近平面向较远平面进行反演时,即使推算深度达到45 m,相对均方根误差仍然在0.21以内;噪声标准差为实际电场最大值的0.1 论文拓展倍时,反演误差仍然小于0.1。该算法抗噪声能力强,精度较高,能较好地用于工程实践。
关键词
Abstract
In order to evaluate the distribution characteristics of submarine underwater corrosion related steady electric field, starting from the mechanism of submarine underwater corrosion related steady electric field, a forward model of submarine underwater corrosion related steady electric field was established based on the equivalent point-type current source theory. And the current intensity of the equivalent point-type current sources was solved by Tikhonov regularization according to the measured electric field data. The effectiveness of the proposed method was verified by commercial finite element software COMSOL. The results show that when inverting the electric field value from depth of 38 m to a water depth of 42.5 m and 33.5 m, the relative root mean square error, maximum relative error, and peak to peak relative error do not exceed 0.06. Even if the calculated depth reaches 45 m, the relative root mean square error is still less than 0.21. When the noise standard deviation is 0.1 times the maximum value of the actual electric field, the inversion error is still less than 0.1. This algorithm has strong noise resistance and high accuracy, can be used in engineering practice.
潜艇水下腐蚀相关稳恒电场是一种重要的物理场,不仅可用于对潜艇进行探测,也可作为水雷等水中武器的引信作用源,对潜艇水下腐蚀相关稳恒电场特征进行评估[1-3]。而由于水下电场测量装置的限制,以及海洋环境的复杂性,在实际测量中往往只能获取潜艇下方某一路径或某个特定区域上的电场数据,无法得到潜艇周围整个海水区域的电场分布特征。虽然有限元法、边界元法等数值方法可以求解出复杂结构潜艇的水下腐蚀相关电场分布,但是需要已知潜艇外形结构以及艇体材料的极化曲线,且要进行复杂的网格剖分,计算量非常大[4-7]。因此,实际情况中往往需要根据有限区域、有限深度的实测电场对潜艇水下腐蚀相关稳恒电场进行反演,得到潜艇周围海水空间不同深度、不同正横距的腐蚀相关稳恒电场分布特性。
目前,对舰船腐蚀相关稳恒电场进行反演主要有物理方程法和等效源法两种思路。文献[8-9]从舰船腐蚀相关稳恒电场的物理特性出发,基于拉普拉斯方程分别提出了舰船腐蚀相关稳恒电场反演的拉氏方程法和格林函数法。该类方法无须对场源进行求解即可实现舰船腐蚀相关稳恒电场深度换算,但只适用于深海环境,且只能由较浅平面向较深平面反演。文献[10]基于舰船腐蚀稳恒电场的物理特性,利用微分递推法对舰船水下电场进行反演,可以实现由近及远和由远及近的舰船水下电场深度换算,但是要测量某一平面的电场值,在实际中难以实现。文献[11]以水平电流线和电偶极子作为舰船水下腐蚀相关电场模拟源,并基于最小二乘法利用实测电场数据对模拟源进行反演,该方法对于深海、浅海环境都适用,但是需要求解的电偶极子未知量较多,导致计算复杂。文献[12]基于点电流源对舰船腐蚀稳恒电场进行建模,并利用先验信息对等效点电流源的位置和数量进行估算,减少了未知量的个数。但该方法需要的先验信息较多,且确定等效点电流源位置的过程较为复杂。
本文从潜艇水下腐蚀相关稳恒电场的产生机理出发,将潜艇表面腐蚀和防腐电流等效为若干个点电流源,并根据分层介质中点电流源在海水区域的电场表达式构建潜艇表面腐蚀电流与海水中电场值之间的线性方程组。然后利用Tikhonov正则化处理潜艇腐蚀相关稳恒电场反演问题的不适定性,基于部分已知电场值求解等效点电流源电流强度。利用数值仿真实验验证了所提方法的有效性,并分析了测量线位置、等效点电流源数量、测点数量和测量噪声对反演精度的影响。
1 腐蚀相关稳恒电场正演模型
假设潜艇水下表面γ的电流密度为Jcorr(γ),根据稳恒电场基本方程,潜艇在海水中场点P处产生的稳恒电场为:
(1)
式中,σw为海水电导率,r 和r′为分别为场点P和源点M的位置矢量。
由于潜艇结构复杂,且不同部位材料和涂层状态各不相同,潜艇表面的电流分布是不均匀的,可将潜艇表面γ划分为n个小面源,当小面源的面积足够小时,每个小面源的电流密度可近似为常数[12-13]。假设第i个小面源的电流密度为J′(γi),面积为Δγi,则小面源的等效电流强度I(γi)=J′(γi)Δγi。根据叠加原理,海水中场点P处的腐蚀相关电场等于所有小面源产生的电场之和,即
(2)
当场点P到潜艇的距离足够远时,可忽略每个小面源的形状,将小面源等效为点电流源,则潜艇腐蚀相关电场的近似表达式为:
(3)
在空气-海水和海水-海床界面处会发生电荷累积,从而对海水中的电场产生影响,因此需要对式(3)进行修正。
由镜像法可知,空气-海水-海床三层介质中,点电流源在海水中的电场可看作无穷多个镜像点电流源的叠加,如图1所示。
图1空气-海水-海床三层介质中镜像点电流源示意图
Fig.1Schematic diagram of mirror point current source in air seawater seabed three-layer medium
因此只要求出每个小面源等效点电流源的等效位置和电流强度,就可以求出潜艇在海水中任意场点Pj腐蚀相关稳恒电场强度。
(4)
式中:
(5)
式中:ri1k、ri2k、ri3k、ri4k为位于海水中(x0i,y0i,z0i)处的点电流源Ii及其在海水和海床中的镜像点电流源的位置矢量,坐标为(x0i,y0i,-2kD+z0i),(x0i,y0i,-2kD-z0i),(x0i,y0i,2(k+1)D-z0i),(x0i,y0i,2(k+1)D+z0i),其中k=0,1,2,···。ri1k当k=0时位于海水中,其余都位于空气中;ri2k均位于空气中;ri3k、ri4k均位于海床中。D为海水深度;c=(σw-σb)/(σw+σb)为反射系数,其中σw和σb分别为海水电导率和海床电导率;k为反射层数,在工程实际中可取其上限10~20[13]。
2 腐蚀相关稳恒电场反演的正则化求解
2.1 反演方程组的建立
若用水下电场采集装置获取潜艇周围某条测线上的腐蚀稳恒电场三分量,则可得到m=3N个电场测量值(N为测点个数)。以电场测量值为已知量,等效点电流源电流强度为未知变量,可得到一个m维线性方程组,写成矩阵形式为:
(6)
式中,b=[E1x,E1y,E1z,E2x,E2y,E2z,···,ENx,ENy,ENz]T为所有测点电场三分量测量值组成的数据观测向量,x=[I1,I2,···,In]T为等效点电流源电流强度组成的待求参数向量,矩阵A为将场点和所有等效源点联系起来的系数矩阵。
(7)
在实际情况中,对潜艇表面进行合理划分后,测量点和等效点电流源坐标是已知的,由式(5)可得系数矩阵A,因此潜艇腐蚀相关稳恒电场反演问题的关键是根据实测电场值求解出等效点电流源的电流强度。
2.2 反演方程组的Tikhonov 正则化求解
在工程实际中,一般选取测量电场值数量m=3N大于待求等效点电流源的数量n,则式(6)为超定方程,存在最小二乘解xLS。
将系数矩阵A进行奇异值分解可得:
(8)
式中:U=(u1,u2,···,um)和V=(v1,v2,···,vn)为正交矩阵;对角矩阵Q=diag(q1,q2,···,qn)为系数矩阵A的奇异值矩阵,元素非负,且按从大到小的顺序排列。可得到式(6)的最小二乘法解向量表达式为:
(9)
然而,在潜艇电场反演问题中,系数矩阵A的条件数非常大,即A存在非常小的奇异值qi,因此较小测量电场误差也会使反演结果远远偏离真实值,即式(6)的最小二乘解是没有意义的。
而Tikhonov正则化技术是最经典和应用最为广泛的处理病态问题的有效手段之一,广泛应用于地球物理反演、磁场延拓各个工程领域的不适定反演问题[14-16]。其基本原理是在最小二乘拟合误差的基础上加上一个对解向量x进行约束的正则项,从而使正则化的解接近真实解并且稳定。正则化后的解向量为:
(10)
式中,L为正则化矩阵,λ为正则化参数,且λ>0。当L=In×n为单位矩阵时,为标准Tikhonov正则化,解向量x具有最小范数解。选择合适的λ,能使反演解在数据拟合最小误差和稳定性之间达到平衡。此时式(10)等价于求解超定方程式(11)的最小二乘解。
(11)
式(11)的解为:
(12)
从式(12)可以看出,Tikhonov正则化方法将对测量误差具有放大作用的小奇异值进行滤波处理,从而抑制测量误差对反演结果的影响。
2.3 正则化参数的选取
选择合适的正则化参数λ是Tikhonov 正则化反演结果好坏的关键。国内外不少学者对正则化参数的选取进行了研究。其中应用最是广泛的是广义交叉验证法(generalized cross validation method,GCV)、Morozov 原理和L曲线法。大量数值试验结果表明,采用L曲线选择的正则化参数在潜艇水下腐蚀稳恒电场反演中的精度和稳定性较好。
L曲线法是以正则化解范数的自然对数为纵坐标,相应的拟合误差范数的自然对数为横坐标的曲线,在拐角处L曲线的曲率最大,解范数和拟合误差范数能达到较好的平衡,如图2所示。
曲率根据式(13)计算[14]。
(13)
式中,,η′为η对正则化参数λ的一阶倒数。
图2L曲线示意图
Fig.2Schematic diagram of the L curve
3 腐蚀相关稳恒电场反演算例
由于有限元法、边界元法已经被广泛应用于潜艇腐蚀相关电场计算、潜艇阴极保护系统设计,且精度较高[4-7],因此利用商业有限元软件COMSOL对某型潜艇水下腐蚀相关稳恒电场分布进行仿真模拟试验,验证基于等效点电流源和Tikhonov正则化对潜艇水下腐蚀相关稳恒电场进行反演的有效性和实用性。
3.1 潜艇腐蚀相关稳恒电场数据
利用商业有限元软件COMSOL 5.4版本电化学模块建立某型潜艇三维模型。以空气-海水界面为xoy 平面,x轴正方向由艇首指向艇尾,y轴正方向由右舷指向左舷,z轴正方向垂直向下,艇首在xoy平面的投影点为原点,建立三维坐标系,如图3所示。
仿真模型包括艇壳、尾舵、螺旋桨、大轴、指挥台,潜艇长L=81.5 m,型宽B=9 m,型深9.2 m,最大高度为12.5 m(含指挥台),螺旋桨裸露,材料为铝青铜,艇体材料采用低合金钢,表面涂有防腐涂层,在潜艇中前部和后部各装有1对辅助阳极,沿潜艇中轴线对称分布。潜艇涂层局部破损,破损率为2.5%,前部辅助阳极对的输出电流为12 A,中后部阳极对的输出电流为16 A。艇壳材料腐蚀电位φK=-0.69 V,阳极极化率PKa=1.5 Ω·m2,阴极极化率PKc=27 Ω·m2,螺旋桨腐蚀电位φB=-0.369 V,螺旋桨阴极极化率PBc=0.23 Ω·m2。潜艇下潜深度为20 m,海水深度为150 m,海水电导率为4 S/m,海床电导率为0.1 S/m。利用COMSOL软件二次电流模块进行仿真可以得到整个海水域的电场分布。水深33.5 m和42.5 m的水下平面电场分布云图如图4所示。
图3潜艇模型和坐标示意图
Fig.3Submarine model and coordinate diagram
由图4可知,潜艇水下腐蚀相关稳恒电场具有明显的分布特性,整体上x分量和z分量沿艇体纵轴对称分布, y分量延潜艇纵轴反对称分布,且在水下33.5 m和42.5 m的深度上分别可达10-3 V/m和10-4 V/m量级。
图4潜艇水下腐蚀稳恒电场COMSOL仿真值
Fig.4COMSOL simulation value of underwater corrosion steady electric field of a submarine
3.2 点电流源位置和数量的选取
由于辅助阳极和螺旋桨面积远小于艇体,将4个阳极和螺旋桨各等效为1个点电流源。同时潜艇结构具有对称性,艇体表面等效点电流源延中轴线对称分布,且沿艇体纵轴等间隔分布,如图5所示。
图5等效点电流源位置示意图
Fig.5Schematic diagram of the location of equivalent point-type current sources
为了保证反演具有较高的精度,应满足潜艇表面小面源的等效半径相对于测量点P到潜艇中轴线的距离d足够小,可取等效半径ri≤d/4。另一方面,等效点电流源到潜艇表面的距离为B/4,当点电流源的纵向间隔l≤B/4时,影响反演精度的主要因素为等效点电流源在截面上的位置偏差,继续减小纵向间隔l会使计算量增大,但对提高反演精度的作用不大。因此,可令B/4≤l≤d2/(2B),即等效点电流源数量n满足:
(14)
式中,nk为螺旋桨和阳极的数量之和。需要说明的是,上述点电流源位置和数量的选取方法并不是为了对点电流源进行精确选择,而是一种满足等效计算条件下合理的假设与估算。并且为了提高等效点电流源的模拟精度,在满足式(14)的前提下可适当增加等效点电流源的数量。
3.3 反演结果与分析
以x轴范围为-100~200 m、y=-13.5 m、z=38 m路径上6 m等间距共51个测点的电场三分量作为测量值,对等效点电流源电流强度进行求解,并利用等效点电流源模型对潜艇水下电场进行模拟。根据式(14)可知,等效点电流源数量n可在21~150之间,取n=109。为了验证基于点电流源和Tikhonov正则化对潜艇腐蚀相关电场反演中由近及远、由远及近[13]和横向距离反演的有效性,以y=-13.5 m,z=42.5 m和y=-13.5 m,z=33.5 m以及y=-18 m,z=38 m这3条路径作为反演路径,3条路径上的电场反演值与对应路径上的COMSOL仿真值对比情况如图6所示。
图6中,Exc、Eyc和Ezc为COMSOL仿真电场值,Exi、Eyi和Ezi为反演电场值。由图6可知,在3条路径上,利用等效点电流源和Tikhonov正则化反演得到的潜艇腐蚀电场各分量与COMSOL仿真值具有相似的变化的趋势,且基本重合。
图6潜艇腐蚀电场反演值与仿真值
Fig.6Inversion value and simulation value of corrosion electric field of a submarine
为了定性分析反演计算误差,采用反演测量线上的电场反演值与COMSOL仿真值的相对均方根误差RE、最大值相对误差REmax和峰峰值相对误差REpp来描述。
(15)
(16)
(17)
其中,breg为Tikhonov正则化求解得到的电场值向量,bcom为COMSOL仿真得到的电场值向量。3条不同路径上的电场反演误差分别如表1、表2和表3所示。
表1y=-13.5 m,z=42.5 m电场反演误差
Tab.1 Electric field inversion error for y=-13.5 m, z=42.5 m
表2y=-13.5 m,z=33.5 m电场反演误差
Tab.2 Electric field inversion error for y=-13.5 m, z=33.5 m
表3y=-18 m,z=38 m电场反演误差
Tab.3 Electric field inversion error for y=-18 m, z=38 m
由表1、表2和表3可知,3条路径上的电场反演值和COMSOL仿真值相对均方根误差、最大值相对误差和峰峰值相对误差的最大值都不超过0.06,精度较高,说明本文反演算法对由远及近、由近及远和横向距离反演都是有效的。
为了进一步说明本文反演方法的有效性,分别以y=-13.5 m,z=29 m和y=-13.5 m,z=74 m为测量值,在二者深度范围内进行由近及远和由远及近的平面电场反演,电场反演相对均方根误差如图7和图8所示,图中ΔH=|z1-z2|表示测量深度z1和反演深度z2之差的绝对值。
图7由近及远不同深度电场反演误差
Fig.7Inversion error of electric field at different depths from near to far
由图7可知,由近及远进行反演时,电场Ex分量和Ey分量随着深度的增加,反演误差增大,当反演深度为45 m时,Ex、Ey分量误差分别增大至0.110 8和0.204 4;Ez分量反演误差先减小后增大,反演误差在0.012 8和0.038 3之间波动。
由图8可知,由远及近进行反演时,电场Ex分量随着深度的增加,反演误差先减小后增大,反演深度为18 m时有最小反演误差0.059 6;Ey分量和Ez分量随着反演深度的增加,反演误差增大,当反演深度为45 m时, Ey 、Ez分量误差分别增大至0.239 0和0.329 5。
文献[10]由近及远进行反演时,反演深度为11 m时反演误差超过0.5;由远及近反演时,反演深度为20 m时反演误差超过0.5。文献[11]由近及远进行反演时,反演深度为14.5 m时反演误差小于0.037;但由远及近进行反演时,反演深度为1.5 m时反演误差达到了0.125。对比来看,本文反演算法具有较高的精度。
图8由远及近不同深度电场反演误差
Fig.8Inversion error of electric field at different depth from far to near
使用本文算法时,总体上相同反演深度下由近及远进行电场反演的精度大于由远及近进行反演的精度,这是由于由远及近进行反演时,螺旋桨和涂层破损区域等效半径相对于场点到潜艇的距离变大,将涂层破损处和螺旋桨等效为点电流源的误差增大。而由近及远进行反演时,螺旋桨和涂层破损区域等效半径相对于场点到潜艇的距离变小,将涂层破损处和螺旋桨等效为点电流源的误差减小。
4 电场反演精度仿真分析
基于点电流源和Tikhonov正则化对潜艇水下电场进行反演时,除反演深度外,等效点电流源数量、测点数量和噪声对反演精度也有较大的影响。在3.1节的基础上,改变等效点电流源数量、测点数量和添加随机噪声来仿真分析相关因素对反演精度的影响。
4.1 等效点电流源数量对反演精度的影响
选取y=-13.5 m,z=33.5 m作为测量线,y=-13.5 m,z=42.5 m作为反演线,x轴范围为-100~200 m,测点间距为6 m。螺旋桨和阳极等效点电流源位置不变,改变艇体等效点电流源间隔以改变等效点电流源数量,等效点电流源数量分别为5、13、21、29、69、109、167、309时反演电场值和COMSOL仿真电场值的相对均方根误差如表4所示。
表4反演误差随等效点电流源数量的变化
Tab.4 Variation of inversion error with the number of equivalent point current sources
由表4可知,当等效点电流源数小于21时,电场Ex、Ey、Ez分量反演误差均在0.095以上,其中电场Ey分量反演误差超过0.15;当等效点电流源数大于29时,电场Ex、Ey、Ez分量反演误差明显减小,分别基本保持在0.04、0.04和0.02以下。点电流源数量为29、69、109时,反演精度差别不大,等效点电流源增加到167和309时,反演精度虽略有增加,但电场Ex、Ey、Ez分量反演误差仅分别减小约0.001 9、0.002 2和0.001 2,且计算量增大。而根据式(14)可知,等效点电流源数量n可在21~150之间,进一步说明采用3.2节的方法选取等效点电流源的位置和数量是合适的。
4.2 测点数量对反演精度的影响
选取y=-13.5 m,z=33.5 m作为测量线,y=-13.5 m,z=42.5 m作为反演线,x轴范围为-100~200 m。等效点电流源数量为109,通过改变测点间隔来改变测点数量,不同测点数量时电场反演的相对均方根误差如表5所示。
表5反演误差随测点数量的变化
Tab.5 Variation of inversion error with the number of measuring points
在反演问题中,为使式(6)为超定方程,测量电场数量m=3N应大于待求等效点电流源的数量n。而由表5可知,当测点数量为51,即测量数量153大于等效点电流源数量109时,随着测点数量的增多,Ex、Ey和Ez分量反演误差变化不大,Ex分量在0.042左右波动,Ey分量在0.046左右波动,Ez分量在0.020左右波动,Ey分量甚至随着测点数量的增多略有增大。考虑到反演精度和计算量,电场值测量数量m可以选择在等效点电流源数量的110%~150%之间。
4.3 噪声对反演精度的影响
为了检验基于Tikhonov正则化对潜艇电场进行反演的抗噪声能力,在COMSOL仿真电场值的基础上加上随机噪声作为测量值,噪声均值为0,标准差分别为电场各方向分量最大值的0.01、0.02、0.05、0.10、0.15和0.20倍时电场反演值和测量值的均方根误差,如表6所示。
表6反演误差随噪声幅值的变化
Tab.6 Variation of inversion error with noise amplitude
由表6可知,虽然随着噪声的增大,反演误差不断增大,但当随机噪声标准差为0.1倍最大电场值时,反演精度仍在0.1以内,即使噪声幅值达到0.2,反演误差仍然小于0.2,这说明本文反演算法抗噪声能力较强,能较好地用于工程实践。
5 结论
1)利用等效点电流源和Tikhonov正则化对潜艇周围海水区域的腐蚀相关稳恒电场进行反演,点电流源位置确定较为简单,且反演精度较高,抗噪声能力较强。
2)由近及远进行反演时,电场Ex分量和Ey分量反演误差随着反演深度的增加而增大,而Ez分量反演误差先减小后增大。由远及近进行反演时,电场Ex分量反演误差随着深度的增加先减小后增大, Ey分量和Ez分量反演误差随着反演深度的增加而增大。
3)优化等效点电流源数量和测点数量可以提高反演精度。适当增大点电流源数量可提高反演精度;兼顾反演精度和计算量,电场测量数量可选择在等效点电流源数量的110%~150%。反演误差随着噪声的增大而增大,但噪声标准差达到电场最大值0.2倍时,反演均方根误差仍然小于0.2,抗噪声能力较强。