基于探地雷达的水库坝基渗漏正演模拟

吴学礼，贾江波，孟凡华，李素康，张浩浩

(1.河北科技大学电气工程学院，河北石家庄 050018；2，河北省生产过程自动化工程技术研究中心，河北石家庄 050018，3.石家庄市佳信电气有限公司，河北石家庄 050018)

基于探地雷达的水库坝基渗漏正演模拟

吴学礼1,2，贾江波1,2，孟凡华1,2，李素康3，张浩浩1,2

(1.河北科技大学电气工程学院，河北石家庄 050018；2，河北省生产过程自动化工程技术研究中心，河北石家庄 050018，3.石家庄市佳信电气有限公司，河北石家庄 050018)

水库坝基渗水将带来溃堤的风险，造成水资源浪费及生命财产的损失，针对现有技术无法确定坝体渗水洞坐标的实际情况，提出采用探地雷达对坝体进行日常检测。利用时域有限差分法(FDTD)对在二维空间进行FDTD按照Yee氏网格离散处理的原理，得到其电磁分布方程，基于水库坝基的实际构造，建立了空气层、混凝土层、黏土层及两个渗水洞的理想模型，并依据理想模型建立了带有石灰岩干扰点的混凝土层干扰模型。利用MATLAB实现系统架构，并进行正演模拟。结果显示：探地雷达可以用于较深目标体探测，通过改变探地雷达电磁波中心频率，并对比50，100和200 MHz三种频率雷达电磁波的检测频谱，又对带干扰模型进行了模拟，验证了100 MHz频率雷达电磁波，用于复杂地形的检测效果更优。

空间物理探测；探地雷达；水库大坝；Yee网格；正演模拟

水库大坝是重要的水利工程，在国家生产和安全中起着关键作用。从整体看，大坝主要包括混凝土坝和土石坝两大类，据资料显示中国95%的水库大坝为土石坝，由大坝所在地附近的土、石块等经过加工建成[1]。20世纪50—70年代，是中国水坝的主要修建时期，由于当时的现实情况，以现在的要求来审视普遍存在着一些问题[2]，坝基填土不实导致坝基以下原始层出现渗漏点，水库中的水通过坝底经过渗漏点进行渗漏，对坝体结构安全和防渗安全带来很大的威胁，坝基一般位于坝体下，无法直接通过肉眼检测，一旦发生自然灾害，极易造成溃堤，带来极大危害[3]。根据中华人民共和国水利部批准的《堤防隐患探测规程》(SL436—2008)所提出的检测方法，探地雷达(ground penetrating radar，GPR)的防患探测方法在实际操作中有很大的可行性和便利性[4]。GPR是根据高频脉冲电磁波在地下不同物质媒介间的反射及绕射等物理波动规律,分析接收到的波形、振幅强度和时间等特征推断地下情况来进行工作的无损检测仪器。

SZERBIAK等[5]研究了探地雷达在模拟三维水库中流体渗透率模型的问题。任爱武等[6]研究了高泉水库的水文地质条件和水库渗透原因，表明探地雷达在岩溶水库地区渗漏原因分析方面具有一定的适用性。郭秀军等[7]研究了探地雷达在检测南水北调东线渠道两种砌混凝土，给出了施工中材料选择的建议，指出GPR是实现大型渠道混凝土砌混凝土厚度快速检测的有效方法。周立刚等[8]研究了高频探地雷达探测土壤表层含水量的测定，探测表面频率在1～2 GHz时对土壤表层含水量的识别度更高。何开胜等[9]以甘肃省民乐山丹地震时受到一定破坏的李桥水库为对象,用探地雷达对水库大坝填土不密实引起的渗漏病害进行了探测研究，探测结果发现了上下游的渗漏通道并对整个水库地形做出了分析。上述的研究表明，探地雷达在探测水库坝基渗漏方面可以发挥很大的作用。

1 二维空间有限时域差分方程法

图1 二维Yee氏网格图及电磁场分布Fig.1 Two-dimensional Yee's mesh and electromagnetic field distribution

