一种求解三维抛物方程的迭代方法_吴晓萍.pdf

吴晓萍，梁志禧，黎子豪，等.一种求解三维抛物方程的迭代方法J.电波科学学报，2023，38（1）：159-163.DOI:10.12265/j.cjors.2022009WU X P,LIANG Z X,LI Z H,et al.An iterative method for solving 3D parabolic equationJ.Chinese journal of radio science，2023，38（1）：159-163.（inChinese）.DOI:10.12265/j.cjors.2022009一种求解三维抛物方程的迭代方法吴晓萍*梁志禧黎子豪龙云亮*（中山大学,广州 510006）摘要 为了高效且高精度地求解三维抛物方程(parabolic equation，PE)，提出了一种求解三维 PE 的迭代方法.该方法采用有限差分（finite difference，FD）近似方法，保留交替方向隐式（alternating direction implicit,ADI）方法的计算模式.推导了该方法对应的相对误差关系式，在推导过程中考虑计算式分裂项丢失带来的误差影响.对比现有的求解三维 PE 的 ADI 方法(记为 ADI-PE)和克兰克尼科森(Crank-Nicolson)求解方法(记为 CN-PE)，该方法能够提高计算精准度，且可通过迭代达到一定计算精度.该方法保留了 ADI 的计算高效特性，是一种求解三维 PE 的高效高精准度的方法.关键词抛物方程(PE)；有限差分(FD)方法；迭代；交替方向隐式(ADI)；相对误差中图分类号TN011文献标志码A文章编号1005-0388（2023）01-0159-05DOI 10.12265/j.cjors.2022009An iterative method for solving 3D parabolic equationWU Xiaoping*LIANG ZhixiLI ZihaoLONG Yunliang*（Sun Yat-Sen University,Guangzhou 510006,China）AbstractTo solve 3D parabolic equation efficiently(PE)and accurately，an iterative method for solving 3D PEis proposed.In this method,finite difference(FD)approximation is applied.The proposed method retains the characterof the alternating-direction-implicit(ADI)method.The numerical error caused by the split term of the proposed methodis considered in the process of derivation.The proposed method can improve the calculation accuracy compared to theADI-PE method and Crank-Nicolson parabolic equation(CN-PE)method.A certain accuracy can be achieved byseveral times of iteration.This method is efficient and accurate for solving 3D PE.Keywordsparabolic equation(PE)；finite difference(FD)method；iteration；lternating-direction-implicit(ADI)；relative error 引言抛物方程(parabolic equation，PE)法是预测无线电波传播最常用的方法之一.该方法能模拟近轴方向，即模拟电磁波传播方向上的圆锥形区域内的电磁波传播.在远距离传播问题上，PE 具有计算高效的特点，且能够同时处理大气结构不均匀和地形不规则引起的折射和绕射效应1-2.求解 PE 的方法主要有两种：有限差分(finite difference,FD)近似方法3和分步傅里叶(split-step Fourier,SSF)方法4.和 SSF 方法对比，FD 方法在处理计算空间内地形边界问题上具有较为灵活的特点，近年来也成为国内外研究者的研究热点.目前，利用 FD 方法求解 PE 的研究中，郭琪等人推导了高阶近似的二维 PE 求解模型以及适用于不规则地形的高阶模型，并分析了模型所能计算的最大传播仰角和最大的坡度角5-6.此外，郭琪等人7还对吸收边界问题进行研究，提出了一种基于全局透射边界条件的宽角PE 电波预测模型，该模型提高了PE求解的计算效率.Lytaev 等人8-9采用非本地边界条 收稿日期：2022-01-12资助项目：国家自然科学基金(41376041)；广东省自然科学基金研究团队项目(2015A030312010)通信作者：吴晓萍 E-mail:；龙云亮 E-mail: 第 38 卷第 1 期电波科学学报Vol.38，No.12023 年 2 月CHINESE JOURNAL OF RADIO SCIENCEFebruary，2023 件处理计算域的上下边界，提出了一种更高阶的求解二维 PE 的 Numerov-Pad方法，对比 SSF 方法，该方法能够减少 PE 所需的计算区域范围，具有一定计算效率与计算稳定性.利用 PE 的一些求解方法进行应用分析的研究中，文献 10 采用时域有限差分(finite-difference time-domain,FDTD)和二维单向/双向 PE 方法，分析了低频短距离高传播角电波传播特性.王丹丹等人11将宽角与窄角 PE 方法应用到地表低频电波传播的问题研究中，分析了复杂地表参数设置对 PE 求解精度的影响.利用三维 PE 进行高效求解的研究中，Martelly 等人12提出一种求解三维PE 的交替方向隐式(alternating direction implicit,ADI)方法(记为 ADI-PE)，对比克兰克尼科森(Crank-Nicolson)求解方法(记为 CN-PE)1，ADI 方法具有计算效率高的优势.何姿等人13推导了双向 ADI-PE，并将该方法应用在高速公路的实际场景中.黄璐等人14将 ADI-PE 应用到机场盲降系统中，并提出了更高效的并行算法.综上所述，近年来国内外研究者在用 FD 方法求解 PE 的研究主要集中在二维 PE 问题，对三维 PE 问题的研究尚少.在有关 FDTD 求解方法的研究中，文献 15 提出了一种迭代方法，该方法能减少 ADI-FDTD 的计算误差，并保留 ADI 的高效计算特性.ADI-PE 方法是求解三维 PE 的方法中应用较广泛方法，因此，考虑通过迭代的方式对该方法进行优化具有一定的意义.本文提出一种求解三维 PE 的迭代方法，其保留了 ADI 方法的高效计算模式.对比现有应用广泛的ADI-PE 及 CN-PE 方法，该方法能够提高计算精准度、减小误差，且可通过迭代达到一定计算精度.1 三维 PE 迭代方法 1.1 求解三维 PE 的 CN-PE 与 ADI 方法由麦克斯韦方程组推导得到三维直角坐标系下的 PE 为2x2+2y2+2z2=k20n2=0.（1）k0 x式中：为水平极化或垂直极化下的电场或磁场；为波数；n 为折射率(当为空气时，n=1).当电波传播方向为 方向时，引入辅助函数(x,y,z)=u(x,y,z)ejk0 x.（2）假设|2x2|k0|x|，（3）将式(2)代入式(1)得ux=12jk0(2uy2+2uz2).（4）对式(4)中的微分算子进行 FD 近似，化简得到CN-PE 求解方法：(1ry4jk0yrz4jk0z)uk+1i,j=(1+ry4jk0y+rz4jk0z)uki,j.（5）式中：|ry=xy2rz=xz2;yuki,j=uki1,j2uki,j+uki+1,jzuki,j=uki,j12uki,j+uki,j+1.（6）式(5)的 CN-PE 方法求解涉及五对角线的稀疏矩阵，求解效率较低.为提高求解三维 PE 的效率，研究者提出了求解三维 PE 的 ADI 方法12，形式如下：(1rz4jk0z)uk+12i,j=(1+ry4jk0y)uki,j;(1ry4jk0y)uk+1i,j=(1+rz4jk0z)uk+12i,j.（7）1.2 三维 PE 的迭代方法推导将式(5)改写为以下形式：(1ry4jk0y)(1rz4jk0z)uk+1i,j=(1+ry4jk0y)(1+rz4jk0z)uki,jry4jk0rz4jk0yz(uk+1i,juki,j).（8）式(8)等号右边最后一项为 ADI-PE 方法所丢弃的分裂误差项，其对求解三维 PE 精准度的影响与网格步进大小有关.考虑以下关系式：M=(1ry4jk0y)(1rz4jk0z);（9）P=(1+ry4jk0y+rz4jk0z)uki,j;（10）N=ry4jk0rz4jk0yz.（11）根据式(9)(11)，式(8)可以写为Muk+1=Nuk+1+P.（12）该式具有以下线性求解模式：Muk+1n+1=Nuk+1n+P.（13）将式(9)(11)代入式(13)中，可得 160电波科学学报第 38 卷(1ry4jk0y)(1rz4jk0z)uk+1n+1=(1+ry4jk0y)(1+rz4jk0z)uk+1ry4jk0rz4jk0yz(uk+1n+1uk+1).（14）式(14)分裂成两步得到求解三维 PE 迭代方法的计算式：(1ry4jk0y)utmpn+1=(1+rz4jk0z)uk12(ry4jk0rz4jk0yz)(uk+1nuk);(1rz4jk0z)uk+1n+1=(1+ry4jk0y)utmp12(ry4jk0rz4jk0yz)(uk+1nuk).（15）观察式(15)可知，其计算式子保留了 ADI 高效的求解模式，可有效快速求解三维 PE.1.3 误差分析本节推导所提出的三维 PE 迭代方法的误差关系式.亥姆霍兹波动方程的精确解为(x,y,z)=Aej(kxx+kyy+kzz).（16）u(x,y,z)式中，A 为常量.表示为u(x,y,z)=Aekxx+kyy+(kzk0)z.（17）k0将式(17)代入到不同 PE 求解方法中，归一化波数得到，不同方向上的波矢分量为|kx=k0cos sin ky=k0sin sin kz=k0cos.（18）将式(18)代入到式(15)可得所提出的求解三维PE 迭代方法的相对误差关系式为：ejk0,1(1cos)N2=1+jsin2k0sin sin 2N1+jsin2k0cos sin 2N(122sin2k0sin sin 2Nsin2k0cos sin 2N)(Fn11)(1+jsin2k0cos sin 2N);ejk0,2(1cos)N2=1+jsin2k0cos sin 2N1+jsin2k0sin sin 2N(122sin2k0sin sin 2Nsin2k0cos sin 2N)(Fn11)(1+jsin2k0sin sin 2N).（19）=x/(k0z2)z=yN=/z=2/(k0z)Fn1式 中：，；为上一迭代式右端结果.2 仿真结果=4N=5=45将本文所提求解三维 PE 的迭代方法与 CN-PE方法、ADI-PE 方法性能进行对比.当与网格步进 大小相关的系数、方位角时，三种方法相对误差与仰角 关系如图 1 所示.可以看出，三种方法的相对误差都随着传播仰角的增大而增大，但三维 PE 迭代方法的相对误差比 CN-PE 方法和ADI-PE 方法低，且能带来更高的计算精度.510152000.010.020.030.040.05相对误差/()CN-PE 方法ADI-PE 方法三维 PE 迭代方法图 1 仰角 不同时三种计算方法的相对误差Fig.1 Relative error vs.the zenith angle by 3 methods =4 =36 =45zz，或，时 相 对 误 差 与 参 数N 关系如图 2 所示.可以看出，三种方法的相对误差随 N 的增大而减小，最后趋向相等的效果.这是因为N 与波长 及有关，在单位波长下，越大，N 越小.这说明较大的横截面网格步进会导致较大的误差，而所提的迭代方法相比 CN-PE 方法和 ADI-PE 方法能带来较低的误差效果.此外，传播仰角角度越大，迭代方法降低误差效果越明显.42681000.010.020.030.040.060.05相对误差NCN-PE 方法(=3)ADI-PE 方法(=3)三维 PE 迭代方法(=3)CN-PE 方法(=6)ADI-PE 方法(=6)三维 PE 迭代方法(=6)图 2 N 不同时三种计算方法的相对误差Fig.2 Relative error vs.N by 3 methods =4 N=5 =6，时不同方位角对三种方法的影响如图 3 所示.可以看出，三种方法在方位角 90和 180上达到相同的计算误差，但在 45和 135方位角上，本文推导的三维 PE 迭代方法能在对应的角度上降低计算误差.第 1 期吴晓萍，等：一种求解三维抛物方程的迭代方法161 450901351803.03.54.04.5相对误差/()CN-PE 方法ADI-PE 方法三维 PE 迭代方法103图 3 方位角 不同时三种计算方法的相对误差Fig.3 Relative error vs.azimuth angle by 3 methods N=5 =45 =6xxx，时相对误差与 关系如图 4所示.可以看出，三维 PE 迭代方法在 较大时，能够减少相对误差.这是因为 与传播方向上的网格步进大小成正比关系，越大，CN-PE 方法和 ADI-PE 方法计算出的误差越大；而三维 PE 迭代方法考虑了分裂项的补偿，分裂项与网格步进成正比关系.6214182610223000.0050.0100.0150.020相对误差CN-PE 方法ADI-PE 方法三维 PE 迭代方法图 4 不同时三种计算方法的相对误差Fig.4 Relative error vs.by 3 methods =15 N=5 =45，时本文提出的求解三维 PE迭代方法的相对误差如图 5 所示.可以看出：迭代次数为 1 时，误差降低效果明显；迭代次数为 2 时，迭代方法能够再次降低误差；当迭代一定次数时，所能降低误差的效果趋向于相等.说明本文所提出的三维 PE 迭代方法在一定迭代次数下能降低计算误差.1023450.020.040.060.08相对误差迭代次数图 5 三维 PE 迭代方法的相对误差Fig.5 Relative error by 3D PE interative method=15N=5利用一个具体算例对本文所提方法性能进行验证.将三维 PE 迭代方法、ADI-PE 方法应用在一个密闭的隧道空间，隧道横截面规格为 80 80.工作频率为 1 GHz、半功率波瓣宽度为 15的水平极化高斯源放置在隧道口正中间位置，离隧道口 50 m 处横截面的场强值仿真结果如图 6 所示.用欧几里得误差范数计算相应的数值误差，分别为15.2%(ADI-PE 方法)、13.6%(1 次迭代)、12.8%(2 次迭代)，可以看出，所提的三维 PE 迭代方法计算场强值更接近参考值.60121824601218240.70.60.50.40.30.20.10z/my/m(a)ADI-PE 方法(a)ADI-PE methodV/m60121824601218240.70.60.50.40.30.20.10z/my/m(c)三维 PE 迭代方法(1 次迭代)(c)3D PE iterative method(1 iteration)V/m60121824601218240.70.60.50.40.30.20.10z/my/m(d)三维 PE 迭代方法(2 次迭代)(d)3D PE iterative method(2 iteration)V/m60121824601218240.70.60.50.40.30.20.10z/my/m(b)理论参考(b)Theoretical referenceV/m图 6 两种方法场强仿真结果Fig.6 Simulation result comparison of field strength by 2methods 此外，三维 PE 迭代方法与 ADI-PE、CN-PE 的计算效率对比如表 1 所示.可以看出，虽然三维 PE 迭代方法比 ADI-PE 方法计算用时较多，但其时间增量在一定的线性增加范围，远低于 CN-PE 计算量.即在少量的时间开销下，三维 PE 迭代方法能够提高PE 求解的准确度.表 1 三种方法计算效率对比Tab.1 Comparison of computational efficiency方法CN-PE方法 ADI-PE方法三维PE迭代方法(1次迭代)三维PE迭代法(2次迭代)计算时间/s585.6281.7283.1514.028 3 结论本文提出一种求解三维 PE 的迭代方法.对比现有的求解三维 PE 的 ADI-PE 和 CN-PE 方法，该方法能够提高计算精准度，且可通过迭代达到一定计算 162电波科学学报第 38 卷精度，是一种高效且计算精准度高的三维 PE 求解方法.但该方法未考虑阻抗边界处理问题，是下一步研究的方向.参考文献 LEVY,M F.Parabolic equation methods for electromag-netic wave propagation M.London:IEE Press,2000.1 APAYDIN G,SEVGI L.Radio wave propagation and para-bolic equation modeling M.Wiley-IEEE Press,2017.2 黎杨.抛物型方程的有限差分解法及其在复杂电磁环境中的应用D.武汉:武汉理工大学,2010.LI Y.The finite-difference method to parabolic equationand its application in complex electromagnetic environ-ments D.Wuhan:Wuhan University of Technology,2010.(in Chinese)3 郭建炎,王剑莹,龙云亮.森林中电波传播的抛物方程法J.电波科学学报，2007，22（6）：1042-1046.GUO J Y,WANG J Y,LONG Y L.Parabolic equationmodel for wave propagation in forest environmentsJ.Chinese journal of radio science，2007，22（6）：1042-1046.（in Chinese）4 GUO Q,ZHOU C,LONG Y L.Greene approximationwide-angle parabolic equation for radio propagationJ.IEEE transactions on antennas and propagation，2017，65（11）：6048-6056.5 GUO Q,LONG Y L.Pade second-order parabolic equa-tion modeling for propagation over irregular terrainJ.IEEE antennas and wireless propagation letters，2017：2852-2855.6 郭琪,黎子豪,朱琼琼,等.基于全局透射边界条件的宽角抛物方程电波预测模型J.电波科学学报，2019，34（4）：473-478.GUO Q,LI Z H,ZHU Q Q,et al.Radio propagation pre-diction model based on wide-angle parabolic equation withnon-local boundary conditionJ.Chinese journal of radioscience，2019，34（4）：473-478.（in Chinese）7 LYTAEV M S,VLADYKO A G.Split-step Pad approx-imations of the Helmholtz equation for radio coverage pre-diction over irregular terrain C/Advances in Wirelessand Optical Communications(RTUWO),2018.8 LYTAEV M S.Numerov-Pad scheme for the one-wayHelmholtz equation in tropospheric radio-wave propaga-9tionJ.IEEE antennas and wireless propagation letters，2020，19（12）：2167-2171.ZHOU L,WANG D,MU Z,et al.LF Radio wave predic-tion at short ranges with high propagation angles over irreg-ular terrainJ.IEEE antennas and wireless propagation let-ters，2016，16：732-735.10 WANG D D,PU Y R,XI X L,et al.An analysis of narrow-angle and wide-angle shift-map parabolic equationsJ.IEEE transactions on antennas and propagation，2020，68（5）：3911-3918.11 MARTELLY R,JANASWAMY R.An ADI-PE approachfor modeling radio transmission loss in tunnelsJ.IEEEtransactions on antennas and propagation，2009，57（6）：1759-1770.12 HE Z,ZENG H,CHEN R S.Two way propagation model-ing of expressway with vehicles by using the three-dimensionalADI-PE methodJ.IEEE transactions on antennas andpropagation，2018，66（4）：2156-2160.13 HUANG L,WU X P,LI Z H,et al.A parallel FDTD/ADI-PE method for ultralarge-scale propagation modeling of ILSsignal analysisJ.IEEE antennas and wireless propaga-tion letters，2020，19（12）：2245-2249.14 WANG S,TEIXEIRA F L,JI C.An iterative ADI-FDTDwith reduced splitting errorJ.IEEE microwave and wire-less components letters，2005，15（2）：92-94.15作者简介吴晓萍（1994），女，广东人，中山大学在读博士生，研究方向为计算电磁学、抛物方程、无线通信应用.E-mail:梁志禧（1989），男，广东人，中山大学副教授，研究方向为天线理论与设计、电波传播理论和电磁数值计算等.E-mail:黎子豪（1991），男，广东人，中山大学在读博士生，研究方向为计算电磁学、时域有限差分方程和电磁应用.E-mail:龙云亮(1963)，男，重庆人，中山大学教授，博士生导师，研究方向为天线理论与设计、电波传播理论和电磁数值计算等.E-mail:第 1 期吴晓萍，等：一种求解三维抛物方程的迭代方法163