Application of portable frequency domain three-dimensional electromagnetic inversion in urban geological survey

JianXin LIU; HaiJian HUA; ZhanHui QING; RongWen GUO; Rong LIU

doi:10.6038/pg2025HH0562

Progress in Geophysics >

2025 , Vol. 40 >Issue 1: 328 - 336

DOI: https://doi.org/10.6038/pg2025HH0562

Application of portable frequency domain three-dimensional electromagnetic inversion in urban geological survey

JianXin LIU ^,¹^,² ,
HaiJian HUA ¹^,² ,
ZhanHui QING ³ ,
RongWen GUO ¹^,² ,
Rong LIU ^,¹^,²^,*

Expand

¹ School of Geosciences and Info-Physics, Central South University, Changsha 410083, China
² Hunan Key Laboratory of Non-Ferrous Resources and Geological Hazard Detection, Changsha 410083, China
³ Guangdong Geological Environment Monitoring Station Guangzhou, Guangdong 510510, China

Received date: 2024-04-13

Online published: 2025-03-13

Copyright

Fold

Abstract

The portable frequency domain electromagnetic detection instrument has been widely used in the field of urban near-surface exploration because of its high sensitivity, strong anti-interference ability, simple operation and high measurement efficiency. Due to the smaller scale of urban geological exploration targets, three-dimensional inversion is more advantageous than one-dimensional or two-dimensional inversion in terms of data interpretation accuracy. In order to study the effectiveness of the portable frequency domain electromagnetic detection instrument in urban shallow surface exploration and the superiority of 3D inversion, this paper conducts 3D inversion of GEM-2 measured data from a test site in Shanghai, and compares the results with the traditional 1D inversion results. The excavation results show the correctness of the 3D inversion algorithm proposed in this paper and the feasibility and superiority of its application in urban underground space detection.

Key words： Urban geological survey; Near-surface exploration; Controllable source electromagnetic method; Three-dimensional inversion

Cite this article

JianXin LIU , HaiJian HUA , ZhanHui QING , RongWen GUO , Rong LIU . Application of portable frequency domain three-dimensional electromagnetic inversion in urban geological survey[J]. Progress in Geophysics, 2025 , 40(1) : 328 -336 . DOI: 10.6038/pg2025HH0562

0 引言

随着科技发展和社会的不断进步，城市化进程快速向前推进，针对城市地下空间的探测活动也愈加频繁，目前城市地下空间的可利用深度，依据城市规模的不同，一般可达20 m至30 m左右(何继善等，2022).合理有效的利用地下空间是城市发展的重要一环，一般而言，城市地下空间探测的主要目的是探明地层结构，地质构造以及包括巨石，岩溶，管线等在内的异常体的大小与位置等，进而为后期的相关工程实施提供可靠的先验资料(柳建新等，2023).由于城市内探测存在建筑物密集，噪声干扰大等不利因素，使得探测难以获取较为精确的地下结构特征和信息，这对地球物理方法技术提出了较高的要求.

城市地下空间探测可归属于浅地表探测，传统的方法主要是通过钻孔或盲挖来获取地下信息，该方法获得信息方式直观，但存在较大的局限性，成本高昂、勘探效率低，对地表的破坏程度高.地球物理方法发展至今，探地雷达、高密度电法以及地震波法均可以实现浅地表的信息采集.探地雷达法采用发射天线或线圈向地下目标体发射高频电磁波，通过反射回波来分析地下介质分布，其在城市地下设施探测(陈思静等，2021)、水利工程(汪魁峰，2014)、考古(李鑫，2020)等领域均有应用.对比传统钻探，探地雷达法成本较低，勘探效率高，但其信号衰减很快，有效探测深度较浅，无法满足深部信息的采集作业(李广才等，2023)；高密度电法通过在地表布设电极并由主机实现对电极的控制与切换，避免了电极的反复布设工作(何继善等，2022).目前，高密度电法在地下水勘探(郑智杰等，2019)、矿产勘查(马凤山等，2016)、滑坡监测(甘建军等，2019)等领域具有广泛的应用.高密度电法勘探效率高，采集速度快，但也存在抗干扰能力差的缺陷；地震波法利用地震波在不同介质中传播速度的不同来探测地下介质的构造，具有探测范围大，勘探精度高等优势，但地震波法探测速度较慢，成本较高，且其产生的冲击会对正常生活造成干扰，因此该方法在城市地质探测中存在较大的局限性.

