0.   引　言

岩溶地质作用是地表水和地下水对可溶性岩进行的以化学溶蚀作用为主的综合地质作用(包括水的化学溶蚀、机械侵蚀、物质的运移和再沉积)。山东临沂地区岩溶发育，长期抽采地下水极易产生地下水降落漏斗^[1]，并伴随发生岩溶塌陷。自1980年以来，临沂城区发生过19次岩溶塌陷，形成的塌陷坑多达30个，岩溶塌陷具有不确定性、隐蔽性和突发性^[2]，对当地居民的生命、财产安全构成极大的威胁^[3−4]，同时，岩溶塌陷还可能连通地表水与地下水，引发地下水穿层，致使地下水遭受污染。因此在工程选址时应重视对易发塌陷等地质灾害的排查。

近年来，地球物理勘探技术在工程地质施工中应用广泛。例如：在高山峡谷的桥址区选址时，采用浅层地震、音频大地电磁、对称四极测深等方法^[5]效果较好；复杂地质条件下，岩溶发育区的选址可利用高密度电法、浅层地震、大地电磁和孔间弹性波CT等方法^[6−7]。三维高密度电法的基本原理与二维高密度电法类似，都是以导电性的物性差异为基础，研究在稳定电流场作用下，地下空间中传导电流的分布规律^[8]。三维高密度电法被广泛应用于探测隐伏构造、海水入侵界面、土石坝渗漏探测等问题^[9−11]，在岩溶探测中也广泛使用^[12−13]；孟凡松等^[14]把二维高密度数据整理转化成三维反演数据，并将反演结果通过三维成像软件实现任意三维地电结构切片，从而获得丰富的地电结构信息；通过三维高密度电法的数值模拟，在空间定位上具有较好的准确度^[15]；王钧等^[16]利用Res3dinv软件进行正反演数值模拟，获得了异常体的响应规律，分析了高密度电法存在的体积边界效应。

溶洞的规模小到几厘米，大到几十公里^[17]，无疑增加了岩溶探测的难度，部分学者利用反演以及数值模拟等技术，对岩溶进行勘测，郑智杰等^[18−19]认为当岩溶管道深度大于15倍电极距时，高密度电法不易探测到该深度的岩溶管道；而小于10~11倍电极距时，能较精准的定位岩溶管道位置，为本次工作提供了参考。朱欣月等^[20]自主研发了的岩溶溶蚀演化模拟软件进行数值模拟，并以三维成果进行展示，韩凯等通过高精度示踪技术和大功率充电方法定位岩溶地下河管道的分布特征和展布方向，成果显著^[21]；徐佳等^[22]在Voxler平台下对高密度电法的数据体进行了三维可视化，更直观地表达了地电信息，通过高密度电法对岩溶进行勘察的实例更多^[23−24]，这些文章为本次研究提供了丰富的理论和实践依据。

作者在查阅文献时发现，目前没有学者在岩溶水降落漏斗中心利用高密度电法对岩（土）层界面以及岩溶塌陷开展系统研究。本次在预测的岩溶水降落漏斗中心开展三维高密度电法工作，划分了第四系与灰岩的界面，并对溶洞进行了推测和验证，验证效果好，解决了当地岩溶发育区工程选址问题，为以后工程建设以及开展此类科学研究提供了依据。

1.   地质环境背景

1.1   区域水文地质背景

研究区位于山东省临沂市城西区，属于临沂单斜水文地质单元，地形平坦，坡度小，是灰岩裂隙水的富水地带，地表建筑物及人类活动密集，经调查，2003年时存在一个以梨杭村—后岗头村为中心的岩溶水降落漏斗(图1)。该岩溶水降落漏斗呈椭圆状，长轴北西走向，长轴长度约17 km，短轴长约11 km，地下水埋深降落高程差为2~20 m。在降落漏斗范围内，自1980－2006年，城区内地下水超采，水位大幅下降，水动力场发生变化，据记载，曾发生过18次岩溶塌陷和多次地下水污染，形成的塌陷坑多达30余处，且岩溶(塌陷)的位置大多位于城区或居民地，对人民群众的生命财产安全造成严重威胁^[25]。本次工作选址区位于该地下水降落漏斗边缘，郭士昌等^[26]利用有限差分算法进行空间和时间的离散，对临沂城区的岩溶地下水流场进行了数值模拟，模拟期为2001年1月至2007年5月，结果认为：模拟期末岩溶水降落漏斗由梨杭—后岗头村转移到齐家庄村附近，漏斗中心水位标高44.42 m，扩展面积1.1 km²。经调查，2012年后，地下水开采受控，水位变化不大，岩溶塌陷明显减少，但可能仍存在以往遗留的溶洞，本次城建选址区位于齐家庄村附近，属于岩溶塌陷易发区，给工程建设带来了困难。

