多层岩土环境的溶洞智能识别方法研究

巴泽群 * 周玉锋 周其勋 ** 周立成 * 刘逸平 * 刘泽佳 * 汤立群 *

*(华南理工大学土木与交通学院，亚热带建筑与城市科学全国重点实验室，广州 510640)

†(广州广佛肇高速公路有限公司，广州 510289)

**(广州市高速公路有限公司，广州 510030)

喀斯特地貌广泛分布在我国西南部和中部地区。这种地貌是由地下水对可溶性岩石进行溶蚀、沉积等作用形成的，因此容易产生地下溶洞等地质现象。在实际工程中，溶洞会影响桩基施工安全，严重时可能导致地面坍塌。因此，在施工前了解地下环境并排除溶洞等安全隐患，有助于确保施工安全并减少经济损失。

目前，传统溶洞探测方法包括钻孔取样法、弹性波CT 法[1]、地震反射波法[2-4]等机械类方法以及电磁波跨孔CT 法[5]、探地雷达法[1]、高密度电法[1]、电阻率层析法[6-8]、瞬态电磁场法[1,9]等电磁类方法。这些方法技术成熟且应用广泛，但仍存在着一些缺点，如主要以二维截面探测为主，数据反演和解释需要依赖专业知识和专家经验等。亦有方法可进行三维溶洞形态反演[10]，但该方法建立在已对溶洞进行定位的前提下，且数据反演较复杂。因此，传统溶洞探测方法在探测效率、经济成本、时间成本和智能化方面亟待改进。

在地质勘探领域，已有部分学者利用智能方法开展了一些探索性的研究工作。例如，Zhang 等[11]使用U-NET 网络实现地球物理勘探信息的智能反演；Di 等[12]基于支持向量机/多层感知机等机器学习方法实现地震数据的智能解释和反演；Zheng 等[13]用主成分分析、K 聚类和线性判别式分析方法，使用声速信号、电阻率、密度、孔径数据等测井相信息进行地层类型的预测和识别；Xiong 等[14]使用深度自编码网络推断矿物储备情况；Huang 等[15]基于交叉小波变换处理振动时序信号获得振动相位谱和功率谱，并通过卷积神经网络来进行爆破位置的预测；Li 等[16]使用地震属性数据结合随机森林、遗传算法–支持向量机等机器学习方法，实现了煤层地下结构的预测；Yang等[17]基于卷积神经网络和支持向量机等算法，通过振动信号的时频特征对地震源深度进行了识别。这些研究使用了人工智能方法代替人工数据反演与解释，表明了人工智能方法在地质勘探领域的可行性和应用潜力。

本文作者团队在前期对高效、智能的溶洞探测技术进行了初步研究，提出了一种适用于单一土层环境的基于机器学习算法的球形溶洞定量探测方法[18]。在此基础上，本文将研究场景推广到更符合实际工程的多层岩土环境，并设计了基于(bidirectional long short-term memory, BiLSTM)神经网络的深度学习模型来实现多层岩土环境下的溶洞三维定量智能探测。本文技术路线如图1所示，主要步骤如下：（1）借鉴浅层地震反射波法，建立含有单个无填充球形孔洞的多层岩土结构模型，并计算桩锤激振下地表的加速度响应信号；（2）基于应力波传播的基本原理对响应信号进行分析以提取特征数据；（3）针对不同溶洞工况进行大量建模，构建用于深度学习的数据集；（4）设计并优化深度学习模型，实现多层岩土环境下溶洞的三维位置和大小的定量智能识别。本文提出的深度学习模型可不依赖于专业知识和专家经验，为高效、智能化的溶洞探测方法的工程应用奠定基础。

图1 技术路线Fig.1 Technical roadmap

1 含溶洞的多层岩土结构有限元建模

在实际工程中，大部分地下岩土体为分层的多层结构，具有多样的岩土层种类和复杂的岩土层分布，这导致多层岩土环境下对溶洞大小和位置的预测变得非常困难。由于获得实测数据需要高昂的经济成本和时间成本，本文采用有限元方法对多层岩土环境下的大量溶洞工况进行建模，从而构建不同溶洞工况下的响应数据库，为基于深度学习的溶洞定量识别方法的探索奠定基础。

