树坪滑坡变形的张力样条空间插值

钟传欣

（贵州省有色金属和核工业地质勘查局核资源地质调查院，贵州 贵阳 550002）

空间插值是基于已知采样点的值模拟估计未知点或表面值的过程，常用于将离散数据转换为连续的数据曲面[1，2]。空间插值广泛应用于气象领域[3-5]，能基于有限且离散的气象站监测数据获得一定区域的空间连续气象预测值。针对滑（边）坡变形空间插值应用，王建民等[4]基于时空Kriging法根据67个监测点数据对某露天坑边坡变形场进行插值计算，分析了边坡整体稳定性。Kriging法是根据变量空间自相似性插值，仅在数据量较大且分布均匀区域插值精度较高，并不十分适用于树坪滑坡变形空间插值。余礼仁等[5]提出了在D-InSAR基础上融合三次样条插值法用于煤矿采动变形监测，取得较好效果，适用于已获取高分辨遥感影像研究区。

不难看出，空间插值法应用于区域变形的研究仍处于起步阶段。由于张力样条插值法为分段逼近真值，插值结果数据面光滑[6]，较好地适用于根据少量反映整体变形趋势的监测数据插值生成滑体区域变形。因此，本文将张力样条插值法引入树坪滑坡(主滑区)变形场研究。

1 树坪滑坡概况

树坪滑坡[7]属于古崩滑堆积体，位于湖北省秭归县沙镇溪镇树坪村一组，地理坐标：110°37′0″E，30°59′37″N。滑坡发育于由三叠系中统巴东组泥岩、粉砂岩夹泥灰岩组成的逆向斜坡地段，地层产状120°～173°∠9°～38°。本文以树坪滑坡主滑区为研究对象，滑坡整体成“长舌”状，后缘高程350m～400m，前缘沉入长江（剪出口高程60m），东侧以山脊为界，东西两侧边界到中间冲沟距离相近。滑体规模巨大，面积约22.34m2×104m2，体积约11504m2×104m3。

滑坡前缘水位每年在145m至175m高程之间升降，树坪滑坡作为众多典型涉水滑坡之一，在30m水位变幅条件下变形持续剧烈增长。2003年6月起，开始对滑坡实施专业监测工程，先后在滑体上共设置了11个GPS地表变形监测点，目前主滑体上仅7个监测点能正常获取数据（如图1）。

图1 树坪滑坡主滑区GPS地表变形监测点（截至2017年累积变形）平面布置图

2 滑坡变形空间插值

2.1 张力样条插值数学原理

张力样条（Tension Spline）插值函数由Schweikert首次提出，为了消除三次样条插值函数有时会出现多余的拐点，能有效兼顾拟合曲面精确通过已知点和曲面光滑，且不出现交叉。

该插值函数基本构思是分段插值函数为直线插值和两个双曲函数shσx和chσx的线性组合[8-10]：

式中，σ为张力系数；式(1)整体具有二阶连续导数，σ的作用是控制拐点的位置和曲线的形态[9]。通过有关复杂的数学推导，可以得到张力样条函数的表达式：

式中，xi≤x≤xi+1，(i=0,1,2,…,n-1)，hi=xi+1-xi

2.2 滑坡变形特征

根据多年监测数据显示[9]，树坪滑坡累积位移曲线呈“阶跃型”增长，即每年长江水位下降至145m高程及低水位运行期间变形陡然增加，而长江水位上升至175m高程及高水位运行期间变形增长缓慢。鉴于树坪滑坡规模巨大且变形加速，2013年对其实施应急治理工程之后，滑坡体累积位移增速明显放缓趋于稳定（如图2）。

图2 树坪滑坡2007～2017年累积变形监测曲线

2.3 空间插值

图3 树坪滑坡变形张力样条空间插值图

树坪滑坡地表位移专业监测点（GPS）相对均匀布置于175m库水位以上滑体。本文截取7个专业监测点2007至2017年十年间累积位移监测数据，基于ArcGIS平台空间分析模块中张力样条（Tension Spline）插值工具，分析编制了树坪滑坡变形空间插值图（如图3）。

