吉 熙，廖 成，張東民，鄧小川，馮 菊

（西南交通大學 電磁場與微波技術研究所，四川 成都 610031）

0 引言

礦井、鐵路及公路等通信需求，使得隧道中的電波傳播特性研究一直深受國內外學者的關注[1]。傳統(tǒng)上，隧道被建模為理想光滑壁，很少考慮到表面隨機粗糙對電波傳播特性的影響。而最近隨著5G 時代的來臨和無人駕駛技術的應用需要，精確的隧道電磁態(tài)勢分布的重要價值愈加凸顯。因此，粗糙表面的隧道電波傳播研究具有重要意義。

在大尺度的隧道電波傳播預測分析中，目前已經(jīng)有許多數(shù)值方法提出，如射線追蹤法[1]、模式理論[2-3]和拋物方程方法。其中，射線追蹤法利用幾何光學原理追蹤由發(fā)射源發(fā)出的每條射線的傳播路徑，計算出其與周圍環(huán)境所發(fā)生的反射作用，根據(jù)給定的條件判斷其是否能夠到達接收天線；計算每一條收發(fā)天線之間的路徑，并存儲其總長度和每一次反射系數(shù)；利用疊加原理計算到達接收天線的所有路徑的信號總能量[4]。在簡單的隧道場景，射線追蹤法非常有效。但是，對于復雜的隧道場景，該方法計算過程十分繁瑣。模式理論是由麥克斯韋方程組在特定的邊界條件下推導而來的，假設電波是以不同模式的組合進行傳播[3]。模式理論是解析解，計算精度高，但僅適用于少數(shù)形狀規(guī)則、簡單的隧道環(huán)境。拋物方程（Parabolic Equation，PE）由波動方程推導而來[5]，采用步進迭代的方式進行求解。引入分布傅里葉變換技術（Split-step Fourier Transform，SSFT）[6]或交替方向隱式差分技術（Alternating Direction Implicit，ADI）[7-8]后，在保證計算精度的同時，也保證了很高的運算效率，非常適用于隧道中的電波傳播預測。目前，人們對拋物方程應用于隧道的研究還主要集中在光滑壁上。雖然R.Martelly 與R.Janaswamy 使用PE 的ADI 解法研究了粗糙隧道的電波傳播，但并未研究半圓拱形截面隧道的電波傳播[7]。Xingqi Zhang 研究了粗糙半圓拱形的隧道電波傳播，采用ADI 方法。為了保證差分格式的穩(wěn)定性和計算精度，通常要求較小的計算步長，影響了求解效率[9]。

針對上述問題，本文提出了一種基于SSFTPE 的粗糙隧道電波傳播計算方法。通過在傳統(tǒng)的SSFT-PE 基礎上引入等效粗糙因子，得到能應用于表面粗糙隧道的修正模型，既能兼顧解算精度又能保證求解效率。本文首先回顧了隧道中模式分析方法，其次詳細介紹了基于SSFT-PE 的表面粗糙隧道電波傳播計算方法，最后結合數(shù)值算例，將本文方法與模式分析方法的計算結果進行對比驗證，并利用本文方法分別研究了矩形隧道結構與半圓拱形隧道結構下，不同粗糙程度和不同電磁波頻段對隧道電波傳播的影響。

1 模式分析方法

對于形狀規(guī)則、結構簡單的隧道環(huán)境，模式理論能夠得到隧道內場分布的解析解。模式理論認為，隧道中的場分布是一系列不同模式EHp,q的加權組合，其中下標p和q表示模式的階數(shù)。矩形隧道如圖1所示，隧道內任意位置的電場可由式（1）求出[1]：

其中，Ap,q是本征模函數(shù)，Et是由發(fā)射功率決定的發(fā)射電場，隨機粗糙衰減常數(shù)和相位常數(shù)分別為：

圖1 矩形隧道

（1）平行極化時，Ap,q、αp,q由式（4）與式（5）求出：

（2）垂直極化時，Ap,q、αp,q由式（6）與式（7）求出：

其中：