1.1 有限元模型

本文的建模过程参考了广东省广州市广佛肇高速公路SG03 段的地质勘探数据，根据其地勘数据对岩土层种类、岩土层组合以及各层岩土的厚度变化范围进行建模。本文建立的有限元模型如图2(a)所示。其中，整个岩土体模型为一个边长30 m 的立方体，分为8 层，按照从上往下的顺序进行编号，1～6 层为土层，7～8 层为岩层。为简化建模过程，土层和岩层均设置为各向同性、均匀密实的材料，且土层采用Mohr–Columb 本构模型，岩层采用Drucker–Prager 本构模型，详细参数见表1。同时，为简化模型，假设土体内含一个无填充、真空的球形孔洞。激振源使用一个高度0.5 m，直径1.5 m，重约6.8 t 的圆柱形桩，其材料采用线弹性本构，密度为7800 kg/m3，弹性模量为200 GPa，泊松比为0.3。激振形式采用桩锤自由下落的方式，撞击地表时的速度设为4.42 m/s。如图2(b)所示，为了获取地表的加速度响应，以激振点为圆心，在距离其10 m 的采样圆上均匀布置6 个加速度传感器，按照顺时针方向依次编号为1～6 号，且传感器的采样频率为1000 Hz，采样时间为0.6 s。值得注意的是，为了防止应力波在模型边界上的反射，土体侧面和底面设置为无限元边界。

表1 岩土层分布示例Table 1 Example of rock-soil structure

图2 有限元模型示意图Fig.2 Diagram of finite element model

1.2 工况设置

如表2 所示，球形孔洞参数包括洞顶与激振线之间的水平偏移距离、洞顶深度和孔洞直径。其中，洞顶水平偏移距离的取值范围为1～10 m，取值间隔为1 m；洞顶深度取值范围为5～21 m，取值间隔为3 m；孔洞直径取值范围为1～10 m，取值间隔为1 m。上述3 个取值变量依照取值要求随机组合，生成大量涵盖不同孔洞工况的模型。同时，如表3 所示，设置9 种典型的多层岩土结构，与表2 的孔洞模型随机组合，即可构建9 种不同岩土结构下含有不同孔洞工况的海量模型。

表2 溶洞工况设置Table 2 Settings of karst cave parameters

表3 岩土层结构设置Table 3 Settings of rock–soil structures

通过统计，每个有限元算例在搭载了Intel Core i5-10400F 型CPU 的计算机上的计算时间平均为5 h 以上。因此，构建数据集需要大量的算力和时间开销。针对有限元算例计算时间长的问题，本文将有限元算例进行等效工况扩充。具体地，如图2(c)所示，虚线圆框代表溶洞，蓝色圆和绿色圆分别代表两组加速度传感器，蓝色对应x – y坐标系，绿色对应x′-y′坐标系。可以看出，以不同坐标系描述的同一个溶洞具有不同的空间位置。因此，对于一个有限元算例，本文在其加速度采样的圆圈上均布了288 个节点，即每隔1.25°方位角布置一个加速度传感器，进而通过旋转变换采样不同的传感器组合（6 个为一组）的方式来扩充溶洞工况数据，即每一个有限元算例可以等效扩充为288 个工况。最终，每种岩土结构对应的溶洞工况数大致相同，平均约包含2880 条数据工况，总数据量为25 344 条。

2 特征工程

为了提高深度学习模型的识别效果，本文基于应力波反射的基本原理，对地表加速度响应信号进行特征工程，力求挖掘与预测目标相关性强的特征并作为模型的输入，从而提高模型的预测精度。

地下多层岩土结构中存在许多波阻抗界面，如岩土层分界面、溶洞表面等，当应力波遇到这些界面时，会发生反射现象。因此，可通过计算特征界面的应力波反射信号的行程与时间，捕捉受溶洞和地下岩土体结构影响的加速度信号，进而进行特征工程。

首先，对于不含水的各向同性致密岩土体，其纵波波速的计算公式为

