基于有限元邊界元法的艦船開域空間靜磁計算研究

左 超，耿 攀，夏建超，徐正喜

（武漢第二船舶設計研究所，武漢 430064）

0 引言

艦船磁場特征分析是艦船消磁的關鍵技術之一，隨著艦船消磁技術的戰術需求從單純的水下防御向著立體化三維防御發展，艦船的遠區磁場計算越來越受到科研人員的重視[1,2]。國外消磁工作者通過磁特征的基礎理論研究并結合磁場模型的仿真技術，使得磁場建模尤其是遠區磁場建模技術得到了重大突破[3]。國內艦船磁場計算方法也逐步從傳統的磁偶極子、磁偶極子陣列、旋轉橢球體陣列等算法向著基于現代計算機技術的微粒群、邊界元、有限元等算法發展[4]。

有限元算法適合處理非線性、多介質問題，但隨著求解域的增大，其計算量成倍增加，對于開域無限空間的求解并不適用。邊界元法由于具有降維、高精度，又能良好適應無限域、半無限域和奇異場的要求，所以可以有效適用于開域問題[5]。本文采用有限元－邊界元混合方法，將整個無限區域劃分為有限元區域和邊界元區域，通過艦船空間靜態磁場的計算仿真，驗證了該方法的正確性、高精度及高效性。

1 艦船無限空間磁場計算方法

1.1 有限元算法

艦船在地磁場作用下會造成空間磁場畸變，可以看作是靜態弱磁場對鐵磁物質的靜態可逆磁化過程。對艦船外部空間磁場分布的研究可以歸結為一個開域空間的靜態磁場問題。

在靜態磁場問題的求解中可以引入標量磁位φ，磁場強度表示為

其中：-gradφ是的無旋分量（ 即rot(-gradφ)＝0），而為有旋分量或為零。有限元法靜磁場標量模型的通用求解方程如下

有限元法就是將求解區域細分為大量的有限單元區域（如三角形），在每個單元區域中建立（2）式中方程。可以通過加勒金法或直接由目標能量尋找法獲得定義域內場的分布函數，此函數根據單元的節點取極值，然后將所有同節點有關的單元積分，給出一個設定節點磁勢（法線方向）多項式變分。這樣，每個節點產生一個等式，從而得到大量的聯立方程。解出這些方程，就可以確定每個節點的磁勢。對于磁特征的求解問題不可避免地要采用三維場的方法，但是采用標量磁勢方程可以大大地減少未知參數的數量，因為每個節點只需一個變量，而不是像矢量磁勢方程那樣需要三個變量[6]。由于艦船幾何模型能很全面的反映艦船實際構成，因此采用有限元法可較準確地計算艦船近區磁場。

艦船近區磁場有限元計算流程依次是物理模型的選取，建立艦船的幾何模型、網格剖分、設置物理屬性，其中包括材料屬性、機械屬性、外部電路、再對求解進行相應的算法選取和參數設置，求解完后就進行相應的后處理。

1.2 邊界元算法

邊界元法的基本原理是根據場源周圍閉合曲面的外法線方向的方向導數求得閉曲面上的標量磁位，進而推算出閉曲面外周圍空間各點的標量磁位分布，進而可以推算出空間各點的三分量磁場值。假設有一閉合面S將磁場之“源”完全包于其內，如下圖1所示艦船包絡面示意圖，稱S為邊界面，S之外的空間為場域V。

圖1 艦船包絡面示意圖

為研究有限場域內的磁場分布，設u為標量磁位，如果在邊界面S上，u或（u的法向導數）是已知的，則場域中任意一點的u都是唯一存在的。由邊界上u或的分布求場域中u的問題稱為磁場邊值問題。設艦船包絡面外的區域V(即遠區換算點)為待換算的區域，艦船包絡面S為V的邊界面，設V中任一點為M，包絡面上任一點為M0，則由格林第二公式，可推導出M0處的磁位值為：

式中r為M0至變點M的距離，為沿邊界面外法線方向的方向導數。

2 仿真算例

本文以電磁場仿真軟件FLUX為基本的仿真平臺進行艦船開域空間磁場計算，法國CEDRAT公司的FLUX軟件被廣泛應用于艦船磁場的計算及消磁繞組布置的設計[8,9]。FLUX封閉域求解基于有限元算法，可進行艦船的復雜結構、薄殼超大縱橫比、多介質界面耦合仿真計算，同時還可以在無限區域耦合多類型的邊界條件，以求解艦船遠區磁場。

第一步，根據計算需求，選擇FLUX軟件的靜態磁場三維求解器Magneto Static 3D進行艦船靜態磁場仿真計算。

第二步，通過造型系統對艦船的主要殼體部件進行幾何模型建模，由于曲面造型難以進行網格剖分，該模型對艦船曲面進行了大量簡化，肋骨、設備、上層建筑也進行了簡化處理，FLUX建模最大的特點是對于艦船殼體部分的簡化處理可以采取蒙皮模型，而非實體模型，該方法的好處是可以減少剖分網格量，也方便求解薄殼模型。如圖2所示，簡化艦船模型長100米，寬20米，高10米。可以看到艦船位于雙層空氣包以內，空氣包外壁邊界用于模擬開域問題來求解無窮遠處的磁場問題。