復介電常數(shù)可由式（9）表示：

其中，εra,b是側壁和上下壁的相對介電常數(shù)，σa,b是側壁和上下壁的電導率。

2 基于SSFT-PE 的表面粗糙隧道電波傳播計算方法

2.1 SSFT-PE 基本理論

在直角坐標系(x,y,z)中，假設電磁場的時諧因子為e-iwt，電場或磁場分量ψ滿足如式（10）形式的三維標量Helmholtz 方程[10]：

其中ψ是任意方向的電磁場分量，為拉普拉斯算子，為介質折射指數(shù)。令沿z軸正方向傳播的波函數(shù)為u(x,z)=e-ik0zψ(x,y,z)，代入式（10）并對其進行因式分解可得：

其中Q為偽微分算子：

這里僅考慮電磁波的前向傳播，得到前向拋物方程為：

將偽微分算子Q進行Feit-Fleck 分解得到：

將式（14）代入式（11），得到如式（15）形式的拋物方程：

對式（11）采用分步傅里葉變換（Split-Step Fourier Transfrom，SSFT）算法進行求解，可以得到如式（16）形式的遞推公式[11]：

其中，?2與分別表示二維傅里葉正變換和逆變換。二維傅里葉變換可以通過二維FFT 技術實現(xiàn)快速計算，因此SSFT 算法是一個高效的解算算法。

2.2 適用于表面粗糙隧道的修正SSFT-PE 方法

2.2.1 介質隧道計算方法

對每種本征模而言，根據(jù)歐拉定理可以將本征模式展開為兩種平面波的疊加。邊界條件為第一類時對應奇對稱，鏡像場為負場，而第二類邊界條件對應偶對稱，鏡像場為正場。因此，將矩形隧道進行2 次鏡像后，原計算域擴展為4 倍。通過計算域的擴展和鏡像場的填充，可以等效考慮隧道的兩類邊界條件。根據(jù)各模式本征函數(shù)可以求出其對應的入射角，由入射角可以算出隧道內對應方向反射一次時對應的距離[12]：

傳播距離為z時，平面波在介質壁上的損耗L可以通過式（18）和式（19）得到：

其中L0代表拋物方程的計算過程中不同模式的本征損耗，ρ⊥,||是反射系數(shù)。kx與ky代表沿x與y方向的傳播常數(shù)，k0是自由空間中電磁波的波數(shù)。因此，介質矩形隧道環(huán)境中的遞推公式為：

在實際隧道工程中，通常半圓拱形截面形狀的隧道十分常見[13]，見圖2。

圖2 半圓拱形隧道截面形狀

文獻[12]給出了介質半圓拱形隧道環(huán)境中的遞推公式：

令：

則有：

2.2.2 表面粗糙隧道計算方法

對于均勻平面波入射到光滑壁的情況，可以由菲涅爾定律計算出反射系數(shù)[14]，如圖3 所示。

圖3 光滑壁對平面波的反射

此時，菲涅爾反射系數(shù)為：

其中，⊥,||分別表示垂直極化和平行極化，α⊥,||是入射角，且有：

一般來說，粗糙表面對電磁波的反射通過統(tǒng)計學方法進行量化。目前，常認為粗糙表面符合均值u為0、標準差為σh的高斯分布。當σh＜＜λ時，由于粗糙表面造成的反射系數(shù)為[1]：

將式（32）與式（28）相乘，可得到修正后粗糙表面的反射系數(shù)：

將式（33）代入式（18）和式（19），即得到計算矩形粗糙隧道的修正SSFT-PE 算法。同理，將式（33）代入式（21）、式（26）和式（27），即得到計算半圓拱形粗糙隧道的修正SSFT-PE 算法。

3 數(shù)值算例分析

由于現(xiàn)代工程施工的精度越來越高，根據(jù)生活實踐經(jīng)驗以及文獻[1]等參考可以得出，現(xiàn)代隧道的σh不會太大。為了不失一般性，本文選取σh滿足如下條件：

