海水包裹空心橢球體的低頻散射場計算及分析

徐龍飛 鄭奎松 張世田 李亞靜

(1.西北工業(yè)大學電子信息學院,西安 710072;

2.中國電波傳播研究所 電波環(huán)境特性及?；夹g重點實驗室,青島 266107)

引 言

隨著水下目標探測和海洋環(huán)境勘探的快速發(fā)展,水下目標的探測研究逐漸成為人們關注的焦點[1-2],特別是開發(fā)一種處理水下應用的準確算法是多年來一直的挑戰(zhàn). 在解決電磁散射問題的眾多計算方法中,時域有限差分(finite-difference time-domain, FDTD)方法由于其直觀的物理過程和簡單的計算格式,已經(jīng)發(fā)展成為計算水下低頻電磁波傳播的有效工具[3-4]. FDTD方法是將微分形式的麥克斯韋旋度方程進行差分離散從而得到一組時域推進公式,并在時間軸上逐步求解空間的電磁場. 對于一個三維仿真計算問題,FDTD方法要求將問題空間離散成許多Yee元胞. 假設Yee元胞的尺寸為δ,則時間間隔Δt的大小與元胞尺寸有關,受到Courant穩(wěn)定性條件的限制[5].

對于傳統(tǒng)的FDTD方法,在計算仿真低頻電磁波時,需要大量的時間步才能將仿真計算完成,計算時間非常長. 并行散射場源FDTD(tolal-field scattered-field source FDTD, TSS-FDTD)方法的提出為解決該問題提供了一種途徑,但有些情況下依舊需要非常長的計算時間[6]. 文獻[7]提出了一種在太赫茲頻段利用色散介質(zhì)特性提高計算效率的新思路. 該方法是通過改變色散介質(zhì)的控制參數(shù)來模擬不同媒質(zhì)材料[7]. 例如,設置模擬空間填充的背景材料為介電常數(shù)為80的海水,時間間隔的取值可大約設置為背景材料為空氣時的9倍,從而計算效率可大約提高9倍.

為了解決水下低頻電磁波波長與殼體目標厚度尺寸在電長度方面相差較大的電磁散射計算問題,本文提出了一種基于移位算子時域有限差分方法的加速計算方法.利用該方法計算了水下橢球殼體目標的電磁散射.數(shù)值結果表明該方法可以提高計算效率,縮短計算時間,為水下殼體目標的探測、識別和跟蹤提供一定的工程應用價值.

1 色散介質(zhì)參數(shù)的設置

本文采用的色散介質(zhì)FDTD方法的控制方程為

(1)

(2)

D(ω)=ε0εr(ω)E(ω)．

(3)

式中:式(1)和(2)是時域麥克斯韋旋度方程;εr(ω)為相對介電系數(shù),是頻率的函數(shù),故式(3)表現(xiàn)為頻域的形式.

考慮有耗Drude模型的色散介質(zhì),其相對介電系數(shù)為

(4)

式中:ωp為等離子體頻率;νc為電子平均碰撞頻率. 按照移位算子方法,式(4)可改寫為

(5)

通過改變式(4)中參數(shù)ε∞、ωp和νc的取值,來模擬出某一頻率下不同的介質(zhì)材料.

2 應用算例

2.1 低頻驗證

考慮整個計算空間填充海水媒質(zhì)的計算模型,海水媒質(zhì)采用Drude色散介質(zhì)模型來模擬. 在工作頻率f=30 Hz、150 Hz時,式(4)中相關參數(shù)設置如下:

νc=1×10-20,ε∞=1×108.

(6)

將式(6)代入式(4)可以得到海水在工作頻率f=30 Hz、150 Hz下的介電常數(shù)都為80. 計算過程中,空間離散間隔為δ=dx=dy=dz=0.05 m. 時間間隔ΔtN=N·dx/(2c),其中c為真空中的光速. 在這里,變量N取9 000,即Δt9000=9000·dx/(2c). 為了模擬開域空間,設置整個計算空間被厚度為8個網(wǎng)格的完全匹配層(convolutional perfect matched layer, CPML)包圍[5].

輻射源設置為x方向極化的磁偶極子源,加源的表達式為

(7)

式中,磁偶極矩pinc(t)為

(8)

式(8)中:脈寬τ=1.5 ms;頻率ωi=2πi,i=1,2,3;系數(shù)β為初始振幅,在這里β=1.

2.1.1 工作頻率f=30 Hz

整個計算空間為長方體,其物理尺寸取為1000 m×2 m×3 m. 偶極子輻射場中Hx和Ez分量在yoz面上振幅相位的快照如圖1所示.

(a) Hx振幅 (b) Hx相位(a) Ampulitude of Hx (b) Phase of Hx

(c) Ez振幅 (d) Ez相位(c) Ampulitude of Ex (d) Phase of Ex 圖1 f=30 Hz時yoz面上振幅相位的快照Fig.1 The amplitude and phase distributions in the yoz plane at the working frequency f=30 Hz