式中，vp为纵波波速，ν为泊松比，E为弹性模量，ρ为密度。基于式（1），结合岩土体结构中各层的厚度，可以估算不同界面应力波反射信号的行程和走时。应力波传播路径示意图如图3所示，其中图3(a)展示了应力波被激发后传播至岩土层交界面被反射回地表并被传感器接收的过程；图3(b)展示了应力波被激发后传播至溶洞顶点被反射回地表并被传感器接收的过程。考察岩土层界面反射情况时，取传感器与激振点连线的中位线与各岩土层交界点的位置作为反射参考点，以计算应力波的行程。考察溶洞表面反射情况时，近似取溶洞洞顶作为反射参考点进行计算。以一列应力波的传播过程分析作为示例，对于任意传感器m，其自激振作用起到接收到岩土层界面反射信号的时间可按式（2）进行计算

图3 应力波路径示意Fig.3 Diagram of stress wave propagation

对于任意传感器m，其自激振作用起到接收到溶洞洞顶反射应力波的时间可按公式(3)进行估算

其中，x，y，z代表溶洞洞顶反射点的空间坐标；xr，yr和zr代表岩土层反射点r 的空间坐标；xm，ym和zm代表传感器m 的空间坐标。hi代表应力波经过岩土层i的垂直距离，vpi代表岩土层i中纵波的波速。基于式（2）和式（3），本文计算了各个岩土层交界面和溶洞洞顶的应力波反射信号的接收时间，如表4 所示。其中包含了洞顶深度为2 m，6 m 和10 m 的溶洞模型的计算结果，所含溶洞的半径均为2 m。图4(a)给出了不同溶洞工况下的多层岩土体结构加速度响应，其中黑色曲线为无溶洞情况下的多层岩土结构的加速度响应时程曲线，红色、蓝色和绿色曲线为溶洞深度递增的多层岩土结构的响应时程曲线，且4 条图线所对应的岩土体分层情况一致。

表4 反射信号走时Table 4 Travel time of reflected stress wave

图4 加速度响应信号对比Fig.4 Contrast of acceleration response signals

综合表4 和图4(a)的结果可知，岩土层界面反射信号和溶洞洞顶反射信号的接收时间与加速度时序信号的波峰波谷时间接近，虽然无法准确说明每个波峰波谷包含哪些反射波成分，但是可以推断在加速度曲线的极值数据中蕴含着与岩土层界面位置和溶洞深度等相关的信息。同时，若溶洞直径增大，其反射面增大，响应信号波峰波谷对应的加速度大小以及所对应的时间会相应发生变化。因此，可以提取加速度信号的极值数据作为输入数据的特征。具体地，将每个工况下的6 个传感器的时序信号的极值点全部提取出来，并按照出现的先后顺序进行排列，形成极值序列数据。通过统计，平均每个加速度时序信号包含19 个极值点，故极值特征的数据格式为[19, 6]。由以上的分析可知，极值特征与溶洞的深度和大小有较强的相关性。

图4(b)展示了对应于一个含有溶洞的多层岩土结构的一组传感器的时序信号，其中溶洞洞顶的水平偏移距离为5 m，溶洞半径为2 m。观察不同传感器的时序信号可知，由于传感器空间方位的不同，传感器接收的响应信号也有一定的差异，这种差异来源于传感器与溶洞的不同距离。可以推断，不同空间方位的传感器采集的加速度时序信号蕴含着溶洞的位置特征信息，即可以提取一组6 个不同空间位置传感器的时序信号，作为输入数据的时空特征。具体地，首先提取每个传感器采样信号的0.15～0.60 s 时间段的加速度数据，然后按照传感器编号顺序即空间位置顺序排列时序信号，得到格式为[6, 451]的时序数据作为后续深度学习模型的输入。

3 多层岩土环境的溶洞智能识别模型

深度学习是人工智能研究领域的重要方法，其基本原理是通过搭建人工神经元组成的神经网络，建立输入与输出的映射关系，解决各种模式识别问题。目前，深度学习方法广泛应用于文本翻译、图像识别、时间序列预测等领域。其中，BiLSTM 作为循环神经网络（recurrent neural network, RNN）的一种改进模型，有着强大的时间序列数据学习能力。考虑到本文提取了以时间序列数据为主的输入特征，因此引入BiLSTM 网络，以BiLSTM 网络为核心模块设计了双数据通道的分离架构网络模型，对多层岩土环境下的不同溶洞工况进行定量识别。

