寬角拋物方程模型處理復(fù)雜地形的一種方法*

俞紅兵 ，胡繪斌 ，儲(chǔ)飛黃

（1.國(guó)防科技大學(xué)電子對(duì)抗學(xué)院，安徽 合肥 230037；2.94860部隊(duì)，江蘇 南京 210018）

0 引 言

在現(xiàn)代電子戰(zhàn)中，根據(jù)戰(zhàn)場(chǎng)環(huán)境的實(shí)際數(shù)據(jù)準(zhǔn)確分析和評(píng)估復(fù)雜地形環(huán)境對(duì)電波傳播和電子戰(zhàn)系統(tǒng)作戰(zhàn)性能的影響，對(duì)于輔助指揮員進(jìn)行正確的戰(zhàn)術(shù)部署和選擇恰當(dāng)?shù)膽?zhàn)略戰(zhàn)術(shù)，以獲取戰(zhàn)場(chǎng)電磁優(yōu)勢(shì)、掌握戰(zhàn)場(chǎng)制電磁權(quán)具有極其重要的意義。迄今為止，人們已經(jīng)建立了多種計(jì)算復(fù)雜地形上電波傳播特性的數(shù)學(xué)模型，其中，拋物方程模型（Parabolic Equation Model，PEM）是一種基于Maxwell方程建立起來(lái)的確定性電波傳播模型，與幾何光學(xué)(GO)或物理光學(xué)（PO）模型相比，其計(jì)算過(guò)程簡(jiǎn)單且精度高，與矩量法（MOM）、有限元法（FEM）等積分方程方法和時(shí)域有限差分法（FDTD）等微分方程方法相比，PEM對(duì)內(nèi)存要求低、計(jì)算速度快，可以進(jìn)行遠(yuǎn)距離電波傳播特性預(yù)測(cè)。由于上述優(yōu)點(diǎn)，近年來(lái)，PEM得到了人們廣泛的關(guān)注和研究[1］。

和求解任何形式的Maxwell方程組都需要初始條件和邊界條件一樣，用PEM計(jì)算電波傳播問(wèn)題時(shí)，主要需解決三個(gè)問(wèn)題：（1）如何確定PEM的初始值？（2）采用何種數(shù)值方法求解拋物方程（PE）？（3）在相應(yīng)的數(shù)值求解方法中，如何處理PE的邊界條件？對(duì)前兩個(gè)問(wèn)題，前人已進(jìn)行過(guò)比較深入的研究[2-3］，在此不再贅述，本文主要探討第（3）個(gè)問(wèn)題——PEM處理復(fù)雜地形邊界條件的方法。為此，首先介紹計(jì)算電波傳播特性的PEM。

1 PEM概述

二維波動(dòng)方程通過(guò)一定的近似處理可以轉(zhuǎn)化為以下的Feit-Fleck型二維寬角PE[2］：

式中，u(x,z)為任意距離-高度點(diǎn)上的場(chǎng)分布，n(x,z)=n0(x,z)+z/Re為隨地形高度變化而修正的大氣折射指數(shù)。n0為發(fā)射天線所在高度處的大氣折射指數(shù)，x、z分別表示傳播距離和高度，Re≈6370 km為地球半徑，k0=2π/λ為自由空間傳播常數(shù)。

式（1）就是寬角PEM。根據(jù)收發(fā)點(diǎn)之間的地形高程數(shù)據(jù)（可通過(guò)數(shù)字地圖來(lái)獲得）和地表電磁特性（主要為介電常數(shù)εr和電導(dǎo)率σ，可按ITU_R 527-4的規(guī)定來(lái)確定或人為確定[4］），結(jié)合收發(fā)點(diǎn)之間路徑上的大氣折射指數(shù)n(x,z)分布（對(duì)于低于1 km的低層大氣可通過(guò)已有成熟模型計(jì)算獲得[5］），采用離散分步混合付立葉變換方法（DMFT）來(lái)數(shù)值求解1)式[3,6］，就可計(jì)算出收發(fā)點(diǎn)之間的電波傳播路徑衰減。本文主要闡述PEM處理由地形高程數(shù)據(jù)描述的不規(guī)則地形邊界的方法和步驟。

2 PEM處理復(fù)雜地形邊界的方法

復(fù)雜地形環(huán)境主要指電波傳播區(qū)域的地勢(shì)是高低起伏的，地表特征則表現(xiàn)為海水、淡水、干燥地、濕地等類型[4］。

為了能夠采用DMFT來(lái)求解PE，需將不規(guī)則地形上的PE轉(zhuǎn)換成平地面上PE的形式，為此，需對(duì)圖1所示的不規(guī)則地表進(jìn)行坐標(biāo)變換。

圖1 不規(guī)則地表的平地面轉(zhuǎn)換

以函數(shù)z=T(x)=t(x)-x2/(2Re)來(lái)描述不規(guī)則地表，其中，t(x)為描述地形起伏邊界的函數(shù)，x2/(2Re)項(xiàng)則是考慮了地球曲率的影響。以不規(guī)則地形上某一點(diǎn)作為坐標(biāo)原點(diǎn)，建立一個(gè)新坐標(biāo)系(x′,z′)，且滿足：

設(shè)u(x,z)=e-iθψ(x,z)，θ為待定變量，ψ(x,z)為任一標(biāo)量場(chǎng)分量，根據(jù)圖1，將式（2）代入式（1），可以得到以下形式的寬角PE[3］：

為了使式（3）與式（1）形式相同，設(shè)

將式（4）、式（5）代入到式（3）中可得：

這就是經(jīng)平地面地形變換后所導(dǎo)出的寬角PE模型，通常將包含n(x′,z′)的項(xiàng)以修正折射指數(shù)表示為：

式（6）和式（1）相比僅僅是折射指數(shù)項(xiàng)不同。在式（7）中，若?T2/?x′2=0（即平地面），則m(x′,z′)=n(x′,z′)-1，這時(shí)式（6）就與式（1）完全一致?？梢?jiàn)，通過(guò)平地面轉(zhuǎn)換，不規(guī)則地形對(duì)寬角PE的影響完全轉(zhuǎn)換為平地面上對(duì)大氣折射指數(shù)分布的影響。

當(dāng)從數(shù)字地圖中抽取傳播地形剖面數(shù)據(jù)時(shí)，這種離散數(shù)據(jù)使得地形剖面是分段線性的，只存在?T/?x′項(xiàng)，?T2/?x′2=0，式（6）符合平地面條件，可采用DMFT求解。以下本文將詳細(xì)闡述對(duì)分段線性地形的處理方法。

為簡(jiǎn)單起見(jiàn)，將 ?2T/?x′2用二階中心差分的形式來(lái)代替，即

假設(shè)地形剖面上A、B兩段是不連續(xù)的，其斜率分別為T(mén)A′和TB′，且TA′≠TB′，A、B兩段地形上方的場(chǎng)分別以ψA和ψB表示，如圖2所示。

圖2 地形剖面中不連續(xù)的兩段斜面

由于A、B的斜率不同，在電磁波由A傳播到B的過(guò)程中，根據(jù)式（4），在A、B段的交點(diǎn)xAB處θ是不連續(xù)的，而u(x,z)=e-iθψ(x,z)，因此ψA和ψB在相位上也是不連續(xù)的。為了反映地形變化對(duì)電磁場(chǎng)相位的影響，本文設(shè)定：

式（9）的思想就是當(dāng)電磁波由A傳播到B時(shí)，根據(jù)兩段地形斜率強(qiáng)行改變B段上方場(chǎng)的相位，這與McArthur和Bebbington提出的扭轉(zhuǎn)波陣面的思想[7］十分相近。

在計(jì)算地表以上各高度網(wǎng)格上的場(chǎng)分布時(shí)，對(duì)于平地面邊界的情況，可采用DMFT方法計(jì)算出任意步進(jìn)上、從0至zmax高度范圍內(nèi)的場(chǎng)分布。對(duì)于不規(guī)則地形，各高度網(wǎng)格上場(chǎng)的計(jì)算也都要從地表算起，但此時(shí)地表高度不一定為零，因此對(duì)于分段線性地形剖面，除了要計(jì)算不連續(xù)地形交界點(diǎn)處的場(chǎng)以外，還需要將不同的地表高度轉(zhuǎn)換為相對(duì)零高度來(lái)計(jì)算。為此，本文采用所謂的“邊界平移法”（Boundary Shift，BS）來(lái)進(jìn)行這樣的轉(zhuǎn)換處理，具體如下。

如圖3所示，設(shè)電波傳播到某一距離點(diǎn)D1處，其地表上方波陣面上的場(chǎng)為u(x,iΔz)，i=1,2,…N，N為DMFT算法的尺度。經(jīng)過(guò)一個(gè)步進(jìn)后，波陣面?zhèn)鞑サ紻2，設(shè)D2與D1處地表的高度差為T(mén)x，則Tx所包含的高度網(wǎng)格數(shù)為

顯然，Tx-NtΔz＜Δz。

圖3 邊界平移（BS）法示意圖

在用BS方法計(jì)算D2處的場(chǎng)分布時(shí)，如果D2處是上升地形，即Tx＞0，則首先由D1處的場(chǎng)分布計(jì)算出D2處的場(chǎng)分布，再將此場(chǎng)分布的網(wǎng)格序號(hào)往下平移Nt這樣，D2處最大高度上所對(duì)應(yīng)的網(wǎng)格序號(hào)就由N-1下降為N-Nt-1，再將超出最大高度的第(N-Nt)～(N-1)個(gè)網(wǎng)格點(diǎn)上的場(chǎng)設(shè)為零，就得到D2處從地表起至最大高度處各網(wǎng)格上的場(chǎng)分布，即

如圖3（a）所示。顯然，此時(shí)D2處的u所包含的元素個(gè)數(shù)還是N-1。

如果D2處是下降地形，即Tx＜0，同樣地，將D2處場(chǎng)分布的序號(hào)往上平移Nt，同時(shí)，不高于NtΔz的各高度網(wǎng)格點(diǎn)上的場(chǎng)也設(shè)為零，即：

如圖3（b）所示。顯然，此時(shí)D2處的u所包含的元素個(gè)數(shù)也還是N-1。

如果D2和D1高度相同，即Tx=0，則維持原來(lái)波陣面上的場(chǎng)分布，不進(jìn)行BS轉(zhuǎn)換。

通過(guò)以上的BS方法處理，在計(jì)算仰角范圍內(nèi)，式（6）就可以計(jì)算不規(guī)則地形上任意距離和高度處的場(chǎng)分布。例如，要計(jì)算距離為D、高度為H的p點(diǎn)處的場(chǎng)值，如圖4所示，圖中以▲表示p點(diǎn)所處的位置。

圖4 不同高度、不同距離點(diǎn)處場(chǎng)的計(jì)算方法

當(dāng)p點(diǎn)恰好位于高度網(wǎng)格上，即mod(H/Δz)=0時(shí)，mod表示取余數(shù)計(jì)算，則第H/Δz網(wǎng)格上的場(chǎng)即為所求；當(dāng)p點(diǎn)不位于任何高度網(wǎng)格上時(shí)，mod(H/Δz)≠0，此時(shí)可以通過(guò)對(duì)p點(diǎn)上下兩個(gè)網(wǎng)格點(diǎn)上的場(chǎng)進(jìn)行插值計(jì)算來(lái)獲得p點(diǎn)的場(chǎng)值。不妨設(shè)

則采用線性內(nèi)插方法時(shí)，p點(diǎn)的場(chǎng)為：

當(dāng)距離點(diǎn)D恰好位于某個(gè)步進(jìn)處時(shí)，該步進(jìn)上的u即為D處的場(chǎng)分布，當(dāng)D不位于任何一個(gè)步進(jìn)上時(shí)，同樣可以采用線性內(nèi)插的方法求得D處的場(chǎng)分布，即

如圖（4）所示。

至此，利用DMFT算法[6］，結(jié)合上述對(duì)不規(guī)則地形的處理方法，就可以用式（6）的寬角PEM計(jì)算不規(guī)則地形條件下任意距離和高度點(diǎn)處的電波傳播衰減。

3 算例

為了便于比較，本文采用文獻(xiàn)[8］中的一個(gè)算例，同樣不考慮自由空間電波傳播損耗。

算例：標(biāo)準(zhǔn)大氣條件下（即?n/?h≈-40×10-6時(shí)的大氣環(huán)境），電波在一個(gè)類似正弦包絡(luò)形狀的不規(guī)則地形上傳播，其中，地形函數(shù)為：

顯然，這是一個(gè)?2T/?x2≠0的連續(xù)函數(shù)。式中，h=229 m為包絡(luò)的最大高度，x1=15 km為包絡(luò)半高程所在的距離，w=10 km為包絡(luò)的半寬度。地形剖面如圖5所示。

圖5 類似正弦包絡(luò)的地形剖面

發(fā)射天線高度為Ht=30.5 m，頻率為f=3.0 GHz，高斯方向圖：A(p)=e-p2ω2/4，p=k0sin(θ-θ0)，θ0為發(fā)射天線的仰角，θ為每一步進(jìn)處的射線仰角，為方向圖的3 dB寬度。電波水平極化。最大傳播距離為rmax=50 km。所有地表均設(shè)為PEC。

圖6（a）顯示了本文計(jì)算結(jié)果的偽彩圖，圖中的不同顏色代表自地面起至600 m高度區(qū)域內(nèi)的電波傳播衰減值分布，其中，圖下方的色帶顯示了不同衰減值所對(duì)應(yīng)的顏色。

圖6（a）和文獻(xiàn)[8］中的圖2基本上是一致的，主要的區(qū)別就是：在0～5 km距離范圍內(nèi)，在一些較高的點(diǎn)上，圖6（a）沒(méi)有給出衰減計(jì)算結(jié)果。這是因?yàn)檫@些點(diǎn)處的仰角超過(guò)了寬角PE模型的計(jì)算仰角，沒(méi)有計(jì)算。文獻(xiàn)[8］中的圖2雖然顯示了0～5 km距離內(nèi)較大仰角處的衰減值，但不排除這是采用其它方法計(jì)算所得的結(jié)果。

圖6（b）則顯示了最大距離處電波傳播衰減值隨高度變化的情況。從圖中可見(jiàn)，從地面起至地形最大高度處，傳播區(qū)域逐漸由陰影區(qū)過(guò)渡到亮區(qū)，在200 m高度以上，由于天線的高度增益，衰減值幾乎沿直線減小?？傮w來(lái)看，電波在類正弦地形上的傳播衰減特性和在單刃峰地形上的繞射效應(yīng)類似，只是其衰減值的振蕩更大。

圖6 電波在類正弦包絡(luò)地形上的傳播衰減

通過(guò)以上算例可以看出，本文方法計(jì)算不規(guī)則地形上的電波傳播特性是比較準(zhǔn)確的。由于算例中所有的地表都設(shè)為PEC，因此，在計(jì)算PE初始場(chǎng)時(shí)，地表反射系數(shù)取為-1。實(shí)際上，對(duì)于非PEC的地表結(jié)構(gòu)，本文方法依然適用，只是需要根據(jù)地表的相對(duì)介電常數(shù)εr和電導(dǎo)率σ，計(jì)算出地表反射系數(shù)，再由反射系數(shù)計(jì)算出新的PE的初始場(chǎng)分布即可，地表類型并不影響本文對(duì)地表的處理方法。

4 結(jié) 語(yǔ)

通過(guò)以上算例可以看出，本文給出的寬角PE模型處理復(fù)雜地形的方法是有效的、準(zhǔn)確的，可以為計(jì)算電磁波在各種復(fù)雜地形上的傳播衰減提供一條行之有效的途徑。