在驗證算例中，選取σh=0 cm 與σh=10 cm；電波傳播特性分析中，研究不同σh對矩形隧道電波傳播的影響時，為了同時滿足σh＜＜λ的條件，故選取σh=0 cm、σh=3 cm、σh=5 cm 與σh=7 cm；研究不同頻率下σh對傳播因子的影響時，選取σh=5 cm。

3.1 驗證算例

為了驗證本文算法處理粗糙隧道的正確性，分別采用修正后PE 法與模式理論法計算矩形粗糙隧道地形下傳播因子隨距離的變化。

矩形隧道參數(shù)見表1，仿真條件設置為高斯反射源頻率0.9 GHz，波瓣寬度為4.34°，垂直極化波與水平極化波，天線安置在隧道中心位置，接收天線高度為3 m，發(fā)射天線高度為3 m。圖4（a）是PE 方法與模式理論垂直極化時傳播因子隨傳播距離變化的對比。當σh=10 cm 時，PE 方法與模式理論吻合得很好。圖4（b）是PE 方法與模式理論水平極化時傳播因子隨傳播距離變化的對比。當σh=10 cm 時，PE 方法與模式理論吻合得很好。圖4證明了修正后的PE 算法計算粗糙隧道的正確性。

表1 矩形隧道參數(shù)

圖4 拋物方程與模式理論的對比

3.2 粗糙隧道的電磁波傳播特性研究

3.2.1 矩形隧道電波傳播特性分析

為了研究不同σh對矩形隧道電波傳播的影響，仿真不同σh時傳播因子隨傳播距離的變化，如圖5所示。為了定量研究不同σh對傳播因子的影響，，其 中F1(i) 表 示σh=0 cm 時的傳播因子，F(xiàn)2(i)表示在給定頻率下σh=3 cm、5 cm、7 cm時沿x軸線方向各網(wǎng)格點的傳播因子。表2 與表3分別給出了2.515 GHz 與2.675 GHz下不同σh的RMSE。隧道參數(shù)如表1 所示，仿真條件設置為高斯反射源頻率分別為5G 頻率2.515 GHz、2.675 GHz，波瓣寬度為4.34°，垂直極化波，天線安置在隧道中心位置，接收天線高度為3 m，發(fā)射天線高度為3 m。從圖5 可以看出，在近場處出現(xiàn)了信號的快衰落，這是由于電波受到隧道壁的作用會產(chǎn)生多次反射，這些反射的分量會產(chǎn)生許多多徑分量形成頻率選擇性信道，造成信道的快衰落。另外，隨著傳播距離的增加，越大的σh傳播因子，電波傳播的損耗越大，這是由于粗糙壁對無線電波的額外散射損耗造成的。

圖5 不同σh 時傳播因子隨距離的變化

表2 f=2.515 GHz 時不同σh 的RMSE

表3 f=2.675 GHz 時不同σh 的RMSE

為了研究在相同σh不同頻率下對傳播因子的影響，定義，其中F1(i)表示σh=0 cm 時的傳播因子，F(xiàn)2(i)表示在給定頻率下σh=5 cm 時沿x軸線方向各網(wǎng)格點的傳播因子。圖6 顯示了垂直極化與水平極化條件下，在σh=5 cm時不同頻率下MSE的變化。其中，2.515～2.675 GHz、3.4～3.5 GHz 與3.5～3.6 GHz 分別屬于中國移動、中國電信與中國聯(lián)通的5G 通信頻段。表4 給出了頻 率 為2.515 GHz、2.675 GHz、3.4 GHz、3.5 GHz與3.6 GHz 的RMSE，表5 給出了0.5～2.4 GHz 的平均RMSE以及這幾個5G 頻點的平均RMSE。從圖6 可以看出，當f從0.5 GHz 開始變化到3.6 GHz，RMSE整體在增加，表明發(fā)射頻率越高，微小的σh將對電波的傳播因子造成顯著影響；從表4 可以看出，0.5～2.4 GHz 時，垂直極化時平均MSE達到0.831 2 dB，而在水平極化時達到了0.521 3 dB；在5G 通信頻段，垂直極化時平均MSE達到了約2.5 dB，水平極化時平均MSE達到了約2.1 dB，明顯高于0.5～2.4 GHz 的平均MSE，說明在5G 通信頻段σh對電波的傳播影響顯著。

