基于运动学与动力学特征探究GPR管道探测模拟

摘要：GPR管道探测数值模拟主要采用基于Maxwell方程组的时域有限差分（FDTD）法与基于几何光学原理的射线追踪法，本文对两种方法进行了阐述。并对地下管道模型GPR数值模拟与物理模拟特征进行了探讨。研究结果表明：管道周围介质导电率与对雷达波的吸收衰减成正比，同相轴幅度随着来自管道顶、底回波能量的减弱而降低。

关键词： GPR管道探测模拟 射线追踪法 探地雷达 时域有限差分模拟

PR（GPR）是利用天线向介质中发射高频电磁脉冲波来探测介质中目标体的一种物理探测技术；GPR探测的物理基础是被探测目标与周围介质之间存在电磁性（导电率、介电系数及磁导率）差异[1-2]。GPR通常采用反射剖面法进行探测，GPR剖面记录来自地下介质中的电磁反（散）射回波的走时及幅度信息，根据GPR剖面上来自目标体回波的运动学与动力学特征来推演地下目标的位置、埋深及属性等相关信息。

通常条件下，由各种不同材质构成的地下管道与其周围土壤介质之间的导电率与介电系数的差异几乎是普遍存在的；另外，各种地下管道一般都具有较平整、光滑的外表面，因而对雷达波而言是良好的散射截面。因此，GPR探测地下管道具有优良的地球物理前提条件。与传统的基于电磁感应原理的管线探测技术相比，由于GPR具有分辨率高、图像直观，而且既能探测金属管道也能探测非金属管道等显著优点[3-4]，因此，GPR目前已成为探测地下管道的主要方法之一。

1、GPR管道探测数值模拟方法介绍

目前，GPR管道探测数值模拟主要采用基于Maxwell方程组的时域有限差分（FDTD）法与闷喊基于几何光学原理的射线追踪法。

1.1 时域有限差分（FDTD）法简介

时域有限差分（FDTD）法由K.S.Yee[12]于1966年提出，因为其计算过程直观、通用性好而得到广泛应用。时域有限差分法是利用离散的二阶中心差分方程近似替代连续的Maxwell微分方程。因此，构建合理的差分格式是保证微分方程解准确性的关键步骤。同时，又由于离散化采样的过程会导致电磁波的频散现象的发生，在FDTD计算过程中，时间步长Δt与空间步长Δl的取值必须满足一定的限定条件，否则差分方程解将会不稳定。当电磁波在介质中传播时，电磁场的边界总被看作是开放的，电磁场分布空间被看作是无限的。而时域有限差分是通过计算机对电磁场的传播及空间分布进行数值模拟，由于计算机内存是有限的，所以其差分网格必然会被截断。为了使电磁波不在差分网格截断处发生明显反射，就必须设置吸收边界条件以确保计算机模拟的有限空间与自然界无限大空间之间的差异尽量的小[13]。

1.2 射线追踪法简介

射线追踪法已广泛应用于波传播的正演模拟研究。射线追踪方法计算速度快，结果直观，在层析成像、偏移、反演及手慧模型试算中都占有很重要的地位。射线追踪方法首先是在地震勘探和研究中发展很快，并取得了实际的进展。由于GPR理论和地震勘探理论的相似性，因而近年来也有学者将射线追踪方法应用于GPR波的传播和正演研究，如Cai（1995）应用射线追踪方法进行了二维介质中雷达波的传播与模拟研究[9]。射线追踪即对射线路径进行准确计算和描述，对于GPR正演模拟，射线路径是指雷达波的传播轨迹。射线追踪法的理论基础是在高频近似条件下电磁波的主能量沿射线轨迹传播。试射法是一种基于斯奈尔定律的射线追踪法，其基本思路是：射线以某一给定的初始发射角从发射点射出，在介质中波的传播遵循斯奈尔定律，记录射线的传播轨迹，直至到达指定的接收点附近满足一定精度的区域内；若射线终止位置不满足该精度，则调整初始出射角，重新计算直至射线满足精度为止。

