[发明专利]一种基于FETD与FDTD耦合的二维模型的探地雷达正演方法

说明书

本发明公开了一种基于FETD与FDTD耦合的二维模型的探地雷达正演方法，包括以下步骤：加载二维模型参数；将探地雷达数值模拟区域分为模型区域和扩展区域；采用非结构化三角形单元对模型区域进行剖分，采用规则网格对扩展区域进行剖分；设置测量点，在与测量点对应的位置加载激励源；设置初始时间步之前所有时间的电场与辅助场分量均为零，根据FDTD辅助场时间步进公式更新辅助场，根据FETD主场时间步进公式更新电场，直到最大时间步为止，得到测量点处各个时间步的电场值；根据各个测量点处所有时间步的电场值得到探地雷达正演剖面图。本发明拟合精度高、计算快速。

技术领域

本发明涉及一种地球物理领域的探地雷达数值模拟方法，特别适合几何特征不规则、物性参数分布复杂模型的快速高效正演。

背景技术

探地雷达数值模拟方法很多，包括有限差分法(Finite Difference TimeDomain，FDTD)^[1-3]、有限单元法^[4-5]、时域伪谱法^[6]、辛分块龙格库塔方法^[7]、时域间断有限元(DGM)^[8]、样条小波有限元等算法^[9]，这些算法各有特色。其中时域有限差分法具有直接时域计算、编程容易、节约存储空间和计算时间等优点而备受青睐(IrvingKnight，2006^[10]；CassidyMillington，2009^[11]；冯德山等，2010^[12])，但FDTD受CFL稳定性条件的限制、不能与非结构化网格结合，对不规则脱空、任意形体的开裂、渗水通道等复杂隧道衬砌病害拟合效果不好。而非结构化网格的时域有限元算法(FETD)能够对典型隧道病害实现贴体剖分，具有较高拟合精度。而在模拟区域边界采用CPML边界条件，它能有效压制来自截断边界处的反射波，同时，对传输波、隐失波、低频波等综合吸收效果较好，能极大提升FETD边界的处理效果。但CPML边界的最初提出是针对FDTD算法(Berenger，1996^[13])，公式体系并不适用于FETD。尽管有学者(MasoudAbdolali,2007^[14]，2009^[15]；DondericiTeixeira，2008^[16]；CorreiaJin，2006^[17])将CFS-PML边界条件应用于FETD的模拟中；但计算公式过于复杂，不便于编程实现。

因此，有必要设计一种拟合精度高、计算快速的二维模型的探地雷达正演方法。

发明内容

本发明所解决的技术问题是，针对现有技术的不足，提供一种基于FETD与FDTD耦合的二维模型的探地雷达正演方法，主场更新时采用时域有限元法求解，辅助场更新时采用时域有限差分法显式计算的耦合方法，充分利用了FETD算法剖分精度高，而FDTD显式递推计算辅助场计算快速的特点。

本发明所提供的技术方案为：

一种基于FETD与FDTD耦合的二维模型的探地雷达正演方法，包括以下步骤：

步骤1、加载二维模型参数；

步骤2、将探地雷达数值模拟区域分为模型区域和扩展区域；模型区域即二维模型所处区域；扩展区域是根据边界物性参数对模型区域进行扩展得到，包括过渡区域和PML区域，其中过渡区域为设置于PML区域和模型区域之间的一个有限宽度的区域；采用非结构化三角形单元对模型区域进行剖分，采用规则网格对扩展区域进行剖分，得到多个离散网格节点；

步骤3、自左向右设置M个测量点，即M个探地雷达接收天线位置，以及与M个探地雷达接收天线位置一一对应的M个探地雷达发射天线位置；初始化测量点序号m＝1；M为大于或等于1的整数；M个测量点分别设置于M个离散网格节点处；

步骤4、在第m个探地雷达发射天线位置加载激励源，计算与其对应的第m个测量点处各个时间步的电场值；