探地雷达的根本原理是研究高频电磁波电场和磁场在空间中的变换，符合Maxwell方程，同时满足电磁波传播的本构关系，在二维平面正演模拟中将每一道电磁波产生的回波称为A扫描，A扫描组合形成的二维平面称为B扫描。在二维空间中采用横电磁波(TM型电磁波)进行模拟，按照Yee提出的Yee氏网格描述[10]对空间进行离散化处理[11-13]。

图1为二维Yee氏网格模型，设二维Yee氏网格的大小分别为Δx和Δy，Δt为时间间隔，观察点(x，y)为Ex，Hx，Hy的节点，表1中标明了节点坐标及对应时刻。

表1 节点坐标及对应时刻

得到TM波FDTD方程：

(1)

(2)

(3)

式(2)、式(3)中电磁参数中标m的取值与其公式中右端电磁或磁场分量的空间位置相同，即：

(4)

2 基于MATGPR的水库坝基二维模型

MATGPR是基于MATLAB的一款跨平台、可定义、可扩展的雷达分析程序包[14]，由雅典大学开发，采用时域有限差分法(FDTD)，由两层软件系统构成，底层包括数据处理、显示的脚本和程序，顶层包括图形用户界面，有数据基本处理、滤波处理和仿真等功能[15]。通过菜单选择输入数据，输入文件之后将输入数据显示出来，进行选择数据处理方式，并显示处理结果。MATGPR3.1可以实现二维及三维雷达实际扫描图的显示，同时可以实现二维模型的仿真得到B扫描以及每一道对应的A扫描波形图[16]，有利于进行正演模拟。在模拟过程中要考虑探地雷达的主要技术参数分辨率和探测深度，并考虑雷达的中心频率、时窗的选择、天线性能的差异、设计角度的不同、参数选择不同等等影响因素[17]。

假设存在一个长10 m，深5 m的混凝土层，混凝土层有一条不规则的空气裂缝，依据常见介质的介电性质特征[18]列出混凝土层和空气的特性如表2所示。

表2 混凝土层和空气裂缝参数选择

按照表2提供的参数，建立了相对应的空气裂缝模型，如图2所示。

对该模型取中心频率f=400MHz的电磁波进行扫描，得到如图3所示的波形图。

图2 空气裂缝模型图Fig.2 Air crack model

图3 空气裂缝400 MHz雷达扫描结果Fig.3 Air crack 400 MHz radar scan results

从以上对空气裂缝的仿真及结果来看，MATGPR可以很好地对模型进行正演模拟得到有一定清晰度的扫描图。下面将介绍基于MATGPR的坝基建模及雷达不同中心频率对应的仿真结果。

3 水库坝基二维模型的构建及扫描分析

在实际情况中，水库坝基以下的地质条件较为复杂，不仅有混凝土，还会存在不同含水量的土壤、石块、动植物遗体、管道等等情况，对模型构建有很大的考验[19]，现假设理想状况：混凝土坝基下湿土壤中有渗透点。为了模型构建简单，假设湿土壤为单一物理状态下的同一物质，其中相对磁导率μ取1。

模型的构成如下：

h1=0.447 92 m的空气(天线与地面间的距离)；

h2=2.020 78 m的湿混凝土(坝基的混凝土层)；

h3=7.531 30 m的土壤层(坝基下的原始土壤层)[20]。

假设有两个渗透点，洞1圆心(4.962 5,4.989 6)半径R1=1.600 0 m，电阻率为ρ1=80 Ω·m，相对介电常数取14，传播速度为γ1=0.080 083 m/ns；洞2圆心(13.962 5,6.239 6)半径R2=1.825 0 m,电阻率为ρ2=30 Ω·m，相对介电常数取25，传播速度为γ2=0.059 891 m/ns，得到坝基双渗漏点二维模型见图4。

对图4所示模型采取中心频率100 MHz的探地雷达B扫描结果如图5所示，可以明显分辨出渗水洞的位置，也可以清楚看出对应的空气层和混凝土层、混凝土层和土壤层的分界线。