圖2 簡化艦船幾何模型

第三步，模型的網格剖分。艦船形狀特殊，結構復雜，因此在剖分的過程中盡量使用只有四面體單元的剖分器，但不宜采用自動剖分模式，手動設定合適的網格大小，是保證求解精度的關鍵。如圖3所示，針對本文模型大小，剖分網格邊長設定為空氣包層 5000 mm，船模部分 1000 mm。為了滿足計算精度要求，剖分的過程中要保證雙層空氣包之間達到有效的網格量（大于 3層剖分單元）[10]。

圖3 艦船網格剖分圖

第四步，設置物理屬性。簡化艦船模型中的船體材料設置為兩種磁導率，殼體部分相對磁導率設定為100；甲板部分相對磁導率設定為300，厚度均為50 mm；其余部分全部設定為真空區域。為了方便分析磁場結果，外加磁場激勵設定為地磁單一分量，即武漢地區水平磁場強度34000 nT，空間均勻分布，艦船航向為北航向。

第五步，求解及后處理。經過計算求解完成之后，即可進行后處理，如圖4所示，可查看艦船表面的磁場分布云圖。此外還可以在遠區空間增加后處理測試線/面，查看開域空間靜磁分布。

圖4 艦船表面磁場分布云圖

3 結果分析

為了驗證有限元邊界元混合方法在求解艦船開域空間磁場過程中的高效性，我們使用同樣的仿真模型，比較單純使用有限元和有限元邊界元混合方法求解過程中的計算資源消耗，如表1所示，有限元邊界元混合方法可將整個求解時間降低81%，極大提高艦船遠區靜態磁場的求解效率。

表1 兩種算法求解效率對比表

同時，為了驗證有限元邊界元混合方法在求解艦船開域空間磁場過程中的精度，我們建立了磁偶極子模型，將上述艦船模型等效為11個縱向均勻分布的磁偶極子陣列，通過100 m高度數據求解偶極子磁矩，然后推算其他高度磁場強度。求解結果對比如下圖5～圖8所示。

圖5 FLUX仿真與磁偶極子結果對比（上方50m）

從圖5結果可以看出，有限元邊界元混合方法能夠在近場區域保證更為豐富的艦船磁性特征，相比磁偶極子算法精度更高。從圖6-圖8結果可以看出，有限元邊界元混合方法同樣能夠在遠區保證可靠的求解精度，與磁偶極子陣列模型相比，誤差控制在2%以內。

圖6 FLUX仿真與磁偶極子結果對比（上方100m）

圖7 FLUX仿真與磁偶極子結果對比（上方500m）

圖8 FLUX仿真與磁偶極子結果對比（上方1000m）

4 結束語

本文分析了有限元邊界元混合法在艦船開域空間靜態磁場求解問題中的適用性，并以電磁場分析軟件FLUX作為仿真計算平臺，通過有限元邊界元混合法對艦船開域空間靜態磁場進行了仿真計算，建立了艦船幾何模型，進行了網格剖分，并深入探討了參數設置方法，得到了仿真計算結果。在艦船開域空間靜態磁場求解問題上，該方法具備高效、高精度等優點，可應用于工程實踐。

[1]粟有鼎.艦船磁防護與聲防護[M].北京: 國防工業出版社,1985,97-99.

[2]劉大明,劉勝道,肖昌漢,周國華,王昭.艦艇閉環消磁技術國內外研究現狀.船電技術,2011,31(10):6-9.

[3]Jeung G,Yang C S,Chung H J,et al.Magnetic dipole modeling combined with material sensitivity analysis for solving an inverse problem of thin ferromagnetic sheet[J].IEEE Transactions on Magnetics,2009,45:4169-4172.

[4]連麗婷,肖昌漢,劉勝道.基于微粒群算法的薄鋼板兩側磁場推算中的位置優化[J].上海交通大學學報,2010,44(7): 975-979.

[5]王占輝,高俊吉.一種開域靜磁場雙標量位混合有限元邊界元法研究[J].船電技術,2013,33(06):19-21.

[6]汪家驊.正則化技術在艦船磁場建模逆問題中的應用[J].水雷戰與艦船防護,2010,18(4): 1-5.

[7]倪光正,楊仕友,錢秀英,等.工程電磁場數值計算[M].北京: 機械工業出版社,2004: 237-241.

[8]Birsan M.User Experience: Defence research for development Canada[J].Flux Magazine,2006,50: 11.

[9]姜智鵬,莊飚,關濤,母海方.基于FLUX軟件的艦船磁場計算[J].船電技術,2010,30(12): 17-20.

[10]張超,高立娥,劉衛東,劉鐸.FLUX(3D)在艦船磁場方面的應用研究[J].聲學技術,2011,30(4):110-113.