基于射线追踪法的GPR数值模拟的核心是射线追踪，射线路径越准确，模拟出的结果才更接近于真实情况。基于试射法的GPR数值模拟实现主要有以下几个步骤[15]：

（1）建立管道模型；（2）用一组射线对模型进行试射，根据斯奈尔定律求出各条射线的终点，找出满足一定精度要求的射线的入射角蚂薯野的大致范围；（3）在区间内应用二分法进行迭代，使射线满足一定的精度要求；（4）根据射线与各个界面的交点和各个层位的物性参数，求出射线的旅行时间和各个界面的反射波幅度；（5）用已求出的反射波幅度的时间序列与脉冲子波进行褶积，即得出单道雷达记录；将各道雷达记录按顺序排列，便得出模拟雷达剖面。因此，从理论上讲，圆柱状管道GPR时—距曲线的张开度（离心率）的大小仅与管道周围介质雷达波速度有关，即波速值越大，双曲线的离心率越大，双曲线张开度越大。由于GPR管道探测总是以地面位置为参考来分析和判断来自管道回波同相轴的变化，因此就有可能把由于管道直径的变化所引起的回波同相轴张开度的变化误判为管道周围介质中雷达速度变化引起的，或者做出相反的判断。

2、地下管道模型GPR数值模拟与物理模拟特征分析

为了对地下管道探测GPR响应特征进行模拟研究，必须建立包含管道的各种地电模型，这必然涉及到各种管道的材质、管道周围介质以及管道内充填介质等。

2.1 不同材质与管内充填不同性质流体的管道GPR响应模拟分析

图1a、b、c和d为采用EKKO100型GPR系统采用主频为225ZHz天线对金属与非金属管道进行物理模拟探测的GPR记录（原始数据由加拿大Sensor&.Software公司提供）。管道埋设在电性不均匀的砂质土壤中，图1a为直径为15.5cm充满空气的金属管道；图2b、c、d分别为充满自来水（电阻率为34.5Ω·m）、空气以及咸水（电阻率为0.212Ω·m）直径为16cm的塑料管道。所有管道顶面埋深均为1.25m。图2a表明，由于金属材质的高反射率与强吸收衰减作用，GPR剖面只能看到来自管道顶面的回波同相轴，管顶回波幅度强，双曲线形同相轴对称性好。图1b显示充满自来水的塑料管道的GPR剖面上不但能看到管道顶、底的清晰回波，还能隐约看到来自管道顶、底之间的多次波同相轴。图1c显示充满空气的塑料管道由于来自管道顶、底回波时差很小（由于电磁波在空气中的速度约为水中速度的9倍，因此在图1c中的GPR记录中不但看不到来自管底的强反射回波，而且由于来自管底回波与管顶回波部分重叠而产生相消干涉而导致GPR记录上管顶回波幅度明显减弱了。图1d显示充满高导电率的咸水管道由于管内咸水的强吸收衰减作用GPR记录上只能看到来自管顶的回波同相轴。

图1采用225MHz天线0.5m天线距对金属与非金属管道进行物理模拟的记录。（a）、充满自来水的塑料管道（b）、充满空气的塑料管道（c）及充满咸水的塑料管道（d）进行物理模拟的GPR记录。图（b）中T为管顶反射，B为管底反射，M为管内顶、底之间的多次波。

2.2 周围介质属性对管道GPR响应特征的影响模拟分析

管道一般埋设在土壤等介质中，管道周围介质属性不同，GPR管道探测的响应特征也会不同。为了考察周围介质属性变化对管道GPR探测响应特征的影响，建立以下管道模型：管道材质为混凝土，管径1.2m，管道周围介质分别为淤泥、粉质黏土和自来水；管顶埋深1.0m，管壁厚度0.05m；雷达激励源子波中心频率为200MHz。采用射线追踪法进行模拟。