图4 坝基双渗漏点二维模型Fig.4 Two-dimensional model of dam foundation double leakage point

图5 中心频率100 MHz B扫描结果Fig.5 Center frequency 100 MHz B-scan results

接着对得到的B扫描进行洞1和洞2的顶端A扫描分析，结果如图6所示，可以看出电磁波经过长距离的扩散，在空气层与混凝土层的界面反射最为明显，经过混凝土层的衰减之后，虽然可以分辨出混凝土和土壤的分界线以及洞顶的位置，但是回波波形衰减严重，振幅较小。

图6 100 MHz雷达波A扫描分析Fig.6 100 MHz radar A-scan analysis

对图4所示模型进行50MHz和200MHz的雷达波扫描，结果分别如图7、图8所示，同100MHz的扫描结果对比可以看出来，50MHz显示结果一般，分辨率较低，整体A扫描道数也较少，但是相对应反射波的振幅衰减较小；200MHz的扫描结果分辨率较高，但是由于衰减严重，探测深度不够，因此对于该模型中目标体可见度较低。

图7 50 MHz雷达波A扫描分析Fig.7 50 MHz radar A-scan analysis

图8 200 MHz雷达波A扫描分析Fig.8 200 MHz radar A-scan analysis

从A扫描和B扫描结果分析来看，200MHz在成像效果上较为细腻，但是从A扫描可看出波形衰减较大，不适用于5～6m深度探测。50MHz的A扫描结果波形衰减较小，可有效识别目标体，但是扫描道数少，导致成像效果一般。综合对比来看，100MHz在有效识别目标体的同时还可以兼具良好的成像效果。结合实际情况，在图4所示的二维模型中混凝土层随机加入半径5～15cm的石灰岩模型10个，电阻率取ρ3=2Ω·m，相对介电常数取4，传播速度为γ3=0.104 373 8m/ns，其他参数不变。新建立的带干扰坝基双渗漏点二维模型如图9所示，对该模型采取中心频率100MHz的探地雷达B扫描，结果如图10所示。

图9 带干扰坝基双渗漏点二维模型Fig.9 Two-dimensional model of dam foundation double leakage point with jamming

图10 100 MHz雷达波B扫描分析Fig.10 100 MHz radar B-scan analysis

从图10的B扫描结果来看，虽然在二维模型的混凝土层加入了随机的石灰岩模型，依然可以清楚地分辨出混凝土层和土壤层的分割线，同时也可以看到半径较大的石灰岩会有很强烈的回波反射，半径较小的则难以发现，说明了100MHz无法识别这一大小的目标体，同时也说明了MATGPR可以运用于较复杂的雷达正演模拟实验中。

4 结 语

本文论证了探地雷达运用于水库坝基渗漏检测的可行性，通过对二维空间进行FDTDYee氏网格离散处理，建立了MATGPR水库坝基二维模型，进行正演模拟及分析。通过正演模拟可以看出，对在5m～6m的深度下两个渗漏点的检测，分别选取50，100和200MHz的天线中心频率进行扫描，使用A扫描和B扫描的结果都得出了100MHz比50MHz有更好的成像效果，比200MHz能更好地识别目标的结论，能够满足对该深度目标体的探测。同时针对带干扰模型的扫描结果，在识别一部分干扰物的同时也很好地反应了地表以下地层结构和渗漏点位置，对实际操作有非常好的指导作用。用MATGPR可以很好地仿真出较为复杂的地形并能够得到合理的仿真对比图。

/References:

[1] 孙继昌. 中国的水库大坝安全管理[J]. 中国水利, 2008(20):10-14. SUN Jichang. Dam safety management in China[J]. China Water Resources, 2008(20):10-14.

[2] 李宗坤, 葛巍, 王娟,等. 中国水库大坝风险标准与应用研究[J]. 水利学报, 2015, 46(5):567-573. LI Zongkun, GE Wei, WANG Juan, et al. Risk criteria and application on reservoir damsin China[J]. Journal of Hydraulic Engineering, 2015, 46(5):567-573.