便携式频率域电磁探测仪器因其具备灵敏度高，抗干扰能力强，操作简单，测量效率高等优点，在浅地表的勘探中受到了广泛的应用(何晓萍，2018)，其发展至今已有较多成熟的探测设备，例如美国GEOPHEX公司生产的GEM系列仪器，以及加拿大Geonic公司生产的EM系列仪器等.以GEM系列仪器为例，其在同一个探测时间段内，可以同时工作于多个频率状态，从而获得不同深度地层的地质信息.

一般而言，便携式频率域电磁探测仪器的采集数据往往包含多组发射与接收系统的采集信息，由于其工作时具备的多测站特点，导致其计算量大，反演效率较低，而城市地下异常体通常具有种类多样，分布不均匀，形状不规则等特点(陈健，2023)，一般的一维或二维反演很难满足实际施工所需要的准确位置等要求，对比之下，三维反演在探测空间形态上较一维或二维反演更具备优势(殷长春，1994).因此，应用三维反演对城市地下异常体进行成像具有重要的意义.三维反演的前提是三维正演，正演的快慢与否很大程度上影响着反演的效率.截断边界有限元方法通过减小计算域的面积，降低了计算成本，因而可以显著提高正演模拟的效率(柳建新等，2014).在此基础上，基于足迹计算域的截断边界有限元方法(Liu et al., 2022b)相比前者拥有更小的计算域，且其初始测站计算所得的格林函数可以运用于其他各组测站的计算中去，因而大大加快了正演计算，提高了三维反演的效率.截止到目前，较为常见的三维反演优化算法有共轭梯度法(He and Bo, 2018；Qin et al., 2016；Zhang et al., 2017)、非线性共轭梯度法(Huo et al., 2019；Hu et al., 2022；Wang et al., 2017)、拟牛顿法(Dai et al., 2017；Loke and Dahlin, 2002；Gance et al., 2012)、高斯牛顿法(Feng et al., 2020；Xu et al., 2023)以及Occam算法(Weng，2010；Patel，2020).而随着计算机计算能力和内存的快速发展，收敛速度更快的高斯牛顿法得到了广泛的应用.

本文首先介绍了便携式频率域电磁探测仪的基本构成和工作原理，其次阐述了基于足迹计算域的截断边界矢量有限元方法及高斯牛顿法的基本理论，并通过建立理论模型进行验证分析，同时依据上海某测试场地的GEM-2实测数据进行三维反演，并与一维反演结果进行了对比，经后期开挖取证，验证结果表明本文所提出的三维反演算法具有正确性及有效性，便携式频率域电磁法三维反演应用于城市地下空间探测中是可行且优越的.

1 便携式频率域电磁探测仪

便携式频率域电磁探测仪一般通过发射装置向目标区域发射一定频率的电磁波形成一次场，激发地下具有良好导电性的异常体产生二次场，接收装置在二次场的激励下产生电信号，通过分析该电信号可以获取不同探测深度的地层信息.以美国GEOPHEX公司研发的手持宽频带磁偶源电磁探测仪GEM-2为例，仪器采用小型线圈作为发射与接收装置，包含一个发射线圈和一个接收线圈，相距约1.66 m，另包含第三个补偿线圈，用于从接收线圈中移除一次场的影响，从而保证接收线圈可以接收由地质体引起的纯二次场信号.所有线圈都以固定的几何形状布置在一块长板中，实际工作时由测量人员手持仪器并保持在离地面1 m左右的高度.该仪器的工作频率范围为300 Hz~48 kHz，测量时每个频率均可以测得同相(实分量)和正交(虚分量)两个分量.