图  1  临沂市城西岩溶水降落漏斗(2003年6月)(据^[25])

Figure  1.  Karst water drawdown funnel in the west of Linyi City (June, 2003)(Liqingchun et al., 2005^[25])

下载: 全尺寸图片 幻灯片

1.2   选址区水文地质特征

选址区地表大部分为回填渣石土，仅在西北部出露新生代第四系临沂组灰黄色黏土质粉砂、含黏土粉砂；北端大部分被混凝土覆盖。选址区地层上部为回填杂土、第四系黄色黏土质粉砂、含黏土粉砂，下部为古生代奥陶系马家沟群地层灰色厚层、中厚层泥晶灰岩、云斑灰岩夹藻凝块灰岩、粉晶白云岩，地层起伏不平，可能存在岩溶。选址区主要含水层为奥陶系马家沟群地层岩溶裂隙含水层和第四系砂砾孔隙含水层，直接受大气降水及地表水补给，此含水层厚度大，岩溶裂隙较发育。

2.   工作基础

2.1   岩石物性特征

选址区地层简单，岩性单一（表1），临沂组第四系黏土、砂土电阻率较低，变化范围在10~50 Ω·m，平均值约17 Ω·m；选址区内回填物为渣石土，其含水量很高，回填物电阻率受回填物及含水量的影响，变化范围较大，电阻率值在15~100 Ω·m，平均值约25 Ω·m；奥陶系灰岩电阻率较高，变化范围在150~3370 Ω·m，平均值约100 Ω·m；溶洞分为两类，一类为含水、含泥的溶洞，电阻率一般小于50 Ω·m，另一类为未充填的溶洞，电阻率一般大于3 000 Ω·m，因此，利用高密度电法划分地层界面和寻找岩溶是具有物性基础的。

表  1  区内岩、土层电性统计表

Table  1.  Electrical properties of rocks and soil layers in the study area

地质代号	岩性	电阻率ρ/Ω·m
		平均值	变化范围
Q	回填土	25	15~100
	临沂组黏土、砂土	17	10~50
Omt	马家沟组土峪段 灰岩、白云岩	1 000	150~3 370
	含水、含泥溶洞	＜50
	未充填的溶洞	＞3 000

下载: 导出CSV 

| 显示表格

2.2   场地溶洞发育条件

选址区地层主要分为两层，上半部分上覆5~12 m第四系松散层，表现为“上土下砂”的二元相结构，下半部分为基岩，岩性为马家沟组土峪段灰岩、白云岩。

3.   高密度电法工作原理及布置

高密度电法与常规电阻率法的工作原理基本相同，都是以岩土体的电性差异为基础的一种勘探方法。根据在施加电场作用下地下传导电流的分布规律，推断地下具有不同电阻率的地质体的赋存情况。高密度电阻率法的物性前提是地下介质间的导电性差异，和常规电阻率法一样，它通过A、B电极向地下供电，然后在M、N极间测量电位差△V，从而可求得该点(M、N之间)的视电阻率值。根据实测的视电阻率剖面，进行计算、分析，便可获得地下地层中的电阻率分布情况，从而可以划分地层、判定异常成因等。

本次在选址区共布设14条高密度电法测线(图2)，测线间距约12 m，测点间距2 m，采用温纳装置，隔离系数最大为25。温纳装置又称对称四极装置，A、M、N、B等间距排列，其中A、B是供电电极，M、N是测量电极，AM=MN=NB为一个电极距，电极间距按隔离系数由小到大的顺序等间隔增加，四个电极之间的间距也均匀拉开。该测量方式为剖面测量，所得数据断面为倒梯形，为保证数据体的完整性，本次采用Kriging网格化方式，将倒梯形数据补齐为矩形三维数据体。各测点使用米尺及RTK实地放点布置，质检率为7.2%，视电阻率均方相对误差M=±2.71%，符合规范要求。

图  2  选址区三维高密度电法测点布置图

Figure  2.  Layout of survey points for 3D high-density electrical method in the study area

下载: 全尺寸图片 幻灯片

4.   正反演数值模拟

4.1   数值模拟原理