[3] 苏怀智, 吴中如, 戴会超. 初探大坝安全智能融合监控体系[J]. 水力发电学报, 2005, 24(1):122-126. SU Huaizhi,WU Zhongru,DAI Huichao. System of the ories and methods monitoring intelligently dam safety[J]. Journal of Hydroelectric Engineering， 2005, 24(1):122-126.

[4] 冷元宝, 任建东, 王锐,等. 堤防隐患探测与监测技术展望[J]. 工程地球物理学报, 2004, 1(1):74-77. LENG Yuanbao,REN Jiandong,WANG Rui,et al. Prospect for the technology of detecting and monitoring embankments’ hiddendanger[J]. Chinese Journalof Engineering Geophysics, 2004, 1(1):74-77.

[5] SZERBIAK R B, MCMECHAN G A, CORBEANU R, et al. 3-D characterization of a clastic reservoir analog: From 3-DGPR data to a 3-D fluid permeability model[J].Geophysics, 2001, 66(4):1026-1037.

[6] 任爱武, 柯柏荣, 程建设,等. 岩溶地区水库渗漏原因分析与无损检测验证[J]. 水利学报, 2014(sup2):119-124. REN Aiwu, KE Bairong, CHENG Jianshe,et al. Analysis of the causes for the leakage of reservoir in karst area and nondestructive testing verification[J]. Journal of Hydraulic Engineering, 2014(sup2):119-124.

[7] 郭秀军, 张盟, 罗辉,等. 南水北调东线渠道衬砌混凝土厚度GPR快速无损检测[J]. 南水北调与水利科技, 2009, 7(6):282-285. GUO Xiujun, ZHANG Meng, LUO Hui，et al, No destructive testing by GPR method on thickness of concrete channel lining in the east route of the South-to-North WaterTransfer Project[J]. South-to-North Water Transfers and Water Science & Technology, 2009, 7(6):282-285.

[8] 周立刚, 于东升, 王玺洋,等. 基于高频探地雷达的土壤表层含水量测定[J]. 土壤学报, 2016, 53(3):621-626. ZHOU Ligang, YU Dongsheng, WANG Xiyang,et al. Determination of top soil water content based on high-frequency ground penetra-ting radar[J].Acta Pedologiga Sinica, 2016, 53(3):621-626.

[9] 何开胜, 王国群. 水库堤坝渗漏的探地雷达探测研究[J]. 防灾减灾工程学报, 2005, 25(1):20-24. HE Kaisheng,WANG Guoqun. Study of seepage through reservoir dam using ground penetrating radar detection[J]. Journal of Disaster Prevention and Mitigation Engineering, 2005, 25(1):20-24.

[10]YEE K. Numerical solution of initial boundary value problems involving Maxwell’s equations in isotropic media[J]. IEEE Transactions on Antennas & Propagation, 1966, 14(3):302-307.

[11]吴奇. 基于混合ADI-FDTD亚网格技术的GPR高效正演暨水利工程中的应用研究[D]. 长沙：中南大学, 2014. WU Qi. Ahybrid ADI-FDTD Sub-gridding Scheme for Efficient GPR Simulation and Study Application in Hydraulic Engineering[D]. Changsha: Central South University, 2014.

[12]KUNZ K S, SIMPSON L. A technique for increasing the resolution of finite-difference solutions of the Maxwell equation[J]. IEEE Transactions on Electromagnetic Compatibility, 1981, EMC-23(4):419-422.

[13]戴前伟, 冯德山, 王启龙,等. 时域有限差分法在地质雷达二维正演模拟中的应用[J]. 地球物理学进展, 2004, 19(4):898-902. DAI Qianwei, FENG Deshan, WANG Qilong, et al.The apply of finite difference time domain method in the Ground Penetrating Radar (GPR) two-dimension forward simulate[J]. Progress in Geophysics, 2004, 19(4):898-902.

[14]TZANIS A. MATGPR: A freeware MATLAB package for the analysis of common-offset GPR data[C]//Geophysical Research Abstracts. Munich: Copernicus Publications,2006, 8(09448).