觀察圖1中Hx及Ez的相位快照可以看出,在海水中工作頻率f=30 Hz的波長大約為λ=385 m. 根據(jù)電磁波理論,這個頻率的電磁波在海水中波長的理論解大約為333 m. 這兩個結果有差別是由于根據(jù)圖1相位變化的方法測得的波長,觀測位置離偶極子源很近,沒有位于偶極子輻射場的遠區(qū)區(qū)域,這就造成與利用平面波特性方法計算得到的理論解結果不一致. 如果圖1所示的相位變化法的觀測位置是離偶極子源更遠的位置,這兩個值之間的差別會逐漸縮小. 這一原因也通過計算工作頻率f=150 Hz在海水中波長的算例得到了驗證.

為了驗證偶極子產(chǎn)生場振幅的變化情況,根據(jù)電磁波理論給出了電偶極子在有耗媒質(zhì)中產(chǎn)生場的解析解為

(9)

(a) 歸一化Hx幅值 (b) 相對誤差(a) Normalized ampulitude of Hx (b) Relative error 圖2 解析解和計算解計算結果的對比Fig.2 Comparison between analytical and computational solutions

2.1.2 工作頻率f=150 Hz

本算例中,計算空間的物理尺寸為450 m×4 m×4 m. 偶極子輻射場中Hx和Ez分量在yoz面上的振幅相位快照如圖3所示.

(a) Hx振幅 (b) Hx相位(a) Ampulitude of Hx (b) Phase of Hx

(c) Ez振幅 (d) Ez相位(c) Ampulitude of Ex (d) Phase of Ex 圖3 f=150 Hz時yoz面上振幅相位的快照Fig.3 The amplitude and phase distributions in the yoz plane at the working frequency f=150 Hz

同理,根據(jù)電磁波理論,海水中工作頻率f=150 Hz波長的理論解約為λ=149 m. 根據(jù)相位變化法從圖3中可以讀出波長的計算解約為λ=157 m. 通過對比圖1中的結果,可以看出隨著工作頻率的升高偶極子產(chǎn)生的場會更快地到達遠區(qū)場. 圖4給出yoz面上輻射場沿y軸的變化情況.

將計算所得結果進行歸一化,與解析解結果進行對比,結果及相對誤差如圖4所示. 通過觀察圖2(b)和圖4(b)可以看出,計算結果和理論結果吻合.

(a) 歸一化Hx幅值 (b) 相對誤差(a) Normalized ampulitude of Hx (b) Relative error 圖4 解析解與歸一化后計算結果的對比Fig.4 Comparison between analytical and computational solutions

2.1.3 加速性能的測試

為了測試本文方法的加速計算性能,設置了本算例.計算環(huán)境為惠普Z620工作站(CPU E5-2620 v2 2.10 GHz,內(nèi)存32 GB). 計算空間的網(wǎng)格尺寸為30δ×30δ×30δ,時間間隔為ΔtN=Nδ/(2c),變量N取正整數(shù). 在本算例中N取三個數(shù)值,即N=1,5,7.加速比Sp和加速效率Ep的定義與文獻[8]中的定義相同,計算結果如表1所示.

表1 加速性能的測試結果Tab.1 Test results of the accelerated performance

在測試程序加速性能算例中,要使三次計算的總時間保持一致,即三次計算需滿足M×Δt一樣. 通過觀察表1可知,隨著N的增加,加速性能接近線性變化,加速效率達到約95%以上.

2.2 金屬殼體目標的電磁散射

2.2.1 起伏海面的建模及影響

運用Monte Carlo方法建立起伏海面[9],基本原理簡述如下:假設x、y方向的長度分別為Lx、Ly,兩個方向的采樣點數(shù)分別為M、N,采樣間距分別為Δx和Δy,則起伏面上每一點(xm,yn)(xm=mΔx,yn=nΔy)的高度表示為

F(km,kn)exp[j(kmxm+knyn)].

(10)

式中,

(11)

式(11)中:S(km,kn)為起伏海面的功率譜;km=2πm/Lx;kn=2πn/Ly;二維高斯粗糙表面的功率譜為

(12)

式(12)中:δ代表粗糙面的均方根高度;lx和ly代表沿x軸和y軸的相關長度.

考慮海平面為三級海況時,可令粗糙面的均方根高度δ=0.8 m,相關長度為lx=ly=0.2 m,得到海面如圖5(a)所示;保持均方根高度不變,改變其相關長度為lx=ly=0.15 m,得到起伏海面如圖5(b)所示.

將圖5所示的起伏海面模型分別疊加到圖6所示的模型中,對頻率f=150 Hz的低頻場進行計算. 為了查看起伏海面對低頻電磁波的影響,將上述兩種起伏海面的計算結果做差值處理,得到海面對低頻電磁波的影響,如圖7所示. 觀察圖7海面可以看出起伏海面影響低頻電磁場的空間分布,特別是對電場分布的影響.

(a) δ=0.8 m,lx=ly=0.2 m