当管道埋设在粉质粘土中时，由于粉质黏土导电率相对较低，模拟记录上管道顶、底回波同相轴幅度强；当管道埋设在在淤泥中时，由于淤泥的导电率较高，对雷达波具有强吸收衰减，管道顶、底回波能量弱，同相轴幅度低（图8b）；当管道位于在自来水中时，由于自来水的导电率介于粉质黏土与淤泥之间，加上水的介电常数值最大，雷达波在其中的传播速度最低，因此，管顶回波同相轴幅度较淤泥介质中的明显，但比粉质黏土中的明显偏低。另外，管道周围介质雷达波速度越低，管顶回波同相轴形状越偏狭，张开度越小。

3、讨论与结论

基于单道反射的GPR剖面上圆柱状管道的回波时—距曲线为双曲线，该曲线的离心率及其渐近线的斜率与管道半径大小无关，只与管道周围介质中雷达波速度有关，渐近线斜率与管道周围介质雷达波速度成反比，即管道周围介质中雷达波速度越高，管道回波时—距曲线张开度越大。因此，管道周围介质的电磁波速度是影响管道回波同相轴形态的主要因素。管径越大，管道回波同相轴双曲线张开度越大，管顶和管底回波时差越大，回波同相轴分离越开；反之，当管道半径足够小时，管顶和管底反射同相轴趋于重合。通过对GPR记录上来自圆柱状管道回波同相轴的拟合可以确定其时—距曲线即双曲线方程，并由此可以反演管道周围介质的平均雷达波速度及管道半径等数，根据波速可以确定管顶埋深。

对不同材质与管内充填不同性质流体的管道GPR响应物理模拟与数值模拟结果表明：由于金属材质的高反射率与强吸收衰减作用，金属管道GPR剖面只有来自管顶的回波同相轴，回波幅度强，双曲线形同相轴对称性好；充满自来水的塑料管道的GPR剖面上不但有来自管道顶、底的清晰回波存在，还能隐约看到来自管道顶、底之间的多次波同相轴；充满空气的塑料管道由于来自管道顶、底回波时差很小，因此GPR记录中不但看不到来自管底的强反射回波，而且由于来自管底回波与管顶回波部分重叠而产生相消干涉而导致GPR记录上管顶回波幅度明显减弱；充满高导电率的咸水管道由于管内咸水的强吸收衰减作用GPR记录上只能看到来自管顶的回波同相轴。在潮湿黏土介质中随管顶埋深的增加，雷达波受吸收衰减及球面扩散的影响，管道回波能量逐渐减弱；反射回波同相轴双曲线张开度也随管顶埋深加大而增大。管道周围介质的电性不但影响管道回波同相轴的形态，也会明显影响其幅度；管道周围介质导电率越高，对雷达波的吸收衰减越强，来自管道顶、底回波能量约弱，同相轴幅度越低。