[15]周奇才, 李炳杰, 何自强,等. 探地雷达数据处理软件接口编制及应用[J]. 无损检测, 2009, 31(4):300-302. ZHOU Qicai, LI Bingjie, HE Ziqiang, et al. A ground penetrating radar data processing software development and application[J].Nondestructive Testing, 2009, 31(4):300-302.

[16]TORRIONE P A, THROCKMORTON C S, COLLINS L M. Performance of anadaptive feature-based processor for awideband ground penetrating radar system[J]. IEEE Transactions on Aerospace and Electronic Systems, 2006, 42(2): 644-658.

[17]DAVIS J L, ANNAN A P. Ground-penetrating radar for high-resolution mapping of soil and rock stratigraphy[J]. Geophysical Prospecting, 1989, 37(5): 531-551.

[18]GUO P, SHI J, GAO B, et al. Evaluation of errors induced by soil dielectric models for soil moisture retrieval at L-band[C]// Geoscience and Remote Sensing Symposium (IGARSS). Washington DC:IEEE, 2016: 1679-1682.

[19]DAFFLON B, LEGER E, SOOM F, et al. Quantification of Arctic soil and permafrost properties using ground penetrating radar[C]//Ground Penetrating Radar (GPR), 2016 16th International Conference.Washington DC:IEEE, 2016: 1-6.

[20]SHAABAN F, HABEEBULLAH T M, MORSY E A, et al. Ground penetrating radar and 2D electric resistivity studies for tracing hydrocarbon leakage site, close to Abha City: A case study[J]. Arabian Journal of Geosciences, 2016, 9(20): 754.

Forward modeling of seepage of reservoir dam based onground penetrating radar

WU Xueli1,2, JIA Jiangbo1,2, MENG Fanhua1,2, LI Sukang3, ZHANG Haohao1,2

(1.School of Electrical Engineering, Hebei University of Science and Technology, Shijiazhuang, Hebei 050018, China; 2.Hebei Provincial Research Center for Technologies in Process Engineering Automation, Shijiazhuang, Hebei 050018, China; 3.Shijiazhuang Jiaxin Electric Company Limited, Shijiazhuang, Hebei 050018, China)

The risk of the reservoir dam seepage will bring the waste of water resources and the loss of life and property. The ground penetrating radar (GPR) is designed as a daily inspection system of dams to improve the existing technology which can't determine the actual situation of the dam seepage tunnel coordinates. The finite difference time domain (FDTD) is used to solve the Yee's grids discreatization in two-dimensional space, and its electromagnetic distribution equation is obtained as well. Based on the actual structure of reservoir dam foundation, the ideal model of air layer, concrete layer, clay layer and two water seepage holes is described in detail, and the concrete layer interference model with limestone interference point is established. The system architecture is implemented by using MATLAB, and the forward modeling is performed. The results indicate that ground penetrating radar can be used for deep target detection. Through comparing the detection spectrum of three kinds of frequency electromagnetic wave by changing the center frequency of the GPR electromagnetic wave of 50 MHz, 100 MHz and 200 MHz, it is concluded that the scanning result is more accurate at 100 MHz. At the same time, the simulation results of the interference model show that this method can be used for the detection of complex terrain.

space physics detection; GPR; reservoir ram; Yee grid space; forward modeling

2017-04-20；

2017-06-08；责任编辑：李 穆

石家庄市科技攻关项目(141131371A)

吴学礼(1961—)，男(满族)，黑龙江齐齐哈尔人，教授，博士，主要从事控制科学与工程方面的研究。

E-mail：wuxueli@hebust.edu.cn

1008-1542(2017)04-0389-06

10.7535/hbkd.2017yx04011

P225.1

A

吴学礼，贾江波，孟凡华，等.基于探地雷达的水库坝基渗漏正演模拟[J].河北科技大学学报,2017,38(4):389-394. WU Xueli, JIA Jiangbo, MENG Fanhua,et al.Forward modeling of seepage of reservoir dam based on ground penetrating radar[J].Journal of Hebei University of Science and Technology,2017,38(4):389-394.