(b) δ=0.8 m,lx=ly=0.15 m圖5 起伏海面示意圖Fig.5 Schematic diagram of rough sea surface

圖6 計算起伏海面下低頻電磁波傳播的三維模型Fig.6 3D model for calculating the low frequency electromagnetic wave under rough sea surface

(a) Hx振幅差 (b) Ez振幅差(a) Amplitude difference of Hx (b) Amplitude difference of Ex圖7 起伏海面對低頻電磁波的影響Fig.7 Impact of rough sea surface on the low-frequency electromagnetic waves

2.2.2 水下金屬殼體目標

圖8給出了計算平滑海面下低頻電磁波傳播的三維模型.在z方向依次建立空氣和海水兩種分層媒質(zhì),其中空氣媒質(zhì)厚度為2.5 m,海水媒質(zhì)厚度為6 m. 整個計算空間的物理尺寸為21 m×6 m×8.5 m. 建立如圖8所示的金屬橢球狀殼體目標,其半長軸為7.5 m,半短軸為1 m. 殼體目標中心點距離海面高度為3 m,殼體厚度為0.1 m,殼體為金屬材質(zhì),殼體內(nèi)部填充空氣.

圖8 計算平滑海面下低頻電磁波傳播的三維模型Fig.8 3D model for calculating the low frequency electromagnetic wave under the smooth sea surface

整個計算空間的媒質(zhì)材料都采用有耗Drude模型進行模擬. 激勵源為x方向極化的磁偶極子源,如式(8)所示. 在工作頻率f=150 Hz時,空氣、海水和金屬的相關設計參數(shù)為:

空氣

vc=1×10-18,ε∞=1×1010;

(13)

海水

vc=1×10-18,ε∞=1×1010;

(14)

金屬

vc=1×10-18,ε∞=1×1010.

(15)

經(jīng)過計算得到的電磁場各分量在yoz面上的振幅快照如圖9所示.

(a) Hx振幅 (b) Hy振幅(a) Ampulitude of Hx (b) Ampulitude of Hy

(c) Hz振幅 (d) Ex振幅(c) Ampulitude of Hz (d) Ampulitude of Ex

(e) Ey振幅 (f) Ez振幅(e) Ampulitude of Ey (f) Ampulitude of Ez圖9 低頻電磁場各分量在yoz面上的振幅快照Fig.9 The amplitude distributions in the yoz plane

觀察圖9可以看出,殼體目標影響低頻電磁波的傳播. 由于殼體目標是金屬導體,所以對低頻電磁波的傳播起到了引導的作用. 由于工作頻率為f=150 Hz,通過計算發(fā)現(xiàn)低頻電磁波透過殼體進入殼體內(nèi)部傳播,但是量級很小. 海面影響低頻電磁波各分量的空間分布,特別是電場分量,這是由于海水的電導率比空氣中的電導率要大. 觀察圖9(f)還可以得到,當?shù)皖l電磁波穿過海面到達空氣中傳播時,其傳播速度要比在海水中的傳播速度快,特別是對電場的空間分布影響比較大,從圖中可以很明顯地看出海面的空間位置. 由圖9(d)可知,由于水下目標的存在導致在海水中存在電場Ex分量.

為了考慮水下目標散射場的空間分布,將圖9的計算結果重新作了處理,即將存在橢球體目標計算的結果減去無橢球體目標的計算結果,具體結果如圖10所示.

(a) Hx振幅 (b) Hy振幅(a) Ampulitude of Hx (b) Ampulitude of Hy

(c) Hz振幅 (d) Ex振幅(c) Ampulitude of Hz (d) Ampulitude of Ex

(e) Ey振幅 (f) Ez振幅(e) Ampulitude of Ey (f) Ampulitude of Ez圖10 水下目標在yoz面上的低頻散射場Fig.10 The amplitude distributions in the yoz plane

觀察圖10可以明顯地看出,在水下殼體目標周圍存在散射場,且殼體目標左側的散射場要比右側散射場強,這是由于低頻源離殼體目標的左側較近. 由于低頻電磁源的照射,在金屬殼體表面產(chǎn)生了感應電流,從而在殼體目標的右端點會產(chǎn)生較強的散射場,例如圖10(b)的Hy分量和圖10(e)的Ey分量. 由于海水電導率的存在,導致電場的衰減要比磁場的衰減快,這可以從圖10(a)～(c)與圖10(d)～(f)的對比中得出. 顯然,海面也影響散射場各分量的空間分布,特別是對電場各分量的影響比較明顯,而且散射場在海面上方的傳播速度要比在海水中的傳播速度快.

3 結 論

利用本文方法計算了水下殼體目標的散射場分布,并對該方法的加速性能進行了測試.計算結果表明海面的存在可以影響殼體目標散射場的空間分布規(guī)律,特別是對電場空間分布影響比較大,并且不同的磁分量具有不同的空間分布規(guī)律.本文的計算結果可以為水下殼體目標的探測、識別和跟蹤提供參考.