图 1展示了便携式频率域电磁探测仪的工作原理图，利用电磁感应原理，通过一个发射线圈向地下发射交变电磁场，即一次场.一次场作用于地下异常体时，由于其与周围介质电性和磁性的差异，会在异常目标体内部激发出相应的感应涡旋电流，而产生的涡流又会激发产生二次电磁场.接收线圈用于测量返回的包含地下介质电性和磁性信息的归一化二次场.采集得到的数据一般表示为二次场与一次场的比值，由于该数值较小，通常会将结果扩大倍，实际输出结果为无量纲单位10^-6(Won et al., 1996).

显示原图|下载原图ZIP|生成PPT

图1 频率域电磁探测仪的工作原理图

Fig 1 Diagram of the working principle of frequency domain electromagnetic detector

2 方法理论

2.1 基于足迹计算域的截断边界矢量有限元方法

在截断边界矢量有限元法中，将计算区域定义为异常体及该异常体外的一层单元厚度包裹层(Liu R et al., 2018)，如图 2a所示的黄色区域，使用矢量有限元规则网格将计算区域进行剖分，计算区域的边界二次电场可以由内部单元散射电流源与格林函数表示：

(1)

$E_{\mathrm{B}}^{\mathrm{S}}=g^{\mathrm{ee}} J_{\mathrm{I}},$

显示原图|下载原图ZIP|生成PPT

图2 (a) 截断边界矢量有限元以及基于足迹计算域的截断边界矢量有限元的计算区域，蓝线、黄线以及红线分别表示矢量有限元；(b)截断边界矢量有限元的计算方法

Fig 2 (a) The compact finite element method, and the footprint-guided compact finite element method. The blue, yellow, and red lines represent the calculation areas of the vector finite element method; (b) Calculation method for the compact finite element method

式中，g^ee表示格林函数，内部单元的散射电流J_Ⅰ可以通过内部单元棱边中心点电场插值得到：

(2)

$J_{\mathrm{I}}=V\left(\sigma-\sigma_*\right) N^{\mathrm{e}}\left(E_{\mathrm{I}}^{\mathrm{S}}+E_{\mathrm{I}}^{\mathrm{p}}\right),$

其中，V表示内部单元的体积，σ和σ_*分别表示异常体和围岩的电导率，N^e表示插值基函数，E_Ⅰ^S和E_Ⅰ^p分别表示内部单元的二次场和一次场.

异常体内部的二次电场方程(Liu J X et al., 2018)可由式(3)表示：

(3)

$\boldsymbol{K}_{\mathrm{II}} E_{\mathrm{I}}^{\mathrm{S}}+\boldsymbol{K}_{\mathrm{IB}} E_{\mathrm{B}}^{\mathrm{S}}=\boldsymbol{S}_{\mathrm{I}},$

其中，K_Ⅱ为对称稀疏矩阵，K_IB为非对称稀疏矩阵，S_Ⅰ为源项.

联立式(1)和式(2)即可将内部单元棱边中心点的电场与计算区域边界的电场联系起来，从而满足边界条件.由于格林函数的计算量较大，截断边界矢量有限元法的计算速度受到较大的限制.

基于足迹计算域的截断边界矢量有限元方法定义的计算区域如图 2a中的红色区域，对比截断边界矢量有限元法，拥有更小的计算域.因为足迹外部电场几乎为零，忽略足迹外部区域的异常体对包裹层外边界的影响.如图 2b所示，电场分量定义在相应方向的棱边的中心点位置(黑色箭头表示足迹边界单元的电场分量，蓝色箭头表示足迹内的电场分量).异常体单元内部中心存在散射电流源(图中的红点所示)，计算区域边界电场可以由内部散射电流源与格林函数表示.内部单元内散射电流源由内部单元棱边中心点电场插值表示.由于第一组测站计算得到的格林函数可以运用于其他各组测站的计算中去，因而显著降低了计算时间，提升了正演效率.我们将计算区域的边界二次电场E_B^S表述为：

(4)