3.1 BiLSTM 的原理

BiLSTM 由两条方向相反的长短期记忆（long short-term memory, LSTM）神经网络构成，其结构如图5(a)所示，可以从正反两个方向学习时间序列的特征。LSTM 网络的基本单元为LSTM 神经元，其基本原理如图5(b)所示。其中，LSTM 神经元由遗忘门ft，输入门it，输出门ot和细胞临时状态Ct∼组成。LSTM 神经元将上一时刻的输出值ht−1和本时刻的输入值xt进行矩阵拼接，并将拼接矩阵激活后得到遗忘门、输入门、细胞临时状态和输出门。遗忘门控制上一时刻细胞状态Ct−1的信息流入；输入门决定了本时刻的细胞临时状态Ct∼对本时刻细胞状态Ct的影响程度；输出门决定了LSTM 神经元的输出值。将LSTM 神经元连接起来，通过细胞状态的传播，即可实现序列数据的长时记忆，从而全面、完整地提取时间序列的特征。

图5 BiLSTM 示意图Fig.5 Diagram of BiLSTM

3.2 网络架构

本文采用BiLSTM 网络作为多层岩土环境下溶洞识别深度学习模型的核心模块，基于该网络设计了如图6(a)所示的双数据通道融合架构、如图6(b)所示的双数据通道分离架构和如图6(c)所示的单数据通道架构。如图6(a)所示，两支通道的输入数据分别为时序信号序列和极值特征，其分别通过BiLSTM 网络层提取特征，随之与多层全连接网络连接，在相邻网络层间设置随机失活层（dropout）以避免网络过拟合，最后进行两个通道特征矩阵的拼接融合，并通过全连接层得到最后的预测结果。预测结果为4 个标签X,Y,Z和D，分别对应了溶洞洞顶的三维坐标和溶洞直径。其网络各层神经元个数及激活函数设置如表5 所示。如图6(b)所示，双通道分离架构的输入数据分为两支，分别为时序信号序列和极值特征，两支数据通道的训练同步进行。且两支数据通道的网络层神经元权重、偏置等网络参数的更新优化互不影响，其中以时序信号序列为输入的分支对应的输出为溶洞洞顶水平坐标X,Y；以极值特征为输入的分支对应的输出为溶洞洞顶深度Z和溶洞的直径D。该网络架构的各个网络层神经元个数与激活函数设置如表6 所示。图6(c)上部展示了以时序信号序列为输入的单数据通道架构，其通过BiLSTM 层提取特征后，随之连接多层全连接层，同时在相邻全连接层间设置随机失活层（dropout），最终输出溶洞洞顶的空间坐标X，Y，Z和溶洞的直径D；图6(c)下部为以极值特征为输入的单数据通道架构，其网络设置与以时序信号序列为输入的单通道架构相同。单通道架构的各层神经元个数与激活函数设置如表7所示。

表5 双通道融合架构参数设置Table 5 Settings of dual-channel fusion architecture

表6 双通道分离架构参数设置Table 6 Settings of dual-channel separated architecture

表7 单通道架构参数设置Table 7 Settings of single-channel architecture

图6 不同网络架构Fig.6 Different model architectures

3.3 模型训练

基于9 种岩土结构下共25 344 个工况的数据集，以8:2 的数据量比例划分训练集与测试集。其中，训练集数据量为20 275 条，测试集数据量为5069 条。从训练集中划分出10%的数据作为验证集，并在训练集上训练数据，在验证集上观察模型训练效果，在测试集上检验预测效果。在训练过程中，选择平均平方误差（meansquare error, MSE）作为损失函数，其表达式为

其中，ypi表示第i个数据的预测值，yti表示第i个数据的真实值，n代表总数据量。网络优化器选择Adam 优化器，其算法由Kingma 等[19]提出，能够实现自适应学习率调节，并且初始学习率设置为0.001。网络训练轮次上限设置为500 次，收敛条件为：若连续20 次训练迭代后，验证集的MSE 减小值小于0.1，则模型达到收敛，训练终止。平均每个网络模型训练10 次，每次训练结果相互独立。