为便于对异常解译，本文通过数值模拟建立形成溶洞条件下的地质模型，并研究其响应特征。本次正演计算采用Res2dmod软件，反演计算采用RES2DINV软件。其中Res2dmod采用有限差分法，RES2DINV采用阻尼最小二乘法。

有限差分和有限元法是使用2D模型是将地下空间划分为一定数量的块体，采用矩阵网格，A.Dey和Morrison针对有限差分公式做了改进，提高了计算视电阻率值的精度^[27−28]，有限差分方法使用的是计算矩形网格节点的电位值，矩形网格节点组成复杂的地质构造建模，最终用以反映电场响应特征^[29]。

阻尼最小二乘法是求解非线性二乘问题的标准方法，最小二乘法主要来源于拟合问题，即实测一组数据，然后通过含参数的方程，计算出一组数据，不断修改完善参数和方程，将计算数据与实测数据进行观测拟合，求取误差最小值，由此确定地电模型，服务于地质工作要求。高密度电法数据反演是建立在正演基础上进行的，首先假定地电模型，将地电模型划分为多个网格单元，待求电阻率参数反演模拟处理方法是建立在圆滑约束最小二乘法的基础上进行的^[30]，计算每个小单元格内的视电阻率值，并与实测视电阻率值进行比较和优化迭代，直至计算值和测量值的均方根误差达到极限要求^[30]。

4.2   地电模型建立

在对选址区物性参数、地层和不良地质体(溶洞)进行分析的基础上，建立了区内的地电模型，如图3a，第四系覆盖层(包括回填土或砂土层)电阻率值为20 Ω·m，基岩为灰岩或白云岩，电阻率为1 000 Ω·m，因现阶段该区水位实测为62.1 m，地表高程为68 m，水位离地表约6 m，而溶洞位置距地表约20 m，因此溶洞为充水充泥溶洞，假定溶洞电阻率为30 Ω·m，在此基础上对地层(界面)和岩溶的识别上进行分析。

图  3  地电模型及正反演电阻率断面图

Figure  3.  Cross sections of geoelectric model and forward and backward resistivity

下载: 全尺寸图片 幻灯片

图3b和图3c断面上由浅到深地层均表现为电阻率逐渐升高的“G”型曲线，电阻率层状反映较明显，仅在溶洞发育区电阻率出现不同程度的扭曲，总体可划分为两个电性层，浅部为第四系砂土层或回填土，深部对应基岩为灰岩。据正演结果，第四系砂土层或回填土电阻率值为75~100 Ω·m，灰岩电阻率值为150~300 Ω·m，电性层界面反映并不明晰，且在岩溶位置电阻率等值线向深部扭曲，这也反映了溶洞与基岩电阻率值相比较偏低，但岩溶位置与等值线扭曲位置吻合程度不高；据反演结果，第四系砂土层或回填土电阻率值小于50~100 Ω·m，灰岩电阻率值为150~1000 Ω·m，等值线密度明显增大，电性层界面反演明显，在岩溶位置电阻率等值线明显向深部延展，是低电阻率溶洞的响应，在等值线变化梯度最陡处，与实际溶洞位置吻合程度高。总之，在地电模型基础上，对其进行的正反演结果表明，正演结果对界面和溶洞的反映并不明显，而反演结果对界面和溶洞有清晰的响应。分析认为反演方法对溶洞和界面的响应效果更好。

4.3   反演结果及验证

在分析地电模型的基础上，本次针对选址区开展的14条高密度电法工作采用阻尼最小二乘法进行反演，各剖面数据经整理，利用Voxler软件成图，综合形成三维立体图，将其分解为14个剖面开展分析（图4），以P10剖面为例展开论述，并开展了钻孔(ZK01)验证工作，其余剖面与P10剖面类似。

图  4  选址区三维高密度电法视电阻率三维切片图

Figure  4.  3D slice diagram of apparent resistivity of 3D high-density electrical method in the study area

下载: 全尺寸图片 幻灯片

P10剖面位于自西向东第10条剖面，剖面长度138 m。剖面范围内低阻区厚度较大（图5），视电阻率曲线杂乱，视电阻率值在10~50 Ω·m，局部有高阻段，推断为后期回填渣石土引起，以140点处高阻隆起为界分为南北两段，南段向南未圈闭，在剖面100点处厚度约32 m，向北缓慢变小，北段宽度较窄，在剖面160点处厚度约25 m，向北急剧变小，其边界约在剖面186点处；剖面北端浅部高阻段与低阻段高低交错，视电阻率值高阻在150 Ω·m以上，低阻低于40 Ω·m，厚度约3 m，推断为第四系黏土与夹石反映；深部高阻层视电阻率值大于400 Ω·m，视电阻率曲线密集，推断为致密灰岩的反映。

