一种基于深部监测变形序列的边坡临滑预报方法、系统及介质

技术领域

本发明涉及大型水利水电工程建设期和运行期内的地质灾害监测与预警技术领域，具体涉及一种基于深部监测变形序列的边坡临滑预报方法、系统及介质。

背景技术

随着我国西南水电资源的持续开发，在高山峡谷区相继新建了众多的大型水利水电工程。坝址区和近坝库区内形成的大量高边坡通常规模巨大，其复杂的外部影响因素以及多变、不确定的地质条件决定了这类边坡失稳具有巨大的影响性和破坏性，严重威胁着近坝区的人民生命财产安全和大坝长期稳定运行。因此，坝址区及近坝库区内高边坡的稳定评价及临滑预警工作是水电工程建设和运行中的关键核心问题。

随着新型边坡监测手段和无线传感技术的蓬勃发展，实现了对不稳定边坡及潜在滑坡的全天候和“天-地-空”多层次监测，可为边坡变形预测和临滑预警提供重要的技术支持。鉴于边坡变形演化的不确定性和诱发机理的复杂性，仅单一分析边坡表观变形与环境影响因素的相关性仍无法实现边坡变形的中短期精确预测，且确定的失稳判据和临滑时间预报具有唯一不确定性和不具普适性的问题。目前已有的边坡失稳预报预警研究多集中在表观测点的变形量预测且预测模型无法直接输出具体滑坡时间，此外，当边坡经历加速变形阶段时，其变形演化特征与历史监测数据中的变形规律显著不同，由于监测信息和模型的局限性，训练的位移预测模型无法有效预测这种位移突变特征。因此，现阶段暂未有具备明确物理意义的高边坡全过程位移预测和临滑失稳预报方法。

发明内容

针对上述问题，本发明提供一种基于深部监测变形序列的边坡临滑预报方法、系统及介质，以解决如何基于岩土体蠕变理论考虑边坡变形的演化过程特征，以便结合边坡变形演化阶段特征和速率倒数法进行滑坡临滑失稳预警，为水利水电工程高边坡的变形演化分析和临滑预报系统构建提供有力支持。

为实现以上目的，本发明通过以下技术方案予以实现：

第一个方面，本发明提供一种基于深部监测变形序列的边坡临滑预报方法，包括如下步骤：

S1、通过深部变形测量装置采集边坡深部目标测点处的累计变形数据，利用Origin绘图软件绘制累计变形-时间(S-t)曲线，采用坐标变化对变形曲线的横-纵坐标进行同量纲化处理成T-t曲线；

S2、根据岩土体蠕变理论，引入分数阶导数元件和非线性粘塑性体改进西原蠕变模型，建立可反映边坡岩土体全阶段蠕变特性的本构模型，结合T-t曲线辨识模型参数；

S3、动态识别T-t曲线的加速变形起始点(AP点)，当边坡深部变形处于匀速变形阶段时，通过非线性蠕变本构模型进行变形量预测；当判定边坡深部变形进入加速变形阶段时，利用AP点后的监测变形序列，融合非线性蠕变模型和速率倒数法求得速率倒数趋势线与T-t曲线在时间轴上的交点，即边坡失稳时间t

S4、根据上述步骤S2及步骤S3的计算结果，分析边坡深部变形曲线的演化趋势、失稳预测时间t

在采用上述技术方案的同时，本发明还可以采用或者组合采用如下技术方案：

作为本发明的一种优选技术方案：步骤S1中，变形曲线无量纲化处理包括：根据S-t曲线特征初步识别并划分出初始变形阶段、等速变形阶段和加速变形阶段，利用等速变形阶段的平均变形速率v对S-t曲线纵坐标的量纲进行变换，得到纵、横坐标同量纲的T-t曲线的计算公式，如公式(1)所示：

式中，T(i)为量纲转换后的纵坐标值，ΔS(i)为某单位时间段内的变形变化量，