图7(a)展示了双通道融合架构训练过程中损失函数的变化情况，可以看出损失函数稳定下降，收敛程度良好；图7(b)展示了双通道分离架构的情况，其损失函数下降速度先快后慢，最终收敛；图7(c)和图7(d)分别展示了单通道架构（时序信号序列）和单通道架构（极值特征）的训练情况，其中单通道架构（极值特征）的收敛状态下的损失函数值约为3，相比其他网络训练效果较差。由图7 可知，4 种网络架构模型都能较好收敛。双通道架构的收敛效果更优，以双通道分离架构网络模型的收敛程度为最好。

图7 各网络训练效果可视化Fig.7 Visualization of training process

3.4 评价标准

本文提出的多层岩土环境下溶洞识别模型的识别目标为溶洞的三维位置和大小，因此评价网络性能的标准是溶洞参数的预测准确率，对溶洞预测准确与否的标准所示为

其中，r表示容许误差，下标pred 和true 分别表示预测值和真实值。当预测结果同时满足式（5）～式（8）时，可判定模型完成了溶洞工况的准确预测；当预测结果满足上述任一公式条件时，可判定该模型对于溶洞的对应参数预测准确。模型的预测准确率计算式为

其中，Ntrue表示预测准确的样本数，Ntotal表示参与预测的总样本数。

4 结果与讨论

本节对基于BiLSTM 的双数据通道分离架构、双数据通道融合架构和单数据通道架构模型开展对比分析，各模型在3 m 容差下的预测结果如表8 所示。其中，双通道分离架构网络在3 m 容差下的综合预测准确率最高，达到了98.23%，高于双通道融合架构网络的95.95%综合预测准确率、以时序信号作为输入的单通道架构网络的95.51%综合预测准确率和以极值特征为输入的单通道架构的17.99%综合预测准确率。对比双通道分离架构和双通道融合架构的预测结果，可知双通道分离架构充分利用了时序信号与溶洞坐标X和Y的相关性以及极值特征与溶洞深度Z和溶洞大小D的相关性，且双通道分离架构无需特征融合，规避了双通道融合架构特征融合困难、依赖调参的缺点。观察双通道分离架构网络模型与以极值特征为输入的单通道架构网络模型的预测结果，可知极值特征与溶洞深度Z和溶洞大小D的相关性较强，与溶洞水平位置X和Y的相关性较弱。观察基于时序信号的单通道架构网络综合预测结果，其预测准确率相较双通道融合架构网络仅有不到0.5%的差距；而与双通道分离架构网络的综合预测准确率相比，存在约2.7%的差距。可以推断时序信号包含了足够的信息，具有不同水平方位、深度和大小的溶洞对溶洞表面反射应力波的行程和走时会产生影响，因此通过对不同空间位置加速度传感器信号进行完整采样，能够从中挖掘溶洞的三维位置和大小特征，从而实现对溶洞的三维位置和大小的较高精度的识别，而加入极值特征有助于提高识别准确率，说明了本文特征工程的有效性。

表8 预测结果Table 8 Results of prediction

综上所述，采用时序信号和极值特征的双通道输入有助于提高模型预测准确率，其中，双通道分离架构模型能充分利用输入特征与相应识别目标的强相关性，因此其能达到最优的预测准确率。

5 结论

本文提出了一种基于深度学习的多层岩土环境的溶洞智能三维定量探测方法，所得到的主要结论如下：（1）地表加速度响应时序信号包含足够的信息，以时序信号为输入的单通道架构模型在3 m 容差下对溶洞位置和大小的识别准确率达到了95.51%，加入极值特征可提高模型识别准确率；（2）基于BiLSTM 的双通道分离架构网络模型在3 m 容差范围内对多层岩土环境下溶洞的三维位置坐标和直径大小的综合预测准确率达到了98.23%，实现了多层岩土环境下溶洞的定量智能预测。本文研究的不足之处在于：（1）岩土结构和溶洞的组合工况可能性极多，本文受到算力的限制，只能在一些有限的组合工况下进行方法的验证；（2）仅进行了基于有限元数据的可行性研究，未在实测数据上对方法进行验证。