$E_{\mathrm{B}}^{\mathrm{S}}=G^{\mathrm{ee}}\left(E_{\mathrm{I}}^{\mathrm{S}}+E_{\mathrm{I}}^{\mathrm{p}}\right),$

其中，G^ee表示：

(5)

$G^{\mathrm{ee}}=g^{\mathrm{ee}} V\left(\sigma-\sigma_*\right) N^{\mathrm{e}} .$

将式(4)代入式(3)即可得到异常体内部的二次电场方程为：

(6)

$\left(\boldsymbol{K}_{\mathrm{II}}+\boldsymbol{K}_{\mathrm{IB}} G^{\mathrm{ee}}\right) E_{\mathrm{I}}^{\mathrm{S}}=\boldsymbol{S}_{\mathrm{I}}-\boldsymbol{K}_{\mathrm{IB}} G^{\mathrm{ee}} E_{\mathrm{I}}^{\mathrm{P}},$

其中，E_P和E_S分别为一次场和二次场；Ⅰ和B则分别表示计算区域的内部单元以及边界单元的电场.

接收点处得到的二次磁场可以由式(7)得到：

(7)

$H^{\mathrm{S}}=\sum\limits_{i=1}^M \Delta \sigma G_{\mathrm{i}}^{\mathrm{em}} E_{\mathrm{i}},$

其中，G_i^em表示电磁场中电场与磁场的格林函数，E_i表示总场.GEM-2采集得到的数据为二次场与一次场的比值，包括同相分量与正交分量.其结果会在初始值的基础上扩大10⁶倍，即：

(8)

$\frac{H^{\mathrm{S}}}{H^1} \times 10^6,$

式中，H¹表示接收点处的一次磁场.具体计算时取该值的模，实际结果输出为无量纲单位10^-6.

2.2 基于高斯牛顿法的最小二乘反演

由于反演问题的不适定性，为获得稳定的最优解，三维反演的目标函数通常采用Tikhonov正则化方式来构建，目标函数可定义为：

(9)

$\boldsymbol{\varphi}(\boldsymbol{m})=\boldsymbol{\varphi}_{\mathrm{d}}(\boldsymbol{m})+\beta \boldsymbol{\varphi}_{\mathrm{m}}(\boldsymbol{m}),$

其中φ_d(m)为数据拟合差函数，φ_m(m)为模型正则化函数，β为Tikhonov正则化参数, 用来平衡数据中的非拟合项，从而使反演稳定.

为保证反演中模型参数的物理意义，模型参数m通常设置为m=lnσ.数据拟合差函数可表示为：

(10)

$\boldsymbol{\varphi}_{\mathrm{d}}(\boldsymbol{m})=\frac{1}{2}\left[\boldsymbol{H}^{\mathrm{S}}-\boldsymbol{d}^{\mathrm{obs}}\right]^{\mathrm{T}} \boldsymbol{C}_{\mathrm{d}}^{-1}\left[\boldsymbol{H}^{\mathrm{S}}-\boldsymbol{d}^{\mathrm{obs}}\right],$

其中 d^obs为观测的电磁场分量，H^S表示正演算子，C_d为数据协方差矩阵.模型正则化函数为：

(11)

$\boldsymbol{\varphi}_{\mathrm{m}}(\boldsymbol{m})=\frac{1}{2}\left[\boldsymbol{m}-\boldsymbol{m}_{\mathrm{ref}}\right]^{\mathrm{T}} \boldsymbol{W}_{\mathrm{m}}^{\mathrm{T}} \boldsymbol{W}_{\mathrm{m}}\left[\boldsymbol{m}-\boldsymbol{m}_{\mathrm{ref}}\right],$

其中m_ref为反演计算中参考模型，包括模型的先验信息.W_m为模型平滑矩阵.对于最小结构模型泛函，W_m为单位矩阵；对于平滑结构泛函，W_m为一阶差分矩阵.

将目标函数在当前模型变量处进行二阶Taylor展开，忽略展开式的高阶项，即可得到高斯-牛顿法(GN)的法方程为：

(12)

$\boldsymbol{h} \boldsymbol{\delta}=\boldsymbol{g},$

