基于邊界積分法的艦艇磁場模型參數(shù)確定

基于邊界積分法的艦艇磁場模型參數(shù)確定

0重點難以確定的磁場范圍難以延拓的歸位方法或原理消磁是提高潛水員耐磁性素質(zhì)的有效手段。評價消磁效果或者評估潛艇的磁性防護能力,一般都是針對潛艇水線以下的某種特定深度上某些特定點的磁感應(yīng)強度或潛艇上方的特定高度上的磁信號來進行的。受實際測量條件(水深,風浪和高度等)的限制,在實際測量中往往是很難直接獲得這類特定的磁場信息,所以需要從已知場點的磁場值來推算考核點上的磁場值。目前,用于艦艇磁場延拓的方法主要有有限元法、積分方程法和磁體模擬法。有限元和積分方程法都需要對艦船進行剖分,其剖分的合理性和疏密程度嚴重影響著延拓結(jié)果。對于形體復(fù)雜且龐大的艦船而言,剖分是非常復(fù)雜困難的,并且得到的方程組也將非常龐大,運算時間較長。磁體模擬法常用的磁性模擬體主要有磁偶極子和磁橢球體,該方法易受人為經(jīng)驗的影響,穩(wěn)定性和精度難以保證。人們在20世紀70年代將邊界積分方程用于對位場的研究,由于該方法只對求解區(qū)域邊界剖分,因而使問題的維數(shù)降低,這就使代數(shù)方程組的維數(shù)大為減少,節(jié)約了計算成本。1磁體擴展數(shù)學模型1.1各區(qū)域v內(nèi)磁場的各磁感應(yīng)強度計算如圖1所示,閉合面S將場域空間分成2部分V1和V,所有的磁源(包括鐵磁物質(zhì)和傳導(dǎo)電流)都被限制在V1內(nèi)。V可以看成是由2個閉合面限定的區(qū)域(從內(nèi)部限定的閉合面S和伸展到無限遠的S∞),該區(qū)域內(nèi)沒有電流分布(J=0),標量磁位φ滿足拉普拉斯方程??2φ=0.由于磁位φP在無窮遠處定義為0,則直接邊界積分方程可以寫成:cPφP=∮S14πr?φ?ndS-14πr?φ?ndS?∮Sφ??nφ??n14πr14πrdS,(1)其中:cP對V內(nèi)的點為1,對光滑邊界S上的點為1/2,對區(qū)域V1內(nèi)的點為0;r是邊界面上任意點到場點的矢徑;S面上的外法線n從區(qū)域V指向V1.記φ*=?φ?nφ?=?φ?n,則從(1)式可以看出:V1內(nèi)的磁源對V中任意位置磁位φP的貢獻可以用邊界上的面積分來表示。若已知邊界面S上的等效源磁性參數(shù)(φ和φ*)后,可以由該式計算V中任何一點的標量磁位φP.在得到該點的φP后,通過求取該點處的負梯度,并結(jié)合該點處的磁導(dǎo)率μ,即可得到該點處的磁感應(yīng)強度。BP=-μ??′φΡ=-μ4π′φP=?μ4π∮S??′1r?φ?n-φ??n1r?φ?n?φ??n1r1rdS,(2)其中,法線方向的偏微分是對邊界源點坐標進行的,而梯度??′是對場點坐標進行的。將邊界面S劃分為N個面積單元,則兩式可離散化為這就是邊界上的等效源磁性參數(shù)與無磁區(qū)域V中的磁感應(yīng)強度的正演關(guān)系。如果在區(qū)域V內(nèi)測量得到一組潛艇磁場感應(yīng)強度,就可以利用(4)式反演出邊界上等效的磁性參數(shù),再利用這組等效的磁性參數(shù)可以計算V區(qū)域內(nèi)任意位置的磁感應(yīng)強度。在反演的過程中,還需要考慮一些其他信息,如邊界面磁通連續(xù),φ和φ*的平衡,以及兩者在邊界上的連續(xù)性。1.2單元中心處的磁位平衡從測量值反演等效磁性參數(shù)時,約束條件為:1)φ和φ*的平衡;2)磁通連續(xù)。如果將場點取在各個單元的中心,則cp=1/2(為了計算方便,一般取光滑的邊界),φ和φ*在邊界上的平衡,可以描述為φi為第i個單元中心處的磁位。由于所有磁源都被限制在閉合面S內(nèi),進出S的磁力線條數(shù)相等。由高斯定理∮SBdS=0,經(jīng)推導(dǎo)并離散化后得到Ν∑j=1∫Sj∑j=1N∫Sj?φj?ndSj=0.(6)?φj?ndSj=0.(6)1.3型磁管磁模型的正演分析邊界形狀和邊界單元的劃分嚴重影響著正演和反演計算的速度和精度。出于單元劃分和積分計算簡單考慮,一般都選取較規(guī)則的邊界(如長方體表面,圖2).由于邊界S上沒有實際電荷和磁源,所以邊界上的φ和φ*連續(xù);同時邊界S存在拐點(如長方體的棱),在這些拐點處φ*的方向未定。線性單元適合于φ和φ*在邊界上連續(xù)的場合,但是如果將邊界全部使用線性單元分割,解在拐點處的誤差就大。為了改善解的精度,本文將邊界元視為混合單元,即在每個單元內(nèi)φ是連續(xù)線性變化的,而φ*是定常的,此時(4)式和(5)式相應(yīng)的變?yōu)槭街?Q為單元節(jié)點數(shù);ωkiki和ωkjkj為加權(quán)系數(shù);φkiki和φkjkj為單元i和j的第k個節(jié)點處的磁位。ωkiki和Q由邊界面的具體形狀和分割方式?jīng)Q定,一般采用簡單的四邊形元素(本文中為長方形,見圖2),它由4個節(jié)點組成,此時Q=4,ωkiki可按下式計算:{ω1i=(1-ε1)(1-ε2)4,ω2i=(1+ε1)(1-ε2)4,ε3i=(1+ε1)(1+ε2)4,ε4i=(1-ε1)(1+ε2)4,-1≤ε1,ε2≤1.(9)?????????ω1i=(1?ε1)(1?ε2)4,ω2i=(1+ε1)(1?ε2)4,ε3i=(1+ε1)(1+ε2)4,ε4i=(1?ε1)(1+ε2)4,?1≤ε1,ε2≤1.(9)(7)式和(8)式中的積分項與場點Pi和邊界單元Sj有關(guān),積分方法與單元的具體形狀密切相關(guān),一般難以求得解析表達式,都采用數(shù)值積分方法。積分項的具體表達式可以參閱文獻,這里不再贅述。對于(8)式中,應(yīng)當注意的是:當i=j時,源點與節(jié)點在同一單元上,出現(xiàn)積分奇異,兩者的積分結(jié)果如下:(11)式只適用于圖2中的長方形邊界單元,其中l(wèi)和d為長方形的長和寬。至此,作為延拓潛艇磁感應(yīng)強度所需的第1個步驟,我們完成了正演問題的分析。延拓的具體過程可以描述為:在潛艇周圍空間測量得到m個不同位置的潛艇磁場三分量,則綜合(6)～(8)式形成方程組,寫成矩陣的形式為:Ax=b;求解方程組Ax=b反演得到邊界上一組等效的磁性參數(shù)后,再利用(7)式計算艦船在其他位置的磁場感應(yīng)強度。2模型求解和正則化參數(shù)的選取由上節(jié)模型得到的方程組都呈現(xiàn)嚴重的病態(tài)特性,磁感應(yīng)強度測量值的微小擾動會導(dǎo)致解的很大變化。而在求解方程組時,通常按線性最小二乘來求解,即求解法方程組ATAx=ATb.(12)若方程本身是病態(tài)的,則法方程會更加的病態(tài)。國內(nèi)外很多學者都對這類病態(tài)問題進行了研究,這些研究表明TSVD和Tikhonov正則化等方法可以較好的處理逆問題中由系數(shù)矩陣病態(tài)特性而造成結(jié)果的不穩(wěn)定性。本文將采用Tikhonov正則化方法對方程組進行求解。2.1Tikhonov正則化基本原理不妨假設(shè)實際磁場測量過程中測量誤差為ξ,理想情況下磁場測量值為c,則有b=c+ξ,Ax=c+ξ.對系數(shù)矩陣進行奇異值分解A=USVT,則(12)式的解為x=k∑i=1∑i=1k(uΤic)σivi+(uTic)σivi+k∑i=1∑i=1k(uΤiξ)σivi,(13)其中ui和vi分別為U和V矩陣的第i列,σ1≥σ2≥…≥σk為A的奇異值,k=rank(A).分析(13)式可知,噪聲對真解的污染主要反映在右端第2項,當奇異值σi異常小時,該項將被嚴重放大,甚至掩蓋真解(右端第1項)。TSVD方法把容易造成不穩(wěn)定的較小奇異值直接截去,反映在物理實際中就是將該單元內(nèi)的磁化參數(shù)置為0,加劇了磁化參數(shù)的不連續(xù)性。而Tikhonov方法的處理思路:在奇異值上疊加一個極小的正實數(shù)α(正則化參數(shù)),抑制異常小的奇異值σi對噪聲的放大作用,達到正則濾波的目的。其相應(yīng)的解為xα=k∑i=1σi(uΤic)σ2i+αvi+k∑i=1σi(uΤiξ)σ2i+αvi.(14)從(14)式可以看出,若正則化參數(shù)α選取較大(過正則化),則xα與真解差別較大,模型的精確性降低;若正則化參數(shù)α選取較小(欠正則化),則無法抑制對噪聲的放大作用,模型的穩(wěn)定性下降。所以恰當?shù)恼齽t化參數(shù)需要在模型的精確性和穩(wěn)定性之間做出一個合適的折中。2.2正則化參數(shù)的選取正則化參數(shù)的選取方法分為兩類:先驗選取法和后驗選取法。先驗選取方法有Morozov偏差原理,該類方法需要知道與測量數(shù)據(jù)相關(guān)的先驗信息(如測量數(shù)據(jù)的誤差水平),然而這些先驗信息在潛艇磁場反演建模中非常少。后驗選取方法主要有L曲線法,交叉校驗(CV)和廣義交叉校驗(GCV).這3種方法直接從測量數(shù)據(jù)出發(fā),通過預(yù)解方程組來確定正則化參數(shù),而不需要先驗信息,適合本文中逆問題的求解。1)L曲線法如果以‖xa‖為橫坐標,‖Axa-b‖為縱坐標繪成的曲線為L形狀,稱為L曲線。L曲線反映了逆問題模型的估計參數(shù)的精確性和穩(wěn)定性隨著正則化參數(shù)α的變化規(guī)律。一般選取L曲線上曲率最大的點(都為L曲線的拐點)所對應(yīng)的α作為正則化參數(shù),但是磁場延拓模型實際的L曲線在拐點附近集聚了大量的α,最優(yōu)的α值難以確定。所以本文借助極小化如下輔助函數(shù)來選取α.ρ(α)=lg(‖Axα-b‖)+lg(‖xα‖),(15)磁場延拓模型的L曲線和輔助函數(shù)ρ(α)曲線見圖3和圖4.2)CV和GCV定義如下的CV和GCV函數(shù):CVα=1nn∑i=1bi-?bi1-Ηii(α)2,(16)GCVα=1nn∑i=1bi-?bi1-tr(Η(α))/n2,(17)式中:Hα=A(ATA+αI)-1AT;■=Hαb;tr(H(α))為H(α)的跡;Hii(α)是矩陣H(α)的對角元。使函數(shù)CV或GCV取得最小值的α就是由交叉校驗或廣義交叉校驗所求得的正則化參數(shù)。潛艇磁場延拓模型的CV函數(shù)和GCV函數(shù)曲線見圖5和圖6.從圖4～6可以看出,CV函數(shù)確定的α最小,其次是GCV函數(shù),最大的是L曲線。大量分析發(fā)現(xiàn):CV函數(shù)確定的模型精度較高,但魯棒性較差;L曲線確定的模型魯棒性較強,但精度較差;綜合比較之下,GCV函數(shù)選取的正則化參數(shù)比較恰當。本文在船模實驗中,采用GCV函數(shù)確定正則化參數(shù),取得了較為滿意的效果。3船模右病z1-12模以一個長308cm、寬34cm的鐵質(zhì)結(jié)構(gòu)船模為實驗對象,采用三分量磁通門傳感器陣列進行測量。在4m×0.4m的測量范圍內(nèi),以等間距分別測量得到21×5點陣列上的磁感應(yīng)強度,如圖7所示(船模中心在測量平面上的投影與平面中心重合,測量平面高度分別為z1=57.8cm和z2=72cm).閉合邊界面取包含船模的長方體表面,測量平面較其在水平面的投影稍大,整個邊界共劃分為280個單元。分別利用z1(z2)平面上的磁感應(yīng)強度三分量的測量數(shù)據(jù),建立磁場模型,并由模型計算得到了z2(z1)平面磁感應(yīng)強度的三分量值。定義相對殘差:‖Bm-Bc‖2/‖Bm,max‖2,推算誤差:Bci/Bm,max×100%,其中Bm和Bc分別為磁感應(yīng)強度的測量值向量和推算值向量,Bm,max為磁感應(yīng)強度測量值的最大值,Bci為第i點磁感應(yīng)強度的推算值。從表1可以看出:X分量的相對殘差最小,Z分量其次,Y分量最大,這是因為在Y方向的磁場信息比X方向和Z方向少造成的。圖8和圖9分別給出了近場到遠場及遠場到近場的推算值與測量值的對比曲線(由于篇幅的原因,只畫出潛艇龍骨下的磁場比較圖,其他位置有同樣的性質(zhì))。從圖中可以看出,由近場推算遠場時,3個分量的推算誤差均可控制在6%以內(nèi),而從遠場推算近場模型誤差稍大,但均可控制在9%以內(nèi)。圖10和圖11分別給出了經(jīng)正則化處理和未經(jīng)正則化處理2種情況下,從近場到遠場和從遠場到近場,船模右舷z分量的推算值與測量值的對比曲線(其他位置和其他分量具有同樣的性質(zhì))。從圖中可以看出,采用GCV函數(shù)正則