CN107341284A - 高精度预测低频电波传播特性的双向抛物方程方法 - Google Patents

高精度预测低频电波传播特性的双向抛物方程方法 Download PDF

Info

Publication number
CN107341284A
CN107341284A CN201710350234.2A CN201710350234A CN107341284A CN 107341284 A CN107341284 A CN 107341284A CN 201710350234 A CN201710350234 A CN 201710350234A CN 107341284 A CN107341284 A CN 107341284A
Authority
CN
China
Prior art keywords
mrow
msub
field
rho
backward
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Granted
Application number
CN201710350234.2A
Other languages
English (en)
Other versions
CN107341284B (zh
Inventor
席晓莉
王丹丹
张金生
蒲玉蓉
李征委
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Xian University of Technology
Original Assignee
Xian University of Technology
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Xian University of Technology filed Critical Xian University of Technology
Priority to CN201710350234.2A priority Critical patent/CN107341284B/zh
Publication of CN107341284A publication Critical patent/CN107341284A/zh
Application granted granted Critical
Publication of CN107341284B publication Critical patent/CN107341284B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Hardware Design (AREA)
  • Evolutionary Computation (AREA)
  • Geometry (AREA)
  • General Engineering & Computer Science (AREA)
  • General Physics & Mathematics (AREA)
  • Geophysics And Detection Of Objects (AREA)
  • Aerials With Secondary Devices (AREA)

Abstract

本发明公开了高精度预测低频电波传播特性的双向抛物方程方法,具体按照以下步骤实施:步骤1:输入模型文件;步骤2:利用平地面公式计算前向初始场;步骤3:结合前向初始场,基于坐标变换模型,采用分布离散傅里叶变换算法,求解计算区域任意位置电波传播的前向场;步骤4:结合前向初始场,基于阶梯近似模型,采用SSFT算法,递归的求解由于地形影响对电波传播产生的后向场,和由多次反射产生的前向场;步骤5:结合步骤3和步骤4中的传播场结果求出电波传播的总磁场。本发明既能解决FDTD方法计算量大、耗时长的缺点,又能解决积分方程方法和抛物方程方法在地形起伏剧烈路径中由于忽略后向电波传播影响引起的误差大的问题。

Description

高精度预测低频电波传播特性的双向抛物方程方法
技术领域
本发明属于电波传播技术领域,具体涉及一种高精度预测低频电波传播特性的双向抛物方程方法。
背景技术
低频电波以其波长较长,信号的传播损耗小,信号的幅度与相位稳定,被广泛应用于授时、导航、通信等领域。为提高低频无线电工程的使用性能和精度,需要深入研究低频电波在各种地域、时域、频域范围内的变化规律及预测技术。目前,预测复杂环境下低频电波传播已有的理论方法有:积分方程方法、抛物方程方法和时域有限差分(Finite-Difference Time-Domain, FDTD)方法。积分方程方法和抛物方程方法都可以考虑地形变换,对复杂路径具有较高精度,但是,这两种方法都忽略了沿路径后向传播的反射场影响,在地形起伏剧烈的传播路径上会引起较大误差。FDTD方法可方便的考虑任何地质和地形路径的电波传播问题,并且精度较高,但其存在过长的计算时间、较多的内存资源占用等问题,不适合大区域的工程化推广。
发明内容
本发明的目的是提供一种高精度预测低频电波传播特性的双向抛物方程方法,既能解决FDTD方法计算量大、耗时长的缺点,又能解决积分方程方法和抛物方程方法在地形起伏剧烈路径中由于忽略后向电波传播影响引起的误差大的问题。
本发明所采用的技术方案是:高精度预测低频电波传播特性的双向抛物方程方法,具体按照以下步骤实施:
步骤1:输入模型文件;
步骤2:利用平地面公式计算初始位置ρ=ρ0处纵向横切面上的磁场 Hf0,z),并通过其求解双向抛物方程的前向初始场uf0,z);
步骤3:结合步骤2计算出来的前向初始场uf0,z),基于坐标变换模型,采用分布离散傅里叶变换算法,求解计算区域任意位置电波传播的前向场uf(ρ,z);
步骤4:结合步骤2计算出来的前向初始场uf0,z),基于阶梯近似模型,采用SSFT算法,递归的求解由于地形影响对电波传播产生的后向场和由多次反射产生的前向场τ=1,2,3,…;
步骤5:结合步骤3和步骤4中的传播场结果uf(ρ,z)、τ=1,2,3,…和求出电波传播的总磁场
本发明的特点还在于:
步骤1具体为:
计算区域大小Nρ×Nz,其中Nρ为ρ方向的网格数,Nz为z方向的网格数;空间网格步长分别为Δρ和Δz,ρ和z分别表示横向和纵向坐标;ρ0为初始场位置;源的参数;地面相对介电常数εr和电导率σ;地形模型参数。
步骤2具体为:
选取二维柱坐标系(ρ,z),其中ρ和z分别表示为横向和纵向坐标,根据实际发射天线尺寸,通过测量得到垂直电偶极子的电荷间距dl,放置在距离地面高度为d的位置,假设时谐因子为e-iωt,利用平地面公式计算初始距离ρ=ρ0处纵向横切面的磁场Hf0,z):
其中,I为电流大小,k0和kg分别为真空和地面的波数,r0表示从源点到观测点的直线距离,r1表示从源的镜像到观测点的直线距离,P2为中间参量:
F(z)是Fresnel积分并定义为
定义沿ρ轴正方向传播的波函数,即前向场uf0,z)为:
步骤3具体为:
在原坐标系(ρ,z)中,地形函数为T(ρ),以地形表面上任意一点为新的坐标原点O′建立新的坐标系
在新坐标系中,波函数已经修正为:
其中,
采用起伏地形SSFT算法求解下一步进处的前向场:
其中,n为大气折射率,为离散混合傅里叶变换对,p=k0sinθ,θ为电波传播的仰角;
步骤2中得到初始距离ρ=ρ0处电波传播的前向初始场为uf0,z),结合式(4)和(5),即可得到新坐标系下初始场
代入式(6)中步进计算,得到新坐标系下任意位置处的前向场
再考虑原坐标系和新坐标系前向场之间的关系,即式(4),求出原坐标系中前向场uf(ρ,z)。
步骤4具体为:
沿ρ轴正方向传播的各个前向场为:
其中,上标τ表示前向场的序号,τ=0表示电波传播的主传播场,而τ=1,2,3…则表示由后向场多次反射所产生的前向场,为前向磁场分量;
当电波沿ρ轴正方向传播时,采用SSFT算法即可求得下一位置处的前向场:
将步骤2中求得的前向初始场uf0,z)赋值给若传播路径不含起伏地形,则计算区域仅存在τ=0时对应的代入式(8)依次迭代计算出各个位置的场若存在起伏地形,需要考虑后向场影响;
定义沿ρ轴负方向传播的各个后向波函数,即后向场:
其中,表示后向场的序号,每当前向场遇到一个相对于其传播方向的上升阶梯面,就会产生一个后向磁场分量
对应于式(9)的SSFT算法为:
假设前向波沿ρ轴正方向传播到ρ=ρe处存在相对于传播方向的上升阶梯面,由于该上升阶梯面对电波传播的反射作用,前向波分为两部分:继续向前传播的前向波和反射作用产生的后向波;
在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向初始场为:
修正前向场:
其中,τ和都是具体的某一个标号,具体过程为标记为τ的前向波,每遇到一个相对于ρ正方向的上升沿就会产生一个标记为的后向波,并且更新为
假设电波沿ρ轴负方向传播到ρ=ρt处存在相对于传播方向的上升阶梯面,后向波分为两部分:继续向后传播的后向波和由后向波反射产生的前向波;
在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向场多次反射产生的前向场为:
修正后向场:
其中,τ和都是具体的某一个标号。具体过程为标记为的后向波,每遇到一个相对于ρ负方向的上升沿就会产生一个标记为τ的前向波,并且τ更新为τ=τ+1;
转到式(7)和式(9)继续循环对所有前向场和后向场迭代计算下一步进处的场强,直至步进到计算区域边界;其中,前向场的边界为ρmax=NρΔρ,后向场的边界为ρ0
步骤4中循环时设置门限终止计算,循环计算的判定依据为:
其中,u(ρ,z)表示τ=1,2,3…或者判决系数ξ设置为0.01,|| ||为矢量求模。
步骤5具体为:
将步骤3和步骤4中的前向场uf(ρ,z)和多次反射产生的前向场τ=1,2,3,…分别代入式(2)和式(7),分别求出所有存在的前向场Hf(ρ,z)和
将步骤4计算的代入式(9)中,求出所有存在的后向磁场分量
电波传播的总磁场为所有前向磁场分量和后向磁场分量的叠加:
本发明的有益效果是:
①本发明高精度预测低频电波传播特性的双向抛物方程方法,是基于柱坐标系实现的,因此,能够更为方便的处理圆柱对称的电波传播问题,对工程指导比较有意义;
②本发明高精度预测低频电波传播特性的双向抛物方程方法,能够准确的预测复杂环境中电波传播的前向和后向传播影响,克服了积分方程方法和传统的抛物方程方法忽略后向波影响引起的误差。此外,与FDTD方法相比,在相同的计算精度前提下,该方法基于高效的SSFT算法步进求解场强,大大的缩短了计算时间,更适合工程推广。
附图说明
图1是本发明方法的流程图;
图2(a)是高度为0.5km时本发明方法与FDTD方法的地面场强比较;
图2(b)是高度为1.5km时本发明方法与FDTD方法的地面场强比较;
图3是在多个高斯山脉地形情况下,本发明方法与FDTD方法的地面场强比较图;
图4(a)是FDTD方法的空间场强伪彩色图比较;
图4(b)是本发明方法的空间场强伪彩色图比较。
具体实施方式
下面结合附图和具体实施方式对本发明进行详细说明。
本发明提供了高精度预测低频电波传播特性的双向抛物方程方法,改进了传统抛物方程方法在复杂地形传播路径上由于忽略后向传播影响引起的预测误差,具体实现过程是:首先利用平地面公式计算出前向初始场;接着,基于坐标变化模型,采用抛物方程方法预测整个区域在忽略后向传播影响情况下的前向场的场值;然后,基于阶梯近似模型,通过双向抛物方程方法预测由于电波在复杂地形作用下引起的反射场(包括后向场和多次反射产生的前向场)影响。最后,通过抛物方程方法计算的前向场和双向抛物方程方法计算的后向场和多次反射产生的前向场叠加求解总磁场,具体实施步骤流程图见图1。
在采用双向抛物方程方法求圆柱对称模型问题中的前向场和后向场时,首先需要推导柱坐标系下的双向抛物方程及其对应的SSFT算法:
假设电磁波的时谐因子是e-iωt,取二维的柱坐标系(ρ,z),其中ρ和z分别表示为横向和纵向坐标,并以标量Φ表示任意标量场分量,假设Φ与方位角无关,则对于水平极化波,而对于垂直极化波,Φ满足二维波动方程
其中k0是真空中的波数,n为大气折射率。首先,作如下变量代换:
其中,ψ(ρ,z)仅考虑了电磁波的柱面传播影响,没有涉及相位修正。将式(18)代入(17)中,可以得到:
然后,对上式作远场近似,即假设ρ>>1且k0ρ>>1,式(19)变为:
当折射指数n随ρ的变化很小时,式(20)可分解为
由于式(21)中的两个因式互不相关,故可以得到关于ρ的两个独立的式子
ψf和ψb分别对应于前向波和后向波。接下来,考虑相位修正,引入合适的相位因子,对ψf和ψb分别作如下修正:
其中,φf和φb分别代表前向和后向场分量。将式(18)代入式(24)和 (25)中,可以分别得到
其中,Φf和Φb分别为前向和后向磁场分量。
下面将式(24)和(25)分别代入(22)和(23),可以得到双向抛物方程:
上面两式分别称为前向和后向抛物方程。采用SSFT算法求解上面两式,可得:
式中为离散混合傅里叶变换对,为大气折射率,Δρ为ρ方向的网格大小。实际计算中,φf和φb可能分别为一系列的其中τ=0,1,2,3,…和
本发明高精度预测低频电波传播特性的双向抛物方程方法,具体按照以下步骤实施:
步骤1:输入模型文件,具体为:
计算区域大小Nρ×Nz,其中Nρ为ρ方向的网格数,Nz为z方向的网格数;空间网格步长分别为Δρ和Δz,ρ和z分别表示横向和纵向坐标;ρ0为初始场位置;源的参数;地面相对介电常数εr和电导率σ;地形模型参数。
步骤2:利用平地面公式计算初始位置ρ=ρ0处纵向横切面上的磁场 Hf0,z),并通过其求解双向抛物方程的前向初始场uf0,z),具体为:
为了更好地解决圆柱对称模型问题,在此选取二维柱坐标系(ρ,z),其中ρ和z分别表示为横向和纵向坐标,根据实际发射天线尺寸,通过测量得到垂直电偶极子的电荷间距dl,放置在距离地面高度为d的位置,假设时谐因子为e-iωt,利用平地面公式计算初始距离ρ=ρ0处纵向横切面的磁场 Hf0,z):
其中,I为电流大小,k0和kg分别为真空和地面的波数,r0表示从源点到观测点的直线距离,r1表示从源的镜像到观测点的直线距离,P2为中间参量:
F(z)是Fresnel积分并定义为
定义沿ρ轴正方向传播的波函数,即前向场uf0,z)为:
步骤3:结合步骤2计算出来的前向初始场uf0,z),基于坐标变换模型,采用分布离散傅里叶变换(Split-Step Fourier Transform,SSFT)算法,求解计算区域任意位置电波传播的前向场uf(ρ,z),具体为:
在原坐标系(ρ,z)中,地形函数为T(ρ),以地形表面上任意一点为新的坐标原点O′建立新的坐标系
在新坐标系中,波函数已经修正为:
其中,
采用起伏地形SSFT算法求解下一步进处的前向场:
其中,n为大气折射率,为离散混合傅里叶变换对,p=k0sinθ,θ为电波传播的仰角;
步骤2中得到初始距离ρ=ρ0处电波传播的前向初始场为uf0,z),结合式(4)和(5),即可得到新坐标系下初始场
代入式(6)中步进计算,得到新坐标系下任意位置处的前向场
再考虑原坐标系和新坐标系前向场之间的关系,即式(4),求出原坐标系中前向场uf(ρ,z)。
步骤4:结合步骤2计算出来的前向场uf0,z),基于阶梯近似模型,采用SSFT算法,递归的求解由于地形影响对电波传播产生的后向场和由多次反射产生的前向场τ=1,2,3,…,具体为:
沿ρ轴正方向传播的各个前向场为:
其中,上标τ表示前向场的序号,τ=0表示电波传播的主传播场,而τ=1,2,3…则表示由后向场多次反射所产生的前向场,为前向磁场分量;与式(2)相比,可以看出,前向场的磁场求解具有相同的表达形式。
当电波沿ρ轴正方向传播时,采用SSFT算法即可求得下一位置处的前向场:
将步骤2中求得的前向初始场uf0,z)赋值给若传播路径不含起伏地形,则计算区域仅存在τ=0时对应的代入式 (8)依次迭代计算出各个位置的场若存在起伏地形,就需要考虑后向场影响;
下面定义沿ρ轴负方向传播的各个后向波函数,即后向场:
其中,表示后向场的序号,每当前向场遇到一个相对于其传播方向的上升阶梯面,就会产生一个后向磁场分量
对应于式(9)的SSFT算法为:
比较式(8)和(10)可以得出,前向场和后向场的SSFT算法迭代形式基本一致。
假设前向波沿ρ轴正方向传播到ρ=ρe处存在相对于传播方向的上升阶梯面,由于该上升阶梯面对电波传播的反射作用,前向波分为两部分:继续向前传播的前向波和反射作用产生的后向波。在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向初始场为:
修正前向场:
这里的τ和都是具体的某一个标号。具体过程为标记为τ的前向波,每遇到一个相对于ρ正方向的上升沿就会产生一个标记为的后向波,并且更新为
同理,假设电波沿ρ轴负方向传播到ρ=ρt处存在相对于传播方向的上升阶梯面,后向波分为两部分:继续向后传播的后向波和由后向波反射产生的前向波。在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向场多次反射产生的前向场为:
修正后向场:
这里的τ和都是具体的某一个标号。具体过程为标记为的后向波,每遇到一个相对于ρ负方向的上升沿就会产生一个标记为τ的前向波,并且τ更新为τ=τ+1。
然后,采用式(7)和(9)继续对所有前向场和后向场迭代计算下一步进处的场强,接着判断是否存相对于传播方向上的上升阶梯面,具体执行过程见上述假设,循环此过程直至步进到计算区域边界(对于前向场,边界为ρmax=NρΔρ,对于后向场,边界为ρ0)。
由此可以看出对于孤立山峰,电波向前传播遇到山峰反射会产生许多后向波,而后向波不存在反射,即不会产生前向波;而对于多个山峰,当电波向前传播遇到第一个山峰时,一部分电波会由于反射作用而产生后向波,另一部分电波则会继续向前传播,当遇到第二个山峰时,又会出现后向波,并且,后向波会在两个山峰之间往复反射,以此类推,从而产生一系列的反射波(其中包括多次反射产生的前向波和后向波)。由于受空间散射及介质的吸收影响,循环计算的反射波场强逐渐减小,继而不会对全域场产生影响。因此,有必要设置门限来终止计算。为了方便计算,给出递归计算的判定依据为:
其中,u(ρ,z)表示τ=1,2,3…或者判决系数ξ设置为0.01,|| ||为矢量求模。
步骤5:结合步骤3和步骤4中的传播场结果uf(ρ,z)、 τ=1,2,3,…,求出电波传播的总磁场具体为:
将步骤3和步骤4中的前向场uf(ρ,z)和多次反射产生的前向场τ=1,2,3,…分别代入式(2)和式(7),分别求出所有存在的前向场Hf(ρ,z)和注意:这里没有采用而是用步骤3求出来的Hf(ρ,z)替代,这是因为采用阶梯近似模型求出来的前向场分量精度在山后不够。
将步骤4计算的代入式(9)中,求出所有存在的后向磁场分量
电波传播的总磁场为所有前向磁场分量和后向磁场分量的叠加:
实施例1
单个高斯山脉地形地面场强预测
垂直电偶极子天线的辐射功率均采用1kW,信号源频率为100kHz。计算区域总大小为ρmax:100km×zmax:102.4km,网格剖分大小均分别为 dρ=200m和dz=100m,初始距离ρ0=10km;地面电参数为εr=13,σ=3×10-3S/m(陆地);在中心位置ρc=50km处有一孤立高斯山峰,其高度函数为
l为山峰宽度取2km,H为山峰高度且分别取0.5km和1.5km。图2(a) 和图2(b)分别为在相同宽度不同高度的单个高斯山脉地形情况下,本发明方法与FDTD方法的地面场强比较。由图2(a)和图2(b)可以看出,两种方法的计算结果吻合一致,并且都可以预测到山前电波反射影响的振荡存在。但是,相比于FDTD方法,双向抛物方程方法计算更快。
实施例2
多个高斯山脉地形地面场强预测
将实施例1中的山峰改为多个高斯山脉地形,第一个山峰的高度为1km,宽度为4km,中心位置为40km;第二个山峰的高度为1.5km,宽度为2km,中心位置为60km,其他参数不变。图3给出了在多个高斯山脉地形情况下,本发明方法与FDTD方法的地面场强比较。由图3可以看出,对于多个高斯山脉地形,双向抛物方程方法也能精确的预测山峰之间波的往复反射影响。并且经统计,FDTD方法的计算时间是双向抛物方程方法的21倍。图4(a) 和图4(b)分别为FDTD方法与本发明方法的空间场强伪彩色图比较,从图中可以直观的看到由于山峰的影响,电波在山前和山峰之间的反射影响。