微波系统的设计越来越复杂,对电路的指标要求越来越高,电路的功能越来越多,电路的尺寸要求越做越小,而设计周期却越来越短。传统的设计方法已经不能满足微波电路设计的需要,使指锋灶用微波EDA 软件工具进行微唯扮波元器件与微波系统的设计已经成为微波电路设计的必然趋势。EDA即Electronic Design Automation, 电子设计自动化；目前,国外各种商业化的微波EDA 软件工具不断涌现，微波射频领域主要的EDA 工具首推Agilent 公司的ADS 软件和Ansoft 公司的HFSS、Designer 软件，其次是比较小型的有MicrowaveOffice, Ansoft Serenade, CST, Zeland, XFDTD, Sonnet 等电路设计软件。下面将会将会简要地介绍一下各个微波EDA 软件的功能特点和使用范围，以期大家有个总体的了解。微波EDA 仿真软件与电磁场的数值算法密切相关，在介绍微波EDA 软件之前先简要的介绍一下微波电磁场理论的数值算法。所有的数值算法都是建立在Maxwell方程基渗组之上的，了解Maxwell方程是学习电磁场数值算法的基础；在频域，数值算法有：有限元法( FEM -- Finite Element Method）、矩量法( MoM -- Method of Moments），差分法( FDM -- Finite Difference Methods），边界元法( BEM -- ），和传输线法( TLM -- Transmission-Line-matrix Method），在时域，数值算法有：时域有限差分法( FDTD – Finite Difference Time Domain ），和有限积分法( FIT – Finite Integration Technology ）。如果想进一步了解各种数值算法的具体实现，可以参阅以下几本书籍：① Microwave Circuit Modeling Using Electromagnetic Field Simulation,②Numerical Techniques in Electromagnetics,③ Electromagmetic Simunation Using the FDTD Method，④ Complex eletromagnetic problems and numerical Simulation Approaches。其中，使用矩量法( MoM ） 的微波EDA软件有ADS，Ansoft Designer，Microwave Office, Zeland IE3D，Ansoft Esemble，Super NEC和FEKO；使用有限元法 ( FEM ）的微波EDA软件有HFSS和ANSYS；使用时域有限差分法( FDTD ） 的微波EDA软件有 EMPIRE和XFDTD，使用有限积分法( FIT ） 的微波EDA软件有CST Microwave Studio和CST Mafia。下面来介绍较流行几种的微波EDA软件的功能和应用。ADS– Advanced Design System，是Agilent公司推出的微波电路和通信系统仿真软件，是国内各大学和研究所使用最多的软件之一。其功能非常强大，仿真手段丰富多样，可实现包括时域和频域、数字与模拟、线性与非线性、噪声等多种仿真分析手段，并可对设计结果进行成品率分析与优化，从而大大提高了复杂电路的设计效率，是非常优秀的微波电路、系统信号链路的设计工具。主要应用于：射频和微波电路的设计，通信系统的设计，DSP设计和向量仿真。现在最新的版本是ADS2005A。Ansoft Designer，是Ansoft公司推出的微波电路和通信系统仿真软件；它采用了最新的视窗技术，是第一个将高频电路系统,版图和电磁场仿真工具无缝地集成到同一个环境的设计工具，这种集成不是简单和界面集成,其关键是Ansoft Designer独有的"按需求解"的技术,它使你能够根据需要选择求解器,从而实现对设计过程的完全控制。AnsoftDesigner实现了“所见即所得”的自动化版图功能，版图与原理图自动同步，大大提高了版图设计效率。同时,Ansoft还能方便地与其他设计软件集成到一起，并可以和测试仪器连接，完成各种设计任务，如频率合成器，锁相环，通信系统，雷达系统以及放大器，混频器，滤波器，移相器，功率分配器，合成器和微带天线等。主要应用于：射频和微波电路的设计，通信系统的设计，电路板和模块设计，部件设计。

电磁学中的算法有很多种，分为时域和频域两类。常见的有时域有限差分法（FDTD）、时域有限积分法（FITD）、矩量法（MoM）、有限元法（FE）、边界元掘胡法（BEM）、 谱域法（SM）、传输线法（TLM）、模式匹配法（MM）、横向谐振法（TRM）、线方法（ML）和野帆解析法等等。频域数值算法有：有限元法 (FEM -- Finite Element Method）、矩量法(MoM -- Method of Moments），差分法(FDM -- Finite Difference Methods），边界元法( BEM --Boundary Element Method），和传输线法(TLM -- Transmission-Line-matrix Method）

最近，又多了一种EIT算法，这是一种基于FDTD算法的革新技术。EIT算法克服了传统FDTD算法在模拟弯曲金属界面和介质界面时的梯形误差，避免精度损失，保持算法精度和效率；解决了共形FDTD算法稳定性要求导致时间步长降低的效率问题，无需减少时间步长，保持计算速度。简单来说，在确保同等精度的情况下，EIT的计算速度更快，且占用内存小。

代表性软件：在传统的国外CAE软件中，CST软件（德国的）用的是FIT有限积分技术颂散雹，HFSS软件（ANSYS公司推出的三维电磁仿真软件）用的是FEM有限元技术。目前，EIT算法主要代表性CAE软件是中望电磁（中望软件推出自主知识产权的全波三维电磁仿真软件），适合求解天线/天线阵列，雷达，微波器件，电磁兼容/电磁干扰，高速互联SI，电磁传播和散射等，以及任意结构电大宽带的电磁问题。

---