滑坡体纵向变形呈后部向中前部逐渐增大，横向变形呈两侧向中间增大。树坪滑坡（主滑区）呈“长舌”状，中轴线附近发育一条由“碟”形渐变至“V”形的拉长反“S”型冲沟。根据滑坡微地形，滑体可看作由3个子变形区组成，含冲沟主体及左右“两肩”。该滑坡滑体厚度较大（约30m～70m），主要受前缘库水影响产生“牵引式”滑动，因此中前部变形明显大于后部；受地形约束，冲沟内覆盖层在滑体前缘下滑或降水汇流引起土体饱和情况下，容易沿沟壁岩土分界面产生下滑位移导致累积变形较大。

3 结果与分析

如图3所示，滑坡大变形区域占滑体面积绝对多数，滑体总面积约22.34m2×104m2，变形量大于3000mm区域面积约17.94m2×104m2，占比高达80.27%。滑体变形在空间上存在明显差异，这与滑体所处斜坡特性有关，其中坡度和分布高程是最常见的两种特性。

3.1 空间变形与坡度

斜坡坡度空间变异性是影响滑坡体变形差异的重要因素。树坪滑坡（主滑体）坡度0°～75°，基于ArcGIS被划分为5类（图4(a)），不同坡度区间各变形量区域面积统计如图4(c)所示。滑体变形量在3000mm之内区域对称分布于各个坡度区间；随着变形量逐渐增大，所在区域坡度相应增加。根据滑体变形量与坡度的相关性研究，可为滑坡削减坡度治理提供理论指导。

3.2 空间变形与高程

滑坡前缘高程70m，后缘高程400m。基于ArcGIS将树坪滑坡（主滑体）数字高程模型按等高距66m重分类为5级（如图4(b)）。通过统计滑坡变形插值图与高程等级区划图交集面积，求得不同变形量区域面积在各级高程范围内分布规律。如图4(c)所示，变形量大于5000mm滑体与高程总体呈负相关，即区域面积越大，高程越小，主要集中于268m高程以下，占滑体总面积的36.32%；变形量为4000mm～5000mm滑体相对均匀分布于334m高程以下；变形量小于4000mm滑体区域面积与高程总体呈正相关，即高程越高，区域面积越大，主要集中于268m高程以上。大变形滑体之所以集中于低高程区域主要受控于滑坡前缘水位介于145m～175m涨落。长江水位上涨，滑体前缘受浮托减重效应产生滑移；库水位消落，滑体前缘受渗透作用形成的动水压力效应再次产生较大滑移。滑体呈“长舌”状，前缘不能牵引滑体后部同步滑移，因此滑体小变形量区域面积主要集中于后部较大高程范围内。

3.3 变形机制分析

图4 树坪滑坡空间插值变形与坡度、高程相关性

滑坡前缘水位升降。根据树坪滑坡累积变形监测曲线[9]，每年4～9月，滑坡前缘长江水位下降或低水位运行阶段曲线斜率显著增大，即变形速率增大；而每年10月至次年3月，长江处于水位上升或高水位运行阶段，变形监测曲线则表现相对平缓，即滑体变形速率减小。说明长江水位下降对树坪滑坡变形影响更大，且滑坡变形滞后于长江水位下降，具有明显的滞后效应。树坪滑坡体由渗透系数较小的含砾粘性土组成，长江水位上升阶段，水体不能及时进入，导致坡外水位高于坡内，产生内向负压差，一定程度上有利于滑体稳定；长江水位下降阶段，坡体内水位不能及时排出，导致坡内水位高于坡外，产生外向负压差，坡内形成渗透水压力加速滑体变形。

大气降水。在降雨条件下，滑坡体从非饱和变为饱和状态，含水率和自重增加，基质吸力减小，宏观上表现为土体软化趋于流塑状态，抗剪强度大幅降低。此外，树坪滑坡中轴部发育一条冲沟，有利于集流汇流，增强了滑体中间的冲刷作用。因此，降水是树坪滑坡变形加剧的另一因素。

4 结论

（1）基于ArcGIS平台和变形监测数据，引入张力样条法插值分析了树坪滑坡变形场，表明通过空间插值方法能有效解决滑坡变形监测数据离散无法全面反映滑坡变形的问题。

（2）根据统计分析，坡度越大、高程越小，树坪滑坡大变形区域分布越广。

（3）基础监测数据较少，滑坡变形空间插值结果存在一定误差，将在后续研究中改进。