作为本发明的一种优选技术方案：步骤S2中，引入常用Riemann-Liouville(RL)型分数阶微积分理论，在恒定剪应力条件(σ(t)＝τ

式中，Γ(·)为伽玛函数，该函数是整数阶运算在实数和虚数范畴的推广，当Re为复数的实部时，如方程(3)所示：

利用分数阶导数元件和非线性粘塑性体改进西原模型的步骤包括：利用分数阶导数元件替换西原模型Kelvin体的牛顿黏壶，利用非线性黏性元件替换Bingham体中的牛顿黏壶。考虑元件的黏滞系数η(t)是与时间相关的非定常数，引入与时间呈指数衰减关系的η(t)黏滞系数后，得到非线性黏性元件的蠕变方程，如下公式(4)所示：

式中，τ为剪应力，τ

S

式中，S

式中，z和β均为函数方程中的变量。

作为本发明的一种优选技术方案：步骤S3中，准确识别T-t曲线的AP点的方法包括：计算短期变形速率平均线(SMA)与长期变形速率平均线(LMA)，利用SMA和LMA的“正向相交点”动态识别AP(t

式中，

作为本发明的一种优选技术方案：步骤S3中，边坡深部变形预测和失稳时间计算包括：当前边坡变形处于等速变形阶段时，利用式(5)拟合T-t曲线，采用最小二乘法和灰狼优化算法(GWO)全局优化辨识模型参数，GWO算法的优化目标函数如公式(9)所示：

式中，a、b、c和β为待辨识模型参数，t

若当前边坡变形处于AP点之后的加速变形阶段时，以AP点为原点建立新坐标系，并将AP点之后的变形序列作为模型输入进行变形曲线拟合，采用最小二乘法和GWO优化算法辨识模型参数，GWO算法的优化目标函数如下公式(10)所示：

式中，M和n为待辨识模型参数，n

在以AP点为原点的变形-时间坐标系下，为进一步简化速率倒数模型，假定速率倒数曲线的趋势为线性关系。利用线性函数拟合式(11)中的1/V～t′曲线，进一步求得线性趋势线与t轴的交点即为边坡失稳时间t

作为本发明的一种优选技术方案：步骤S4中，基于非线性蠕变模型计算得到未来一段时间的位移预测值和滑坡失稳时间，综合分析边坡深部变形趋势、临滑失稳时间t

第二个方面，本发明提供一种基于深部监测变形序列的边坡临滑预报系统，包括以下子模块：

数据采集模块，所述数据采集模块用于采集目标点的深部监测变形量并绘制累计深部变形-时间(S-t)曲线；

数据预处理模块，所述数据预处理模块用于动态识别深部变形S-t曲线的各变形阶段，计算等速变形阶段的平均速率

蠕变模型构建模块，所述蠕变模型构建模块用于引入分数阶导数元件和非线性粘塑性体至西原蠕变本构模型，建立边坡监测点处深部变形与时间的非线性关系；

边坡失稳时间预测模块，所述边坡失稳时间预测模块用于根据边坡的变形阶段预测各变形阶段的变形量，当边坡变形进入加速变形阶段后，利用AP点后的滑坡变形曲线，结合速率倒数法预测边坡失稳时间t

预警信号发布模块，所述预警信号发布模块用于综合分析变形趋势、失稳时间和宏观裂缝时空特征，及时发布对应的边坡预警等级和失稳时间。

还有一个方面，本发明提供一种包括计算机可读指令的计算机可读储存介质，所述计算机可读指令在被执行时使处理器执行如前文所述的基于深部监测变形序列的边坡临滑预报方法步骤。

本发明提供一种基于深部监测变形序列的边坡临滑预报方法、系统及介质，所述方法包括：由变形监测装置获取边坡深部监测点处的累积位移-时间(S-t)曲线，通过坐标变换对S-t曲线的横、纵坐标同量纲化处理，得到T-t曲线；引入岩土体蠕变理论，结合分数阶导数元件、非线性粘塑性体和西原蠕变模型，建立可反映边坡岩土体全阶段蠕变特性的非线性蠕变模型；当判定边坡岩土体处于等速变形阶段时，利用非线性蠕变模型进行位移预测；判定边坡岩土体进入加速变形阶段后，以加速变形起始点(AP)后的位移序列为输入，利用非线性蠕变模型进行变形预测，并结合速率倒数法预测边坡失稳时间t

本发明以边坡地质体的深部监测变形序列为基础，有效避免了表观变形序列中因人类工程活动和环境因素突变产生的无效数据，因此，边坡深部的变形监测数据相较于地表变形信息更具有预警分析价值。此外，本发明充分融合岩土体非线性蠕变理论和临滑分析方法，利用分数阶导数元件和非线性粘塑性体建立岩土体非线性蠕变关系，能够较好描述边坡岩土体的全阶段蠕变特性且模型参数物理意义明确。结合深部变形测量装置、非线性蠕变模型、速率倒数法和边坡宏观裂缝分布特征，实现边坡深部变形信息采集、演化过程分析和临滑预报需求一体化，可为水利水电工程高边坡的变形演化分析和临滑预报系统构建提供有力支持。

附图说明

图1为本发明所提供的基于深部监测变形序列的边坡临滑预报方法的流程图。

图2为某边坡典型监测点处的累计变形-时间(S-t)曲线。

图3为纵、横坐标同量纲转换后的变形-时间(T-t)曲线。

图4为西原蠕变本构模型结构示意图。

图5为非线性蠕变模型结构示意图。

图6为边坡进入加速变形阶段后利用非线性蠕变模型进行曲线拟合和模型参数辨识图。

图7为基于非线性蠕变模型的边坡位移预测及失稳预测时间关系图。

具体实施方式

为了更清楚地说明本发明的技术方案，下面结合附图和具体实施例，进一步阐明本发明的方法和系统，应理解这些实施例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，本领域技术人员对本发明的各种等价形式的修改均落于本申请所附权利要求所限定的范围。

众多学者研究指出边坡的变形演化和失稳特征由其内部岩土体的蠕变力学特性控制，其变形演化过程均包括初始变形、等速变形和加速变形三阶段，可视为滑坡预测预报的标准。因此，以边坡深部变形监测序列为基础，引入岩土体蠕变理论建立参数物理意义明确、概念清晰的边坡变形预测和失稳时间预报模型。请参阅图1，本实施例示意性地以某露天煤矿边坡为计算分析案例，采用10min/次监测周期对边坡地质体的变形进行持续监测，并选取某典型边坡位置处监测点的位移序列来构建位于预测模型和计算临滑失稳时间。

S1、根据该典型监测点的变形监测数据，采用Origin绘图软件绘制得到其累计变形-时间(S-t)曲线图，参见图2，结合S-t曲线特征将其划分为初始变形阶段、等速变形阶段和加速变形阶段。利用等速变形阶段的平均变形速率v对S-t曲线纵坐标的量纲进行变换，得到纵、横坐标同量纲的T-t曲线的计算公式，如下：

式中，T(i)为量纲转换后的纵坐标值，ΔS(i)为某单位时间段内的变形变化量，

利用计算公式(1)将边坡变形-时间(S-t)曲线转换为纵、横坐标同量纲的变形-时间(T-t)关系曲线，参见图3；

S2、西原模型是描述岩土体流变性质的经典元件模型，由Hooke体、Kelvin体和Bingham体串联而成，参见图4，能够较好地描述岩土体的线性非衰减蠕变行为。然而，由于西原模型中的基本元件均为理想线性元件，组合元件模型难以反映岩土材料的非线性加速蠕变行为。为了更好描述岩土体的全过程蠕变特性，在西原模型的基础上，分别引入分数阶导数软体元件和非线性黏性体替换Kelvin体和Bingham体中的牛顿黏壶，构建如下非线性蠕变模型，参见图5：

S

式中，S

根据转换后的T-t曲线分别计算短期变形速率平均线(SMA)与长期变形速率平均线(LMA)，再利用SMA和LMA的“正向相交点”动态识别AP(t

式中，

若边坡变形处于AP点之前的等速变形阶段时，拟合T-t曲线，采用最小二乘法和灰狼优化算法(GWO)优化辨识模型参数，GWO算法的优化目标函数如下公式所示：

式中，a、b、c和β为待辨识模型参数，t

若当前边坡变形处于AP点之后的加速变形阶段，以AP点为原点建立新坐标系，并将AP点之后的变形序列作为模型输入进行变形曲线拟合，参见图6，采用最小二乘法和GWO优化算法辨识模型参数，GWO算法的优化目标函数如下所示：

式中，M和n为待辨识模型参数，S

以边坡进入加速变形阶段为例阐述技术方案，将位移曲线进入加速变形(AP)后的位移序列(t

将辨识的模型参数M和n后带入上式中，对时间t求导后建立变形速率倒数-时间(1/V～t)关系曲线，如下所示：

以AP点为原点建立变形-时间坐标系，并利用线性函数拟合1/V～t曲线，进一步求得拟合线与t轴的交点即为边坡失稳时间t

在AP后的位移数据上应用GWO算法优化式进行模型参数辨识，随着应用越来越接近滑坡时间的监测数据对模型进行逐步更新和修正，得到位移预测值和失稳时间，参见图7，实施案例中模型参数和预测误差如表1所示。

表1模型辨识参数及预测滑坡时间统计表

边坡变形监测数据表明真实滑坡时间为94.23h，表1所示随着输入模型的监测数据逐渐趋近剧滑时刻，非线性蠕变模型的预测时间越来越接近真实滑坡时间，预测时间逐渐趋于稳定，在距离真实滑坡时刻前10h的预测时间误差率在5％以内。结合蠕变模型和速率倒数法最终计算的预测时间为94.33h，与实际滑坡时间基本一致，验证了本发明中基于非线性蠕变模型的滑坡短临预测的有效性。

本发明还提供一种基于深部监测变形序列的边坡临滑预报系统，包括以下子模块：

数据采集模块，用于自动采集边坡目标监测点处的深部变形量，并调用Origin绘制累计深部变形-时间(S-t)曲线；

数据预处理模块，用于自动识别监测变形-时间(S-t)曲线的各变形阶段，计算等速变形阶段的平均速率v，将S-t曲线同量纲化处理得到T-t曲线；

蠕变模型构建模块，用于引入分数阶导数元件和非线性粘塑性体到西原蠕变本构模型，建立边坡测点处的深部变形与时间非线性关系，结合实际边坡变形-时间(T-t)曲线，调用GWO优化算法来辨识非线性蠕变模型的模型参数；

边坡失稳时间预测模块，用于计算预测各变形阶段的变形量，对于在边坡变形处于AP点后(加速变形阶段)，结合速率倒数法计算预测边坡的失稳时间t

预警信号发布模块，用于综合分析变形演化趋势、失稳时间和宏观裂缝的时空特征，及时发布对应的边坡预警等级和失稳时间。

本发明还提供一种包括计算机可读指令的计算机可读储存介质，所述计算机可读指令在被执行时使处理器执行如前文所述的基于深部监测变形序列的边坡临滑预报方法步骤。

上述阐述的基于深部监测变形序列的边坡临滑预报系统、包括计算机可读指令的计算机可读储存介质与上述方法具有相同的发明构思，因此，此处未描述之处，可以参见前述方法实施例中的相关描述，此处不再赘述。本发明基于边坡深部岩土体的监测变形序列，充分融合岩土体非线性蠕变理论和边坡临滑预报方法，能够实现边坡深部变形信息自动化采集、演化过程动态分析、变形量预测及失稳时间预报需求一体化，可为水利水电工程坝区高边坡的变形演化分析和临滑预警预报系统构建提供有力支持，具有重要的工程意义和良好的应用前景。

上述具体实施方式用来解释说明本发明，仅为本发明的优选实施例，而不是对本发明进行限制，在本发明的精神和权利要求的保护范围内，对本发明做出的任何修改、等同替换、改进等，都落入本发明的保护范围。