滑坡预报与有限元实例分析对比研究 (1).pdf

收稿日期:2023-03-10;修订日期:2023-03-29作者简介:曹 伟(1981),男,硕士研究生,工程师,主要从事岩土工程设计及地质灾害研究。第 41 卷 第 4 期2023 年 8 月江 西 科 学JIANGXI SCIENCEVol.41 No.4Aug.2023 doi:10.13990/j.issn1001-3679.2023.04.032滑坡预报与有限元实例分析对比研究曹 伟,陈仲如(国家林业和草原局西南调查规划院,650031,昆明)摘要:为对滑坡的分析预报工作总结一套行之有效的方法,以三峡库区内滑坡为例展开研究。采用基于有限元原理的 ABAQUS 软件及摩尔库伦本构模型,建立三维地质模型并进行各工况下数值模拟计算。运用MATLAB 自编程序滑坡预测集成系统对滑坡既有监测数据进行分析与预测。最后,将数值模拟计算结果、实际监测数据及预测结果进行对比分析,结果显示该分析预报方法具有较好的可靠性和实用性。对滑坡稳定性及流固耦合分析具有积极作用,对滑坡预测及其他领域的预测预报研究也具有广泛的参考意义。关键词:涉水滑坡;监测预报;有限元模拟;稳定性分析中图分类号:P642 文献标识码:A 文章编号:1001-3679(2023)04-794-07Comparative Study of Landslide Forecasting and Finite Element Case AnalysisCAO Wei,CHEN Zhongru(Southwest Survey and Planning Institute of National Forestry and Grassland Administration,650031,Kunming,PRC)Abstract:In order to summarize a set of effective methods for the analysis and forecasting of land-slides,landslide in the Three Gorges Reservoir area is taken as an example.Based on the finite element theory of ABAQUS software and the Moll Coulomb constitutive model are used.The three-dimensional geological model is established,and the numerical simulation is carried out under vari-ous working conditions.Then,the landslide prediction integrated system programmed by MATLAB is used to analyze and forecast landslide.Finally,the numerical simulation results,actual moni-toring data and prediction results are compared and analyzed.The analysis and prediction method is proved to be reliable and practical by these conclusions.This paper plays an active role in landslide stability and fluid-structure coupling analysis.At the same time,this paper also has extensive refer-ence significance for landslide prediction and prediction research in other fields.Key words:fording landslide;monitoring and prediction;finite element simulation;stability analysis0 引言目前,如何对滑坡进行准确的预测预报是世界性难题。根据我国目前大多数已运营水库滑坡的发育情况来看,水库高水位运营期间,特别是蓄水期及降落期间极大概率发生大量的滑坡。三峡库区现处高水位运行期,极可能发生大量的滑坡,这些滑坡一旦下滑入江,除了对崩滑体上的居民造成危害之外,入江涌浪将产生派生灾难。因此,对复活型崩塌滑坡预测预报开展深入研究,是十分必要的,且具有重要的工程意义。1988 年,陈明东、王兰生1两位教授合作,首次在滑坡位移-时间曲线拟合中将灰色系统理论的 GM(1,1)模型纳入,并将拟合结果深入后推,进而提出了一个利用滤波灰色系统分析法模型来实现 对 滑 坡 的 中 期 预 报 的 方 法。廖 小 平2(1994)以弹塑性力学为理论基础提出了用于滑坡预报的功率模型。王尚庆3提出了一套完整的预报模型,包括安全系数、可靠性、变形速率和综合信息,且在三峡库区的实践中,更加系统地提出了一种基于滑坡变形破坏迹象、诱发因素等的宏观和综合预报方法。过去几十年,虽然我国在滑坡监测预报方面取得了显著进展,但大多数滑坡仍然发生在未被发现的地区,这一现象令人担忧4-5。本文以三峡库区内典型滑坡为例,筛选出适合该滑坡的预报模型,建立变形破坏的数学模型;综合分析数学模型和数值模拟试验的结果,并与实际变形破坏状况及滑坡现场监测资料验证分析,建立滑坡的预报模型和预报判据。1 预测模型及判据1.1 预测模型简介回归分析6就是通过对观察对象数据的统计分析和处理,研究认识和确定事物之间内在关系及形式的一种有效方法。按照自变量和因变量之间的关系类型,回归分析可分为线性回归分析和非线性回归分析。通过非线性回归技术,构建出滑坡位移-时间关系数学模型,从而更好地探索和预测滑坡未来的变形特征。1.2 预测判据滑坡时间预报和临界判据7主要系指滑坡出现快速变形滑动或滑坡开始时的征兆、状态及主要诱发滑坡因素,统计和研究分析提出滑坡宏观破裂长度值和宽度值以及滑坡边界长度之间的定量关系值作为判定滑坡临界变形滑动的依据。具体归纳概括出有下列几种时间预报的临界判据。1)滑坡的总体稳定性系数等于 1,仅局部小于 1;2)在滑坡前缘区域分布有带状的地下泉水涌出区;3)滑坡体裂隙带处于扩张变形状态,后缘裂隙带宽度应超过滑坡体长度的 1/2;4)滑坡蠕滑移动的速率持续增大,且变形显著,或滑动平均速率明显超过 10 cm/月,且受季节变化影响显著;5)雨季或单点暴雨,滑动面处饱和状态,且地表水冲蚀处出现明显崩塌、滑塌,或局部破裂变形加剧滑动、塌陷变形;6)地震动峰值加速度为 0.15 g,或滑坡场地地震烈度达到;7)滑坡周界轨迹连续不间断,轨迹分维值在1.401.60 区间。从以上滑坡的评判标准8和确定原则来看,要建立一种恒定的、又普遍适用的评判标准是不大可能。滑坡变形是其自身稳定状态最为直观的表征,加上位移监测信息相对容易获得,故位移变化信息也可用来监测。以与位移变形相关的数据信息作为判据包括以下几个因素:位移速率(mm/d)、最大位移量(mm)、位移加速度(mm/d2)、变形曲线切线角和位移矢量角变化。本文通过采用 MATLAB 自编应用程序作为位移数据信息的整理、预处理和预测分析工具,自动计算生成预测位移、位移速率、位移加速度和变形曲线切线角等位移判据信息的数据、图像文件。软件主界面如图 1 所示。清除图导入监测表格线性图柱形图导出表格和图片批量导出更新网格查看监测表格开始预测生成Word文档批量生成加密曲线回归分析模型时间位移图位移/mm02004006008001 000 1 200 1 400 1 600 1 800 2 000050100150200250观测期数图 1 预测分析集成系统主界面2 流固耦合求解流固耦合的方法一般主要有如下 2种9。1)解析法:由于岩土体材料在天然状态下很难运用解析法对流固耦合问题求得解析解。因此,流固耦合问题通常不采用该方法求解。2)数值分析:针对流固耦合问题的有限元数值分析方法一般主流方法包括有限差分法和有限元法。总的来说,从数学理论这一方向来看,有限差分法的近似化程度要比有限元法的高,但若从597第 4 期 曹 伟等:滑坡预报与有限元实例分析对比研究计算实际和应用层面来看,有限元法往往要更简便、灵活。因此,本文选取有限元法计算模拟。3 工程实例3.1 工程概述滑坡位于秭归县台子湾村,地处长江右岸,距离三峡大坝坝址约 53 km,地理坐标为:经度1103331.8,纬度 310031.1。滑坡位于百福坪背斜北侧,地层为三迭系中统巴东组泥岩、灰岩及砂岩,地层产状 2532,自然斜坡走向北西西,位于长江凹岸,滑坡全貌如图 2 所示。图 2 滑坡全貌图滑坡体前缘高程 6065 m,宽 570 m,后缘高程 560584 m,宽 440 m,纵向长 1 250 m,平面上总体呈长条形,面积约为 0.7 km2,估算方量为2 300104 m3,滑坡平面如图 3 所示。1ZG286ZG285ZG284ZG283ZG282ZG281长0100200 mN剖面及编号图例滑坡边界GPS监测点及编号江1图 3 滑坡平面图3.2 滑坡物质组成及结构特征1)滑体:据钻孔揭露,滑坡体物质平均厚度3035 m,滑坡中部最厚达71 m,滑体组成物质以残积和崩坡积物为主,滑体土层以块石土、碎石土为主,块石、碎石母岩成分以泥岩、灰岩及砂岩为主,滑坡地层如图 4 所示。125175225275325375425475高程/m174路公路长江175 m145 mZG286ZG285QdelZG284T2b3225路高程/m间距/m50100 150 200750700650600550500450400350300250900850800950 1000 1050 1100125175225275325375425475图例Qdel第四系滑坡堆积物碎石土T2b三叠系中统巴东组泥岩、砂质、灰岩滑面岩层产状GPS监测点编号3225ZG285图 4 滑坡剖面图2)滑床:高程在 440 m 以下地层为紫红色粉砂岩及泥岩,高程 440 m 以上地层为灰黄色泥灰岩。滑床总体形态为一斜坡地形,总体坡度为 2022,中间凹陷,侧缘突起,呈凹形。3)滑带:地表未见明显的滑带露头,上覆崩坡积物与下伏基岩接触带为滑带。4 库水作用下滑坡的稳定性分析4.1 计算模型滑坡计算模型采用有限元划分网格,共有126 000 个节点,116 127 个单元。模型单元类型采用减缩积分单元 C3D8RP,该类型单元可进行渗流-应力耦合,可以将所有节点孔隙压力作为一个自由度,保证了模型计算精度及速度。模型计算网格如图 5 所示。图 5 模型计算网格模型中考虑滑坡内的 3 种典型土体:滑体、滑带和滑床,各种材料均假定为各向同性体,综合对滑坡试验数据的统计、类比与反算分析,其抗剪强度参数取值及滑坡各相关参数取值见表 1。由于滑坡中的岩土体物理性质及地质条件较为复杂,其各项强度指标具有高度的非线性分布特征,因此无法采用线性本构模型对其进行模拟,而应采用非线性模型对其进行有限元模拟分析。经多次对比分析论证后,滑坡采用 Mohr-coulumb 模型来研究分析。697江 西 科 学2023 年第 41 卷表 1 滑坡典型土层相关参数典型土层渗透系数K/cms-1弹性模量E/Mpa泊松量 孔隙比 eC/KpaC/Kpa/r/KNm-3r/KNm-3滑体6.010-312.218.5 0.250.540.6614.826.517.518.416.218.912.516.520.122.522.523.4滑带2.210-62.713.40.200.450.6716.417.813.515.415.816.912.113.816.419.621.123.5滑床3.210-67.21040.180.11 1501 25060070029.535.126.130.323.225.326.428.4 计算模型应力分析的边界条件为:底部为三向约束,四周均采用法向约束;渗流分析的边界条件为:底部为不透水边界,四周和坡面均为自由透水面。4.2 计算工况滑坡位于三峡库区内,对滑坡稳定性产生影响的外部因素主要是库水位变化,因此本章主要对以下 3 种工况进行数值模拟,以此分析滑坡在库水位变化影响下的稳定性。工况 1:自重+地表荷载+初始水位(135 m 静水位);工况 2:自重+地表荷载+145 m 水位;工况 3:自重+地表荷载+175 m 水位。4.3 模型计算结果及分析以滑坡及滑坡的纵断面 1-1断面为研究对象,分别从应力场、孔隙水压力、位移及塑性区等方面论述 3 种工况下滑坡在库水位变化情况下应力场与渗流场的变化状态。通过软件计算结果表明如图 6 图 9 所示10,当库水位为初始状态(135 m)时,滑坡的塑性区范围很小;当水库水位从初始水位(135 m)上升到 145 m 时,应力场应力值呈增大趋势,塑性区等效塑性应变值增大的同时,其范围表现出从初始向局部性发展的趋势;当库水位继续上升到175 m 时,应力值继续增大,塑性区范围呈现出继续发展扩大的趋势,但总体来说范围并不大;当库水位从 175 m 降到 145 m 直至降到 135 m 时,渗流场变化显著,水力梯度明显增大;滑坡体位移增大显著,特别是滑坡前缘区域。从上述模拟结果可以得出,库水水位的急剧升降对滑坡的稳定性产生较大的不利影响。4.4 监测数据处理滑坡 6 个位移监测点 ZG281、ZG282、ZG283、ZG284、ZG285、ZG286,分别位于滑坡的滑体东上部、滑体东中部、滑体东下部、滑体西上部、滑体西中部、滑 体 西 下 部。受 篇 幅 限 制,本 节 只 将ZG281、ZG283、ZG285 的监测数据导入预测分析（a）初始水位最大主应力场分布图S,Mises(平均:75%)+1.032e+07+9.464e+06+8.604e+06+7.744e+06+6.885e+06+6.025e+06+5.166e+06+4.306e+06+3.446e+06+2.587e+06+1.727e+06+8.676e+05+7.959e+03（b）145 m水位最大主应力场分布图S,Mises(平均:75%)+1.025e+07+9.396e+06+8.543e+06+7.690e+06+6.837e+06+5.984e+06+5.131e+06+4.278e+06+3.425e+06+2.572e+06+1.719e+06+8.657e+05+1.268e+04（c）175 m水位最大主应力场分布图+1.012e+07+9.277e+06+8.434e+06+7.592e+06+6.750e+06+5.907e+06+5.065e+06+4.223e+06+3.380e+06+2.538e+06+1.696e+06+8.532e+05+1.090e+04S,Mises(平均:75%)图 6 各工况最大主应力场分布图+2.624e+06+2.406e+06+2.187e+06+1.968e+06+1.750e+06+1.531e+06+1.312e+06+1.094e+06+8.748e+05+6.561e+05+4.374e+05+2.187e+05+0.000e+00-4.001e+06?（a）初始水位孔隙水压力分布图POR(平均:75%)（b）145 m水位孔隙水压力分布图+3.000e+06+2.750e+06+2.500e+06+2.250e+06+2.000e+06+1.750e+06+1.500e+06+1.250e+06+1.000e+06+7.500e+05+5.000e+05+2.500e+05+0.000e+00-3.812e+06POR（c）175 m水位孔隙水压力分布图+3.000e+06+2.750e+06+2.500e+06+2.250e+06+2.000e+06+1.750e+06+1.500e+06+1.250e+06+1.000e+06+7.500e+05+5.000e+05+2.500e+05+0.000e+00-3.812e+06POR图 7 各工况孔隙水压力分布图797第 4 期 曹 伟等:滑坡预报与有限元实例分析对比研究+2.624e+06+2.406e+06+2.187e+06+1.968e+06+1.750e+06+1.531e+06+1.312e+06+1.094e+06+8.748e+05+6.561e+05+4.374e+05+2.187e+05+0.000e+00-4.001e+06?（a）初始水位位移图POR(平均:75%)（b）145 m水位位移图（c）175 m水位位移图+1,184e-01+1,086e-01+9,871e-02+8,884e-02+7,897e-02+6,910e-02+5,922e-02+4,935e-02+3,948e-02+2,961e-02+1,974e-02+9,871e-03+0,000e+00U,Magnitude+1,184e-01+1,086e-01+9,871e-02+8,884e-02+7,897e-02+6,910e-02+5,922e-02+4,935e-02+3,948e-02+2,961e-02+1,974e-02+9,871e-03+0,000e+00U,Magnitude+1,449e-01+1,328e-01+1,207e-01+1,087e-01+9,659e-02+8,452e-02+7,245e-02+6,037e-02+4,830e-02+3,622e-02+2,415e-02+1,207e-02+0,000e+00U,Magnitude图 8 各工况位移图（a）135 m水位塑性区分布图（b）145 m水位塑性区分布图（c）175 m水位塑性区分布图POR(平均:75%)+4,927e-03+4,516e-03+4,106e-03+3,695e-03+3,285e-03+2,874e-03+2,464e-03+2,053e-03+1,642e-03+1,232e-03+8,212e-04+4,106e-04+0,000e+00PEMAG(平均:75%)+4,927e-03+4,516e-03+4,106e-03+3,695e-03+3,285e-03+2,874e-03+2,464e-03+2,053e-03+1,642e-03+1,232e-03+8,212e-04+4,106e-04+0,000e+00PEMAG(平均:75%)+1.464e-02+1.342e-02+1.220e-02+1.098e-02+9.758e-03+8.538e-03+7.318e-03+6.099e-03+4.879e-03+3.659e-03+2.439e-03+1.220e-03+0.000e+00PEMAG(平均:75%)图 9 各工况塑性区分布图集成系统进行分析。分别选取 ZG281、ZG283 和 ZG285 3 个监测点,位移计算结果如图 10。通过预测分析集成系统分析已有监测数据,自动生成含有预测位移、位移加速度、位移速率和变形曲线切线角判据的数据、曲线图像文件(本文只列出位移图像文件),通过图像来判断滑坡的未来走势,及对滑动时间进行预报。（a）ZG281时间-位移(b)ZG283时间-位移（c）ZG285时间-位移050100150200250位移/mm曲线回归分析模型时间-位移图0200400600800 1 0001 2001 4001 6001 8002 000观测期数020406080100位移/mm曲线回归分析模型时间-位移图120140160180200050100150200250位移/mm曲线回归分析模型时间-位移图-503003500200400600800 1 0001 2001 4001 6001 8002 000观测期数0200400600800 1 0001 2001 4001 6001 8002 000观测期数图 10 时间-位移图运用中期预报模型中的曲线回归分析模型对6 个监测点 ZG281、ZG282、ZG283、ZG284、ZG285、ZG286 的监测数据进行了位移预测,预测结果为:ZG281、ZG284 监测点附近坡体没有明显的变形,ZG282、ZG283、ZG285、ZG286 监测点附近坡体都处于初始变形阶段,部分预则结果如图 10 所示,各监测点位置如图 11 所示。通过对以上各种数ZG285ZG284ZG283ZG282ZG281ZG286图 11 监测点在三维计算模型中的位置897江 西 科 学2023 年第 41 卷据的研究分析,目前滑坡处于基本稳定欠稳定状态,但是后期在诸多不利因素下,特别是在水位和降雨的双重耦合作用下,可能发生整体失稳。要想进一步地进行准确预报则需要更详细的监测数据、降雨、水位信息及宏观现象的记录。目前从监测信息及现场实际情况来看,预测结果与实际较为吻合,进一步说明了利用曲线回归分析模型处理该滑坡的监测数据是合适的,对该滑坡进行时间预报也是合理的。4.5 Abaqus 计算结果与监测值对比本节以滑坡的地质特征和大量的监测数据为基础,利用 MATLAB 程序对监测数据进行预处理后预测滑坡后续趋势。在对滑坡进行预测的结果(图 12)中可以看出,该处滑坡依然处于相对稳定状态,这与利用 ABAQUS 软件模拟的结果是一致的,并且与实际情况相吻合。对应状态:初始状态145m 175m降到145m降到135m工况01234对应状态:初始状态145m 175m降到145m降到135m工况计算结果监测结果02468101214位移/cm计算结果监测结果01234工况012345678910111213工况工况计算结果监测结果01234工况0123456789101112（a）ZG281的计算结果与监测数据对比(b)ZG282的计算结果与监测数据对比(c)ZG283的计算结果与监测数据对比(d)ZG284的计算结果与监测数据对比(e)ZG285的计算结果与监测数据对比(f)ZG286的计算结果与监测数据对比计算结果监测结果01234012345678910111213位移/cm对应状态:初始状态145m 175m降到145m降到135m对应状态:初始状态145m 175m降到145m降到135m01234计算结果监测结果02468101214位移/cm位移/cm位移/cm计算结果监测结果01234012345678910111213位移/cm对应状态:初始状态145m 175m降到145m降到135m对应状态:初始状态145m 175m降到145m降到135m图 12 计算结果与监测数据对比5 结论滑坡地质灾害预测和预报研究是近年来一项极具理论研究价值和工程应用前景的重要学术课题。本文采用理论分析与数值模拟相结合的方法,(下转第 807 页)997第 4 期 曹 伟等:滑坡预报与有限元实例分析对比研究京市游客与居民签到差异研究J.地理与地理信息科学,2022,38(1):37-44.18 谢永俊,彭霞,黄舟,等.基于微博数据的北京市热点区域意象感知J.地理科学进展,2017,36(9):1099-1110.19 WANG F.Quantitative Methods and Applications in GISM.CRC Press,2006.20 SILVERMAN B W.Density Estimation for Statistics and Data AnalysisM.Routledge,2018.21 张家旗,刘晏男,宋斌玢.基于 POI 数据的郑州市主城区生活服务业空间分布特征研究J.世界地理研究,2022,31(2):399-409.22 郝黎仁.SPSS 实用统计分析M.北京:中国水利水电出版社,2002.23 陈蔚珊,柳林,梁育填.基于 POI 数据的广州零售商业中心热点识别与业态集聚特征分析J.地理研究,2016,35(4):703-716.(上接第 799 页)详细研究了在应力场-渗流场耦合作用下,库岸滑坡在库水位升降条件下的稳定性,并对其后续的变形趋势做了预测,最后将该方法应用到滑坡的稳定性评价和监测预报中。通过理论研究和实例分析,得到的成果和结论主要如下。1)采用数值模拟对滑坡进行稳定性分析,第一步要明确其地质和库水位情况。合理确定地质模型和参数的选择是滑坡稳定性分析的关键。影响库岸滑坡稳定性的因素主要有地质情况和库水位升降。2)库岸滑坡产生的主要因素包括蓄水初期库水位上升和水位骤降。水的渗入使库岸滑坡的摩擦角和内聚力都有一定程度的衰减,同时会增加滑体的质量,进而增大了滑动面上的滑动力。库水的软化作用降低了滑动面土体抗剪强度参数,地下水的入渗也会对土体产生动水压力和静水压力,对滑坡体的稳定性产生不利的影响。3)本文以滑坡的地质结构特征和大量的监测数据为基础,利用 MATLAB 程序,在对滑坡进行预测的结果中可以看出,该滑坡依然处于相对稳定状态,这与利用 ABAQUS 软件模拟的结果是一致的,且与实际情况相吻合。由此说明该预报方法具有较好的可靠性和实用性。参考文献:1 陈明东,王兰生,边坡变形破坏的灰色预报方法C.全国第三次工程地质大会论文选集(下),成都科技大学出版社,1988:1226-1232.2廖小平.滑坡破坏时间预报新理论探讨J.地质灾害与环境保护,1994,5(3):25-29.3王尚庆.长江三峡滑坡监测预报M.北京:地质出版社,2008.4董建辉,唐然,李大鑫,等.四川汉源县二道坪滑坡成因及稳定性分析J.科学技术与工程,2021,21(22):9258-9263.5李效萌,李甜,刘金辉,等.赣东地区小型牵引式残坡积土质滑坡稳定性及破坏模式分析J.科学技术与工程,2021,21(22):9236-9242.6魏小楠.曲线回归分析模型在滑坡预测中的应用J.华东公路,2008(3):65-67.7李天斌.岩质工程高边坡稳定性及其控制的系统研究D.成都:成都理工大学,2002.8丁岩.三峡库区八字门滑坡预报判据研究D.西安:长安大学,2008.9曹伟.论库水位变化对土质滑坡稳定性的影响J.山西师范大学学报(自然科学版),2017,31(4):123-128.10 张欣.基于 ABAQUS 流固耦合理论的库岸滑坡稳定性分析D.济南:山东大学,2005.708第 4 期 董思宁等:基于社交媒体数据的“网红打卡地”空间布局及其相关因素分析