其中g为目标泛函的梯度，h为目标泛函的海森矩阵.它们分别可以表示为：

(13)

$\boldsymbol{g}=\boldsymbol{J}^{\mathrm{T}} \boldsymbol{C}_{\mathrm{d}}^{-1}\left[\boldsymbol{H}^{\mathrm{S}}-\boldsymbol{d}^{\mathrm{obs}}\right]+\beta \boldsymbol{W}_{\mathrm{m}}^{\mathrm{T}} \boldsymbol{W}_{\mathrm{m}}\left[\boldsymbol{m}-\boldsymbol{m}_{\mathrm{ref}}\right],$

(14)

$\boldsymbol{h}=\boldsymbol{J}^{\mathrm{T}} \boldsymbol{C}_{\mathrm{d}}^{-1} \boldsymbol{J}+\beta \boldsymbol{W}_{\mathrm{m}}^{\mathrm{T}} \boldsymbol{W}_{\mathrm{m}},$

其中 J为灵敏度矩阵.模型可以通过式(15)进行更新：

(15)

$\boldsymbol{m}_{k+1}=\boldsymbol{m}_k+a \boldsymbol{\delta},$

其中，a表示步长.

2.3 理论模型验证

为验证本文反演程序的正确性，建立理论模型进行验证，两个尺寸大小为17.5 m×10 m×3 m的六面体状高阻体，电阻率为100 Ω·m，其中心点坐标分别为(36.25，27.5，2.5)m和(36.25，45，2.5)m，嵌于电阻率为33 Ω·m的均匀半空间中.沿y轴布设7条测线，测线间距5 m.测线的观测间隔为2.5 m，选取18325 Hz和5325 Hz这两个频率下的观测数据分别进行一维反演和三维反演.一维反演结果为依据Huang和Won(2000)所提出的半空间模型得到的视电阻率，如图 3所示，视电阻率描绘出了高阻异常体在水平方向上的分布.

显示原图|下载原图ZIP|生成PPT

图3 18325 Hz和5325 Hz下的理论模型一维反演结果

Fig 3 The 1D inversion results of the theoretical model at 18325 Hz and 5325 Hz

在进行三维反演时，依据Liu等(2022a)计算所得，设定足迹的大小为42.5 m×42.5 m×13 m.用2.5 m×2.5 m×1 m的均匀离散单元离散建模体积，选取33 Ω·m均匀半空间作为初始模型.异常体的三维结构可以通过垂直切片与图 4所示的电阻率为60 Ω·m包络面相结合来描绘.研究表明，三维反演结果显示出了比一维反演结果更高的电阻率对比度.特别是60 Ω·m包络面的三维反演结果与合成模型边界吻合较好.

显示原图|下载原图ZIP|生成PPT

图4 理论模型三维反演结果的垂直切片图和包络面效果图

Fig 4 The 3D inversion results of theoretical model by vertical slices and envelope surface

3 现场测试及结果

为验证本文三维反演算法在实际应用中的可靠性，我们将其运用到上海某地质调查中，该调查的动机为将某地块复垦为正常的农用地，为此需要对该地块进行有关的建筑垃圾进行探测，以确保地块内无大型的建筑垃圾存在，从而为后期的工程实施提供依据.为此，我们选取现场的GEM-2实测数据进行反演，并将结果与传统一维反演结果进行对比，本次现场数据采集中，测量人员手持GEM-2仪器，使仪器始终保持在离地面约1 m高的地方且尽可能保持水平.GEM-2的大致测量航迹图如图 5所示.

显示原图|下载原图ZIP|生成PPT

图5 GEM-2测量航迹图

Fig 5 The survey track map of GEM-2

在现场测试中，我们同时采集了5种频率下的数据，分别为425 Hz，1525 Hz，5325 Hz，18325 Hz以及63025 Hz.由于在63025 Hz下采集得到的高频数据不理想，因此我们仅选取425 Hz，1525 Hz，5325 Hz以及18325 Hz下的正交分量数据进行研究.其一维反演结果如图 6所示.由结果可见，一维反演结果仅能大致描绘高阻异常体在水平方向上的分布，且其电阻率对比度不高，在圈定相对高阻地带时存在困难.