图  5  P10剖面反演电阻率断面图

Figure  5.  Cross section of inversion resistivity of P10 profile

下载: 全尺寸图片 幻灯片

在剖面125~134点、垂深28 m以下，视电阻率等值线稀疏且呈竖直状，该部位处于灰岩地层内，推断为落水洞充填含水致视电阻率降低引起，推断落水洞宽度约9 m、厚度大于14 m，上部与回填区连通，已被回填。

在剖面的127点位置布设钻孔（图5），经地质编录，经验证在24~28.1 m位置，见厚度为4.1 m溶洞，与推测结果基本吻合。

5.   分析讨论

5.1   对界面的响应

区内电场分布总体呈低—高两层层状分布，南段浅部低阻层厚度一般在20 m左右，由北向南厚度逐渐加大，视电阻率曲线较杂乱，视电阻率值在10~50 Ω·m之间，局部有高阻异常杂乱分布，推断为回填区回填渣石土与夹石的综合反映，根据其与下伏灰岩不同电性分布特征，划分回填界面，圈定回填等厚线(图6，图7)；北段浅部低阻层厚度一般在3 m左右，视电阻率曲线宽缓，视电阻率值在20 Ω·m左右，推断为第四系黏土、回填杂土的反映；深部高阻层视电阻率曲线密集，分布较均匀，视电阻率值在200 Ω·m以上，推断为奥陶系灰岩引起，南段埋深较深，深度在20 m以下，为开挖后回填所致，北端为原始地层，埋深在3 m左右。经施工钻孔验证，推测结果与钻孔揭露基本一致。

图  6  高密度电法对界面和岩溶的响应

Figure  6.  Response of high-density electrical method to interface and karst

下载: 全尺寸图片 幻灯片

图  7  高密度电法推测界面响应的三维图

Figure  7.  3D diagram of interface response inferred by high-density electrical method

下载: 全尺寸图片 幻灯片

5.2   对岩溶的响应

各剖面异常电性在空间分布上具有连续性，局部因岩层开挖断开，根据其分布特征，剖面中部异常处于回填区下伏奥陶系灰岩地层内，异常呈竖直低阻或倒三角形，埋深一般超过24 m，向下延深较深，底部未封闭，推断为溶洞充水引起，由于上部岩石开挖出露形成落水洞，各剖面异常在空间分布上具有连续性，根据其分布特征，将推测结果统计整理为表2，圈定一个深部岩溶发育区及一个落水洞(图8)：深部岩溶发育区位于测区中部，横贯整个场地，总体呈东西向，西侧转折为北西向，向西未圈闭，总体长度超过130 m，最宽处约8 m；落水洞位于测区东部，向东未圈闭，长度超过10 m，宽度约7 m。经施工钻孔验证，在推测岩溶发育区，有9处钻孔见溶洞(图8)，溶洞大小不一，约0.5~8.4 m。因此，三维高密度电法对寻找岩溶效果显著。

表  2  推测溶洞位置

Table  2.  Inferred location of the karst cave

线号	点位	顶部距地 表距离/m	底部距地 表距离/m	线号	点位	顶部距地 表距离/m	底部距地 表距离/m	线号	点位	顶部距地 表距离/m	底部距地 表距离/m
P1	143	−20	−27	P5	130	−42	−48	P10	131	−28	−48
P1	145	−20	−27	P5	132	−42	−48	P10	133	−28	−48
P1	147	−20	−27	P6	123	−25	−48	P10	134	−28	−48
P1	149	−20	−27	P6	125	−25	−48	P11	125	−28	−48
P1	151	−20	−27	P6	127	−25	−48	P11	127	−28	−48
P1	154	−20	−27	P6	129	−25	−48	P11	129	−28	−48
P2	126	−24	−48	P7	126	−26	−48	P11	131	−28	−48
P2	128	−24	−48	P7	128	−26	−48	P11	133	−28	−48
P2	130	−24	−48	P7	130	−26	−48	P11	135	−28	−48
P2	132	−24	−48	P7	132	−26	−48	P12	124	−25	−48
P3	126	−24	−48	P8	130	−25	−48	P12	126	−25	−48
P3	128	−24	−48	P8	132	−25	−48	P12	128	−25	−48
P3	130	−24	−48	P8	134	−25	−48	P12	130	−25	−48
P3	132	−24	−48	P8	136	−25	−48	P12	132	−25	−48
P4	126	−36	−48	P9	127	−28	−48	P12	134	−25	−48
P4	128	−36	−48	P9	129	−28	−48	P12	136	−25	−48
P4	130	−36	−48	P9	131	−28	−48	P12	137	−25	−48
P4	132	−36	−48	P9	133	−28	−48	P14	144	−22	−48
P4	133	−36	−48	P10	125	−28	−48	P14	146	−22	−48
P5	126	−42	−48	P10	127	−28	−48	P14	148	−22	−48
P5	128	−42	−48	P10	129	−28	−48	P14	150	−22	−48