表4 不同頻點下的平均RMSE

表5 不同頻段下的平均RMSE

圖6 不同頻率下RMSE 的變化

3.2.2 半圓拱形隧道電波傳播特性分析

為了研究不同σh對半圓拱形隧道電波傳播的影響，仿真不同σh時傳播因子隨傳播距離的變化。如圖7 所示，為了定量研究不同σh對傳播因子的影響，表6 與表7 分別給出了2.515 GHz 與2.675 GHz 下不同σh的RMSE。半圓拱形隧道的參數(shù)為W=8.43 m、H=5.99 m，仿真條件設置為高斯反射源頻率分別為2.515 GHz、2.675 GHz，波瓣寬度為4.34°，εr=5.5，σ=0.1 S/m。天線安置在隧道中心位置，發(fā)射天線高度為3 m，垂直極化波，接收天線高度為3 m。從圖7 可以看出，隨著傳播距離的增加，越大的σh傳播因子，電波傳播的損耗增加，與不同σh對矩形隧道電波傳播的影響的結論相同。

圖7 不同σh 時傳播因子隨距離的變化

表6 f=2.515 GHz 時不同σh 的RMSE

表7 f=2.675 GHz 時不同σh 的RMSE

圖8 是垂直極化與水平極化條件下半圓拱形隧道中不同頻率下σh對傳播因子的影響。為了與矩形隧道中電波傳播特性進行比較，設半圓拱形隧道的參數(shù)也為W=8.43 m、H=5.9 m，仿真條件設置為高斯反射源頻率，波瓣寬度為4.34°，εr=5.5，σ=0.1 S/m。發(fā)射天線高度為3 m，接收天線高度為3 m。表8給 出 了2.515 GHz、2.675 GHz、3.4 GHz、3.5 GHz與3.6 GHz 的RMSE，表9 給出了0.5～2.4 GHz 的平均RMSE和這幾個5G 頻點的平均RMSE。從圖8 可以看出，當f從0.5 GHz 開始變化到3.6 GHz，RMSE整體在增加，表明發(fā)射頻率越高，微小的σh將對電波的傳播因子造成顯著影響；從表8 可以看出，0.5～2.4 GHz 時，垂直極化時平均RMSE達到了0.858 9 dB，而在水平極化時達到了0.744 9 dB；在5G通信頻段，垂直極化時平均RMSE達到了約2.6 dB，水平極化時平均RMSE達到了約2.1 dB，明顯高于其他頻段，說明在5G 通信頻段時σh對電波的傳播影響顯著。另外，比較圖6 與圖8 可以看出，矩形隧道的RMSE在2.0～2.4 GHz 時波動明顯，而半圓隧道的RMSE在2.8～3.6 GHz 時波動明顯；垂直極化的額外粗糙損耗大于水平極化時。

表8 不同頻點下的平均RMSE

表9 不同頻段下的平均RMSE

圖8 不同頻率下RMSE 的變化

4 結語

本文提出了一種基于SSFT-PE 的粗糙隧道電波傳播模型。根據(jù)粗糙壁服從高斯分布統(tǒng)計下的反射系數(shù)模型，將該模型與傳統(tǒng)的SSFT-PE 相結合，得到修正后的粗糙隧道傳播模型，為預測真實隧道環(huán)境中的電波傳播特性提供了一種快速準確的方法。在8.43 m×5.99 m 的典型矩形和拱形隧道中進行仿真分析，結果表明：當隨機粗糙標準差σh等于5 cm 時，在0.5～2 GHz 頻段，隨機粗糙帶來額外的約0.5 dB 損耗；在2～3.6 GHz 頻段，隨機粗糙帶來額外的2.1 dB 損耗，在5G 頻段的損耗大于其他低頻段；在相同仿真頻率下，額外的損耗隨著粗糙程度的增加而增加。