显示原图|下载原图ZIP|生成PPT

图6 425 Hz, 1525 Hz, 5325 Hz以及18325 Hz下的实测数据一维反演结果

Fig 6 The 1D inversion results of the field data at 425 Hz, 1525 Hz, 5325 Hz and 18325 Hz

图 7分别显示了z=1.0 m，z=2.2 m，z=3.4 m以及z=4.6 m四种深度下三维反演结果的平面图.结果显示，x=40 m，y=20 m处在3 m至5 m的深度内显示出较为明显的高阻地带.如图 8所示，其三维结构可以通过垂直切片来描绘.对比两种反演结果可以看出三维反演结果比一维反演结果具有更好的横向特征和更高的电阻率对比度，且随着深度的增加，三维反演得到的相对高阻地带面积逐步加大，与垂直切面图显示的高阻地带区域更加吻合，其在划定异常体大小及位置上具备明显的优势.

显示原图|下载原图ZIP|生成PPT

图7 不同深度下的三维反演结果

Fig 7 3D inversion results at different depths

显示原图|下载原图ZIP|生成PPT

图8 三维反演结果的垂直切片图

Fig 8 Vertical slice plot of the 3D inversion results

为验证其精确性，我们针对反演结果中的相对高阻地带(如图 7所示的黑色边框圈定区域)进行开挖验证，开挖结果如图 9所示，结果显示该高阻异常表征为遗留的房屋基础.本文的三维反演算法在实际应用中的可靠性得到验证.

显示原图|下载原图ZIP|生成PPT

图9 挖掘结果显示的高阻异常体

Fig 9 The high-resistance anomalies revealed by excavation results

4 结论

城市地下空间探测是城市未来发展中的重要环节，城市电磁法三维反演在城市地下空间探测中具备优势.本文介绍了便携式频率域电磁探测仪的基本构成和工作原理，阐述了基于足迹计算域的截断边界矢量有限元方法及高斯牛顿法的基本理论，并通过理论模型测试本文三维反演算法的正确性.依据上海某测试场地的GEM-2实测数据进行三维反演，并与一维反演结果进行了对比，经后期开挖取证，验证结果表明本文所提出的三维反演算法具有正确性及有效性，因而能满足城市地下空间探测以及实际施工所需要的准确位置等要求.便携式频率域电磁法三维反演应用于城市地下空间探测中具备可行性和优越性.

本文三维反演所用的GEM-2实测数据由上海元易勘探设计有限公司提供, 在此表示感谢.

References

Publishing order | Descend order by publishing year | Descend order by cited within

Chen J. 2023. Early signal acquisition and shallow surface imaging of urban towed transient electromagnetic method[Ph. D. thesis] (in Chinese). Changchun: Jilin University.

Chen

S J

, Hu

X Y

, Peng

R H

. Review of urban underground pipeline detection. Progress in Geophysics, 2021, 36 (3): 1236- 1247.

模态框（Modal）标题

Abstract

Cite this article

0 引言

1 便携式频率域电磁探测仪

图1 频率域电磁探测仪的工作原理图

2 方法理论

2.1 基于足迹计算域的截断边界矢量有限元方法

图2 (a) 截断边界矢量有限元以及基于足迹计算域的截断边界矢量有限元的计算区域，蓝线、黄线以及红线分别表示矢量有限元；(b)截断边界矢量有限元的计算方法

2.2 基于高斯牛顿法的最小二乘反演

2.3 理论模型验证

图3 18325 Hz和5325 Hz下的理论模型一维反演结果

图4 理论模型三维反演结果的垂直切片图和包络面效果图

3 现场测试及结果

图5 GEM-2测量航迹图

图6 425 Hz, 1525 Hz, 5325 Hz以及18325 Hz下的实测数据一维反演结果

图7 不同深度下的三维反演结果

图8 三维反演结果的垂直切片图

图9 挖掘结果显示的高阻异常体

4 结论

References