下载: 导出CSV 

| 显示表格

图  8  高密度电法推测岩溶(灰白色区域)响应三维图

Figure  8.  3D diagram of karst (gray-white area) response inferred by high-density electrical method

下载: 全尺寸图片 幻灯片

5.3   岩溶塌陷发育讨论

根据岩溶形成过程中地下水的走向，岩溶分为垂直岩溶和水平岩溶两大类，垂直岩溶主要有石芽、溶沟、溶槽、岩溶洼地、岩溶水降落漏斗、落水洞、岩溶竖井等，水平岩溶为溶洞。不同的岩溶形态，其异常形态也各自不同。同时异常形态还受其上覆地层厚薄变化的影响而有所不同。工程选址所研究的岩溶深度一般不超过30 m，在此深度范围内，所遇到的岩溶现象大多数为垂直型岩溶，只有少部分是水平岩溶。而作为不良地质作用的岩土工程勘察，岩溶勘察主要的勘察对象是岩溶水降落漏斗、落水洞、岩溶竖井和溶洞。

一般来说，落水洞和岩溶竖井一般为直立或倾斜型带状低阻异常，而溶洞则为封闭的或是深部的低阻异常或高阻异常，但任何岩溶现象都不会孤立存在，而通常是相互关联的，如落水洞局部的水平溶蚀就会形成溶洞，从而引起异常的不同形态。

临沂市城区西南部岩溶塌陷多发生在岩溶地下水降落漏斗内，区内多为大面积隐伏奥陶系灰岩，上覆5~12 m的第四系松散层(多为黏土和砂层)，二元相结构明显，加之20世纪80年代后，对地下水的开采不规范，超采严重，导致地下水位在基岩面附近起伏变化（图9），局部产生渗透潜蚀和真空吸蚀，并且上覆第四系松散层伴随地下水流动，形成空洞，在受到地下水水流携带、机械振动和人类活动等多因素引发下^[31]，塌陷溶洞逐渐扩大，并可形成多孔相连现象，最终将因顶板失稳产生岩溶塌陷。课题组通过对自1991－2021年水位进行调查研究，调查结果发现，原先在齐家庄附近的超采地下水井已关闭，在2012年和2013年6月份时突发大雨，导致地下水位升高，2013年之后，地下水位一直保持在基岩面以下，岩溶发育条件并不充足，并未发生新的岩溶塌陷，但以往存在的溶洞仍保留至今，因此，城市建设时仍需详细勘察。

图  9  齐家庄村1991－2021年6月份地下水位曲线图

Figure  9.  Curves of groundwater levels in Qijiazhuang village from 1991 to June 2021.

下载: 全尺寸图片 幻灯片

6.   结　论

通过对临沂市城区西南部地下（岩溶）水降落漏斗进行分析，发现选址区位于岩溶发育区，在具备物性工作的前提下，对选址区开展了三维高密度电法工作，首先建立了区内含岩溶的地电模型；其次对地电模型进行了数值模拟，认为针对地层界面和岩溶的反演结果可信；最后将反演结果进行三维成图，在钻孔验证的基础上，认为三维高密度电法对地层界面和岩溶分布响应较好，为以后开展此类工作提供了经验。研究结果表明：

(1)在正反演数值模拟相互结合基础上开展分析，认为反演结果可信；

(2)三维高密度电法数据量大，反映的地电信息丰富，对第四系地层和灰岩界面的反映明显，对溶洞的显示直观，应用于工程选址效果较好；

(3)经分析，临沂市岩溶降落漏斗由2003年时的以梨杭—后岗头村为中心，在2007时转移至齐家庄村附近，期间岩溶塌陷频发，而2012年之后关闭了附近的超采岩溶地下水水井，并未发生过岩溶塌陷，但前期的岩溶发育将遗留部分溶洞，在城市建设时应做好勘察工作。

致谢：感谢编辑和匿名审稿专家的辛苦付出！