Claims (7)

1.高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,具体按照以下步骤实施:
步骤1:输入模型文件;
步骤2:利用平地面公式计算初始位置ρ=ρ0处纵向横切面上的磁场Hf0,z),并通过其求解双向抛物方程的前向初始场uf0,z);
步骤3:结合步骤2计算出来的前向初始场uf0,z),基于坐标变换模型,采用分布离散傅里叶变换算法,求解计算区域任意位置电波传播的前向场uf(ρ,z);
步骤4:结合步骤2计算出来的前向初始场uf0,z),基于阶梯近似模型,采用SSFT算法,递归的求解由于地形影响对电波传播产生的后向场和由多次反射产生的前向场
步骤5:结合步骤3和步骤4中的传播场结果uf(ρ,z)、 求出电波传播的总磁场
2.根据权利要求1所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤1具体为:
计算区域大小Nρ×Nz,其中Nρ为ρ方向的网格数,Nz为z方向的网格数;空间网格步长分别为Δρ和Δz,ρ和z分别表示横向和纵向坐标;ρ0为初始场位置;源的参数;地面相对介电常数εr和电导率σ;地形模型参数。
3.根据权利要求1所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤2具体为:
选取二维柱坐标系(ρ,z),其中ρ和z分别表示为横向和纵向坐标,根据实际发射天线尺寸,通过测量得到垂直电偶极子的电荷间距dl,放置在距离地面高度为d的位置,假设时谐因子为e-iωt,利用平地面公式计算初始距离ρ=ρ0处纵向横切面的磁场Hf0,z):
<mrow> <mtable> <mtr> <mtd> <mrow> <msub> <mi>H</mi> <mi>f</mi> </msub> <mrow> <mo>(</mo> <mrow> <msub> <mi>&amp;rho;</mi> <mn>0</mn> </msub> <mo>,</mo> <mi>z</mi> </mrow> <mo>)</mo> </mrow> <mo>=</mo> <mo>-</mo> <mfrac> <mrow> <mi>I</mi> <mi>d</mi> <mi>l</mi> </mrow> <mrow> <mn>2</mn> <mi>&amp;pi;</mi> </mrow> </mfrac> <mo>&amp;lsqb;</mo> <mfrac> <msup> <mi>e</mi> <mrow> <msub> <mi>ik</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>0</mn> </msub> </mrow> </msup> <mn>2</mn> </mfrac> <mfrac> <msub> <mi>&amp;rho;</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>0</mn> </msub> </mfrac> <mrow> <mo>(</mo> <mrow> <mfrac> <mrow> <msub> <mi>ik</mi> <mn>0</mn> </msub> </mrow> <msub> <mi>r</mi> <mn>0</mn> </msub> </mfrac> <mo>-</mo> <mfrac> <mn>1</mn> <msubsup> <mi>r</mi> <mn>0</mn> <mn>2</mn> </msubsup> </mfrac> </mrow> <mo>)</mo> </mrow> <mo>+</mo> <mfrac> <msup> <mi>e</mi> <mrow> <msub> <mi>ik</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>1</mn> </msub> </mrow> </msup> <mn>2</mn> </mfrac> <mfrac> <msub> <mi>&amp;rho;</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>1</mn> </msub> </mfrac> <mrow> <mo>(</mo> <mrow> <mfrac> <mrow> <msub> <mi>ik</mi> <mn>0</mn> </msub> </mrow> <msub> <mi>r</mi> <mn>1</mn> </msub> </mfrac> <mo>-</mo> <mfrac> <mn>1</mn> <msubsup> <mi>r</mi> <mn>1</mn> <mn>2</mn> </msubsup> </mfrac> </mrow> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mo>-</mo> <mfrac> <mrow> <msup> <mi>e</mi> <mrow> <msub> <mi>ik</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>1</mn> </msub> </mrow> </msup> <msubsup> <mi>k</mi> <mn>0</mn> <mn>3</mn> </msubsup> </mrow> <msub> <mi>k</mi> <mi>g</mi> </msub> </mfrac> <msqrt> <mfrac> <mi>&amp;pi;</mi> <mrow> <msub> <mi>k</mi> <mn>0</mn> </msub> <msub> <mi>r</mi> <mn>1</mn> </msub> </mrow> </mfrac> </msqrt> <msup> <mi>e</mi> <mrow> <mo>-</mo> <msub> <mi>iP</mi> <mn>2</mn> </msub> </mrow> </msup> <mi>F</mi> <mrow> <mo>(</mo> <msub> <mi>P</mi> <mn>2</mn> </msub> <mo>)</mo> </mrow> <mo>&amp;rsqb;</mo> </mrow> </mtd> </mtr> </mtable> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>1</mn> <mo>)</mo> </mrow> </mrow>
其中,I为电流大小,k0和kg分别为真空和地面的波数,r0表示从源点到观测点的直线距离,r1表示从源的镜像到观测点的直线距离,P2为中间参量:F(z)是Fresnel积分并定义为
定义沿ρ轴正方向传播的波函数,即前向场uf0,z)为:
<mrow> <msub> <mi>u</mi> <mi>f</mi> </msub> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <msqrt> <mrow> <msub> <mi>k</mi> <mn>0</mn> </msub> <mi>&amp;rho;</mi> </mrow> </msqrt> <msup> <mi>e</mi> <mrow> <mo>-</mo> <msub> <mi>ik</mi> <mn>0</mn> </msub> <mi>&amp;rho;</mi> </mrow> </msup> <msub> <mi>H</mi> <mi>f</mi> </msub> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>2</mn> <mo>)</mo> </mrow> <mo>.</mo> </mrow>
4.根据权利要求3所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤3具体为:
在原坐标系(ρ,z)中,地形函数为T(ρ),以地形表面上任意一点为新的坐标原点O′建立新的坐标系
在新坐标系中,波函数已经修正为:
其中,
采用起伏地形SSFT算法求解下一步进处的前向场:
其中,n为大气折射率,为离散混合傅里叶变换对, 为电波传播的仰角;
步骤2中得到初始距离ρ=ρ0处电波传播的前向初始场为uf0,z),结合式(4)和(5),即可得到新坐标系下初始场
代入式(6)中步进计算,得到新坐标系下任意位置处的前向场
再考虑原坐标系和新坐标系前向场之间的关系,即式(4),求出原坐标系中前向场uf(ρ,z)。
5.根据权利要求4所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤4具体为:
沿ρ轴正方向传播的各个前向场为:
<mrow> <msubsup> <mi>u</mi> <mi>f</mi> <mi>&amp;tau;</mi> </msubsup> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <msqrt> <mrow> <msub> <mi>k</mi> <mn>0</mn> </msub> <mi>&amp;rho;</mi> </mrow> </msqrt> <msup> <mi>e</mi> <mrow> <mo>-</mo> <msub> <mi>ik</mi> <mn>0</mn> </msub> <mi>&amp;rho;</mi> </mrow> </msup> <msubsup> <mi>H</mi> <mi>f</mi> <mi>&amp;tau;</mi> </msubsup> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>7</mn> <mo>)</mo> </mrow> </mrow>
其中,上标τ表示前向场的序号,τ=0表示电波传播的主传播场,而τ=1,2,3…则表示由后向场多次反射所产生的前向场,为前向磁场分量;
当电波沿ρ轴正方向传播时,采用SSFT算法即可求得下一位置处的前向场:
将步骤2中求得的前向初始场uf0,z)赋值给若传播路径不含起伏地形,则计算区域仅存在τ=0时对应的代入式(8)依次迭代计算出各个位置的场若存在起伏地形,需要考虑后向场影响;
定义沿ρ轴负方向传播的各个后向波函数,即后向场:
其中,表示后向场的序号,每当前向场遇到一个相对于其传播方向的上升阶梯面,就会产生一个后向磁场分量
对应于式(9)的SSFT算法为:
假设前向波沿ρ轴正方向传播到ρ=ρe处存在相对于传播方向的上升阶梯面,由于该上升阶梯面对电波传播的反射作用,前向波分为两部分:继续向前传播的前向波和反射作用产生的后向波;
在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向初始场为:
修正前向场:
<mrow> <msubsup> <mi>u</mi> <mi>f</mi> <mi>&amp;tau;</mi> </msubsup> <mrow> <mo>(</mo> <msub> <mi>&amp;rho;</mi> <mi>e</mi> </msub> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>=</mo> <mfenced open = "{" close = ""> <mtable> <mtr> <mtd> <mrow> <msubsup> <mi>u</mi> <mi>f</mi> <mi>&amp;tau;</mi> </msubsup> <mrow> <mo>(</mo> <msub> <mi>&amp;rho;</mi> <mi>e</mi> </msub> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>,</mo> </mrow> </mtd> <mtd> <mrow> <mi>z</mi> <mo>&gt;</mo> <mi>t</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;rho;</mi> <mi>e</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> <mtr> <mtd> <mrow> <mn>0</mn> <mo>,</mo> </mrow> </mtd> <mtd> <mrow> <mi>z</mi> <mo>&lt;</mo> <mi>t</mi> <mrow> <mo>(</mo> <msub> <mi>&amp;rho;</mi> <mi>e</mi> </msub> <mo>)</mo> </mrow> </mrow> </mtd> </mtr> </mtable> </mfenced> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>12</mn> <mo>)</mo> </mrow> </mrow>
其中,τ和都是具体的某一个标号,具体过程为标记为τ的前向波,每遇到一个相对于ρ正方向的上升沿就会产生一个标记为的后向波,并且更新为
假设电波沿ρ轴负方向传播到ρ=ρt处存在相对于传播方向的上升阶梯面,后向波分为两部分:继续向后传播的后向波和由后向波反射产生的前向波;
在上升阶梯面上,根据磁场所满足切向边界条件可以得到后向场多次反射产生的前向场为:
修正后向场:
其中,τ和都是具体的某一个标号。具体过程为标记为的后向波,每遇到一个相对于ρ负方向的上升沿就会产生一个标记为τ的前向波,并且τ更新为τ=τ+1;
转到式(7)和式(9)继续循环对所有前向场和后向场迭代计算下一步进处的场强,直至步进到计算区域边界;其中,前向场的边界为ρmax=NρΔρ,后向场的边界为ρ0
6.根据权利要求5所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤4中循环时设置门限终止计算,循环计算的判定依据为:
<mrow> <mo>|</mo> <mo>|</mo> <mi>u</mi> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>|</mo> <mo>|</mo> <mo>/</mo> <mo>|</mo> <mo>|</mo> <msubsup> <mi>u</mi> <mi>f</mi> <mn>0</mn> </msubsup> <mrow> <mo>(</mo> <mi>&amp;rho;</mi> <mo>,</mo> <mi>z</mi> <mo>)</mo> </mrow> <mo>|</mo> <mo>|</mo> <mo>&lt;</mo> <mi>&amp;xi;</mi> <mo>-</mo> <mo>-</mo> <mo>-</mo> <mrow> <mo>(</mo> <mn>15</mn> <mo>)</mo> </mrow> </mrow>
其中,u(ρ,z)表示或者判决系数ξ设置为0.01,||||为矢量求模。
7.根据权利要求5所述的高精度预测低频电波传播特性的双向抛物方程方法,其特征在于,所述步骤5具体为:
将步骤3和步骤4中的前向场uf(ρ,z)和多次反射产生的前向场分别代入式(2)和式(7),分别求出所有存在的前向场Hf(ρ,z)和
将步骤4计算的代入式(9)中,求出所有存在的后向磁场分量
电波传播的总磁场为所有前向磁场分量和后向磁场分量的叠加:
CN201710350234.2A 2017-05-18 2017-05-18 高精度预测低频电波传播特性的双向抛物方程方法 Active CN107341284B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN201710350234.2A CN107341284B (zh) 2017-05-18 2017-05-18 高精度预测低频电波传播特性的双向抛物方程方法

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN201710350234.2A CN107341284B (zh) 2017-05-18 2017-05-18 高精度预测低频电波传播特性的双向抛物方程方法

Publications (2)

Publication Number Publication Date
CN107341284A true CN107341284A (zh) 2017-11-10
CN107341284B CN107341284B (zh) 2020-11-17

Family

ID=60221116

Family Applications (1)

Application Number Title Priority Date Filing Date
CN201710350234.2A Active CN107341284B (zh) 2017-05-18 2017-05-18 高精度预测低频电波传播特性的双向抛物方程方法

Country Status (1)

Country Link
CN (1) CN107341284B (zh)

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN111027173A (zh) * 2019-10-31 2020-04-17 中铁二院工程集团有限责任公司 基于改进ssft算法的弯曲隧道电磁建模及仿真方法
CN113272812A (zh) * 2019-01-10 2021-08-17 X开发有限责任公司 用于优化电磁设备的物理特性的系统和方法
CN114741839A (zh) * 2022-03-02 2022-07-12 西北工业大学 一种分析甚低频电磁波在地-电离层中传播的fdtd方法
CN116992192A (zh) * 2023-09-28 2023-11-03 山东科技大学 基于抛物方程预测海冰混合路径中的低频电波传播方法
CN117992707A (zh) * 2024-04-01 2024-05-07 山东科技大学 基于积分方程预测复杂路径中低频电波传播特性的方法

Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5864374A (en) * 1995-04-04 1999-01-26 Mitsubishi Denki Kabushiki Kaisha Method and apparatus for image generation
CN103018719A (zh) * 2012-11-29 2013-04-03 电子科技大学 一种oth雷达发射波形的生成方法
KR20130041682A (ko) * 2011-10-17 2013-04-25 한국과학기술원 지향성이 향상된 다차원 편파 안테나
CN104142908A (zh) * 2013-05-07 2014-11-12 中国人民解放军海军航空工程学院 电波传播抛物方程分步傅里叶变换解的上边界处理方法
CN105740204A (zh) * 2016-03-14 2016-07-06 西安理工大学 不规则地形下低频段大地电导率快速反演方法

Patent Citations (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US5864374A (en) * 1995-04-04 1999-01-26 Mitsubishi Denki Kabushiki Kaisha Method and apparatus for image generation
KR20130041682A (ko) * 2011-10-17 2013-04-25 한국과학기술원 지향성이 향상된 다차원 편파 안테나
CN103018719A (zh) * 2012-11-29 2013-04-03 电子科技大学 一种oth雷达发射波形的生成方法
CN104142908A (zh) * 2013-05-07 2014-11-12 中国人民解放军海军航空工程学院 电波传播抛物方程分步傅里叶变换解的上边界处理方法
CN105740204A (zh) * 2016-03-14 2016-07-06 西安理工大学 不规则地形下低频段大地电导率快速反演方法

Non-Patent Citations (3)

* Cited by examiner, † Cited by third party
Title
DAN-DAN WANG等: "Parabolic Equation Method for Loran-C ASF Prediction Over Irregular Terrain", 《IEEE ANTENNAS AND WIRELESS PROPAGATION LETTERS》 *
GOKHAN APAYDM等: "Groundwave Propagation at Short Ranges and Accurate Source Modeling", 《IEEE ANTENNAS AND PROPAGATION MAGAZINE》 *
张志禹等: "基于宽角双向抛物线方程的高频电波传播预测方法", 《微型机与应用》 *

Cited By (8)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN113272812A (zh) * 2019-01-10 2021-08-17 X开发有限责任公司 用于优化电磁设备的物理特性的系统和方法
CN113272812B (zh) * 2019-01-10 2024-01-16 X开发有限责任公司 用于优化电磁设备的物理特性的系统和方法
CN111027173A (zh) * 2019-10-31 2020-04-17 中铁二院工程集团有限责任公司 基于改进ssft算法的弯曲隧道电磁建模及仿真方法
CN114741839A (zh) * 2022-03-02 2022-07-12 西北工业大学 一种分析甚低频电磁波在地-电离层中传播的fdtd方法
CN114741839B (zh) * 2022-03-02 2024-04-30 西北工业大学 一种分析甚低频电磁波在地-电离层中传播的fdtd方法
CN116992192A (zh) * 2023-09-28 2023-11-03 山东科技大学 基于抛物方程预测海冰混合路径中的低频电波传播方法
CN116992192B (zh) * 2023-09-28 2023-12-12 山东科技大学 基于抛物方程预测海冰混合路径中的低频电波传播方法
CN117992707A (zh) * 2024-04-01 2024-05-07 山东科技大学 基于积分方程预测复杂路径中低频电波传播特性的方法

Also Published As

Publication number Publication date
CN107341284B (zh) 2020-11-17

Similar Documents

Publication Publication Date Title
CN107341284B (zh) 高精度预测低频电波传播特性的双向抛物方程方法
CN106874549B (zh) 一种高精度预测asf的窄带离散分布抛物方程方法
CN107545104A (zh) 基于三维抛物方程的不规则地形电波传播因子预测方法
CN112949134B (zh) 基于非结构有限元方法的地-井瞬变电磁反演方法
Zelley et al. A three-dimensional parabolic equation applied to VHF/UHF propagation over irregular terrain
Azzarone et al. IONORT: A Windows software tool to calculate the HF ray tracing in the ionosphere
CN108180918B (zh) 一种点云测地路径正向跟踪生成方法及装置
Wang et al. Parabolic equation method for Loran-C ASF prediction over irregular terrain
CN108267781B (zh) 任意曲面非均匀介质快速行进程函方程求解射线追踪算法
CN105868571B (zh) 一种“m(2,4)fdtd+fdtd”的低色散低频地波传播时延预测方法
CN105740204A (zh) 不规则地形下低频段大地电导率快速反演方法
CN105573963B (zh) 一种电离层水平不均匀结构重构方法
CN104539340A (zh) 一种基于稀疏表示和协方差拟合的稳健波达角估计方法
He et al. Two-way propagation modeling of expressway with vehicles by using the 3-D ADI-PE method
CN116992192B (zh) 基于抛物方程预测海冰混合路径中的低频电波传播方法
Gherm et al. HF propagation in a wideband ionospheric fluctuating reflection channel: Physically based software simulator of the channel
CN117909626A (zh) 一种分层有耗介质中辐射源层电磁场的高效求解方法
CN117374985A (zh) 一种电力系统最优潮流问题的量子算法及包括其的装置
CN101937026B (zh) 高精度预测地波传播衰减因子的方法
CN106485071B (zh) 一种多层分组结构快速近远场转换方法
Altun et al. Electromagnetic propagation modeling over irregular terrain using a new hybrid method
CN105549031A (zh) 一种卫星信号的电离层传播时延的时域数值计算方法
CN110245316A (zh) 一种电离层参数的反演方法
Yin et al. Study on calculation and verification of radiowave propagation using parabolic equation for the antenna near the ground
CN118013172B (zh) 基于梯形公式求解积分方程的低频地波传播特性预测方法

Legal Events

Date Code Title Description
PB01 Publication
PB01 Publication
SE01 Entry into force of request for substantive examination
SE01 Entry into force of request for substantive examination
GR01 Patent grant
GR01 Patent grant