CN111737871A - 一种结合岩土材料特性的落石轨迹三维预测分析方法 - Google Patents

Abstract

本发明公开了一种结合岩土材料特性的落石轨迹三维预测分析方法，属于地质灾害防治领域与数值计算领域。该方法包括如下步骤：S1：将拟预测区域典型陡崖坡面几何形态转换为坐标数据，将坐标点排序后采用多项式对相邻点直线段进行描述并以张量表达；S2：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按岩土体的不同采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵；S3：将落石初始位置及状态参数装载为初始矩阵；S4：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；S5：在剖面预测概率模型基础上引入包含“横移比”的落石三维预测算法对落石在平面内的分布规律预测；S6：生成满足特定空间规律的防护结构设置区域。

Description

一种结合岩土材料特性的落石轨迹三维预测分析方法

技术领域

本发明属于地质灾害防治工程与数值计算领域，具体涉及一种结合岩土材料特性的落石 轨迹三维预测分析方法。

背景技术

崩塌落石在重庆市乃至我国山区分布广泛、暴发频繁、危害严重，是山区交通、城建、 地质环境保护领域三大地质灾害类型之一。

为工程防治崩塌落石灾害，常常采用拦石网、拦石墙、拦石堤、拦石栅栏、棚洞等被动 防护系统拦截、遮蔽落石体，通过阻止崩塌落石冲击居民区、公路、铁路等被保护对象而达 到减灾防灾的目的。

由于崩塌落石灾害威胁范围即为其运动所达区域，威胁能力来自于其冲击毁损能力，所 以，在研究与工程实践中，预测和计算其可能的运动路径、冲击动能、威胁区域便成了合理 确定被动防治方案，乃至防治成败的关键。国内外技术文献中均规定或建议：拦石网等被动 防治系统布置的位置、走向和范围均需要依据崩塌落石运动计算的结果优选，应优选弹跳高 度小、运动动能小的区段布设，以充分发挥被动防治系统的拦截性能，减少翻越和系统毁损 风险，进而提升防治效益。

但是，由于缺乏可靠的、便于工程人员使用的辅助设计软件，实际工程中按规范方法进 行运动和冲击计算，不仅繁琐和效率低下，而且所得运动路径、弹跳高度、运动动能等关键 设计参数，均为在确定性计算模型、确定性恢复系数、摩擦系数等参数设定下所得，难以做 到可重复、大规模和随机模拟，实际上，以确定性模型和参数来模拟崩塌落石运动随机过程， 并以所得单一运动路径和运动参数作为设计依据，其可靠性是有疑问的。甚至，由于计算繁 琐或不可靠，存在不做运动计算而仅凭经验进行被动防护系统布设的情形。以上问题，轻则 导致被动拦截系统的性能难以有效发挥，重则导致安全储备不足，以致出现落石飞越或拦截 系统被击毁等实例。

为了解决大规模重复计算和敏感参数随机取值的问题，从根本上实现落石运动路径随机 性的模拟，国外开发了Rocfall、STONE等计算软件。

此类软件优点在于简单易用，但也存在以下问题：一是只能做二维分析，可以确定被动 防护系统的位置，但无法确定沿陡崖走向方向的横向布设长度，也就是说无法预测落石三维 威胁范围和进行三维辅助设计。二是没有拦石墙、拦石堤、落石槽、拦石网等不同类型被动 防护系统，软件分析中难以反映不同类型系统的拦截特性。三是法向恢复系数、切向恢复系 数、滚动摩擦系数的取值及随机模型未考虑接触面岩土体自身物理力学参数指标的影响，某 些情形弹跳高度和范围远超实际崩落案例。

发明内容

针对上述现有技术中存在的不足之处，本发明提供了一种结合岩土材料特性的落石轨迹 三维预测分析方法。

一种结合岩土材料特性的落石轨迹三维预测分析方法，该方法包括以下步骤：

S1：将拟预测区域典型陡崖坡面几何形态转换为坐标数据，将坐标点排序后采用多项式 对相邻点直线段进行描述并以张量表达；

S2：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按岩土体的不同 采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵；

S3：将落石初始位置及状态参数装载为初始矩阵；

S4：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；

S5：在剖面预测概率模型基础上引入包含横移比的落石三维预测算法对落石在平面内的 分布规律预测；

S6：生成满足特定空间规律的防护结构设置区域。

进一步，所述步骤S1具体为：将拟预测区域典型陡崖坡面几何形态转换为坐标数据对 落石运动分析问题所需的地形及防护结构外形边界轮廓进行确定；

其中在建模中，应满足以下建模准则：

(1)模型边界：用于区分轨迹预测分析工作区与不参与信息拾取的非工作区；

数据的最低点高程的0.9倍高程作为底界；以图形界面绘制或导入数据的最高点高程的 1.1倍作为顶界；数据的最小平距的0.9倍作为左界；以图形界面绘制或导入数据的最大平 距的1.1倍作为右界；

(2)几何组成：用于绘制具体的建模信息；

①线段：将两个坐标点进行连线，以代表坡面或支护结构的非封闭区域；

②在陡崖及其下部斜坡面边界插值加密线段以提高计算精度；

陡崖及其下部斜坡面与防护结构形成含N点的集合采用N-1条多段线表示，将表达多段 线的多项式以矩阵形式表达。

进一步，所述步骤S2具体为：根据实体落石试验得到不同岩土体接触参数，对构成坡 面的多段线按岩土体的不同采用在材料赋值区间内随机赋值接触参数的方法形成材料参数 矩阵，对包含预设防护结构的模型，同样采用在材料赋值区间内随机赋值接触参数的方法形 成材料参数矩阵；

包括以下步骤：

①将实体落石试验得到的不同岩土体接触参数数据组成数据库，按地层岩性、地质构造 对岩土体的接触参数进行分类和整理，形成不同类别的材料矩阵；

②将陡崖及其下部斜坡轮廓按节点划分为多段线组连续曲线，对不同岩土体属性区域内 多段线按其与数据库对应的接触参数进行赋值，实现对落石接触区域岩土物理力学特性的准 确表达。

进一步，所述步骤S3具体为：预测计算模型中的坐标节点及形态，对初始状态赋值， 初始状态包括水平初始速度、垂直初始速度、水平初始加速度、垂直初始加速度、落石质量 和落石形态，提供落石拟迭代计算次数，形成初始状态矩阵；

包括以下步骤：

①通过数据生成落石位置初始矩阵，仅有坐标节点，其余初值为0；

封闭多线段：将坐标点视为几何形心，代表边界轮廓即表征落石外边界轮廓；

②根据需求确定落石形态：

提供两种形态：一种是质点，另一种是球形；根据不同选择装载形态部分参数；

③块体尺寸分布存在冲突时，默认优先级顺序为：1)落石位置；2)陡崖及其下部斜坡 面边界线段；3)防护结构边界；

④对落石初始状态进行赋值；将落石坐标节点信息、形态信息、拟迭代计算次数，形成 初始状态矩阵。

进一步，所述步骤S4具体为：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代 计算；

包括以下步骤：

①根据落石初始条件采用ODE轨迹算法对落石在第一次碰撞地面前的轨迹进行求解，ODE 轨迹算法指采用变步长4阶与5阶Runge-Kutta算法的非刚性常微分方程；

②碰撞时刻，采用改进GJK算法将拟接触地面多段线视作单纯形，保证明科夫斯基解的 存在，避免普通GJK算法可能出现的无解情况；在碰撞发生后，落石初始入射状态受接触材 料参数影响并以不同出射状态弹出，弹出状态作为第二次落石轨迹计算的初始状态进行下一 步计算；

③基于不同坡面岩土体接触参数随机区间赋值，对落石轨迹按迭代次数进行计算；将每 次迭代计算得到的落石轨迹按笛卡尔坐标系点集形式储存于预测矩阵。

进一步，所述步骤S5具体为：根据剖面范围内ODE轨迹算法和GJK改进算法得到的落 石轨迹矩阵为基础，引入包含横移比的落石三维预测算法对落石在三维空间内的分布规律预 测；

包括以下步骤：

①根据ODE轨迹算法和GJK改进算法得到的落石轨迹矩阵在剖面上的轨迹分布图，按高 斯过程对落石在剖面上落点分布特征进行分析；

②引入包含横移比的落石三维预测算法对落石在三维空间内的分布规律预测。

进一步，所述步骤S6具体为：生成满足特定空间规律的防护结构设置区域；

包括以下步骤：

①按S5步骤得到的落石在三维空间内的分布规律过滤出满足特定空间规律的防护结构 设置区域；

②将上述区域采用高亮显示，设置被动防护结构的位置和类型。

与现有技术相比，本发明的有益效果在于：

1)根据工程实际情况生成能准确表达落石所在区域陡崖及其下部斜边坡形态的典型剖 面，生成落石形态及其初始位置、初始物理力学参数；可按设计需求提供防护结构试算位置 及其结构形式；降低预测模型与实际工程落石形态和剖面形态的差异性；解决现有软件无法 实现危岩凹腔形态下的落石轨迹预测分析。

2)结合落石冲击碰撞实体模型试验，岩土体接触参数真实、可靠，计算过程中采用给 定赋值区间后在区间内随机取值进行迭代计算的方法，既体现了岩土体类别的不同，又对岩 土体物理力学特征的离散性进行了完整表达；从根本上解决落石轨迹预测中地面参数取值依 据不充分、经验参数取值的缺陷，为提高落石轨迹预测准确性的提高提供支撑。

3)落石预测分析模型采用ODE轨迹算法与改进GJK算法，弥补了现有商用软件碰撞过 程算法可能出现接触过程穿越接触面的缺陷。

4)引入包含“横移比”的落石三维预测算法，弥补了现有落石计算方法和程序无法对 防护结构进行计算与设置位置的缺陷。

5)本发明分析方法可根据预测结果满足防治要求的防护结构区域；并推荐设置被动防 护结构的位置和类型。

附图说明

图1为一种结合岩土材料特性的落石轨迹三维预测分析方法的流程图；

图2为某危岩落石及其所在陡崖及下部斜边坡建模示意图；

图3为落石轨迹示意图(单一轨迹计算)；

图4为落石轨迹示意图(多次轨迹计算)；

图5为含防护结构的落石轨迹示意图；

图6为落石轨迹预测平面分布图；

图7为运动体间的距离示意图；

图8为落石计算剖面轨迹分布图；

图9为落石威胁区域示意图。

具体实施方式

下面结合附图和具体实施方式对本发明作进一步详细地描述。

一种结合岩土材料特性的落石轨迹三维预测分析方法，如图1-6所示，包括以下步骤：

S1：将拟预测区域典型陡崖坡面几何形态转换为坐标数据，将坐标点排序后采用多项式 对相邻点直线段进行描述并以张量表达；

S2：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按岩土体的不同 采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵(对包含预设防护结构的 模型，同样采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵)；

S3：将落石初始位置及状态参数装载为初始矩阵；

S4：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；

S5：在剖面预测概率模型基础上引入包含“横移比”的落石三维预测算法对落石在平面 内的分布规律预测；

S6：生成满足特定空间规律的防护结构设置区域(包含防护结构的初始模型则根据区域 更新原模型中防护结构位置)。

步骤S1具体为：采将拟预测区域典型陡崖坡面几何形态转换为坐标数据，将坐标点排 序后采用多项式对相邻点直线段进行描述并以矩阵表达；

其中在建模中，应满足以下建模准则：

(1)数据的最低点高程的0.9倍高程作为底界；以图形界面绘制或导入数据的最高点 高程的1.1倍作为顶界；数据的最小平距的0.9倍作为左界；以图形界面绘制或导入数据的 最大平距的1.1倍作为右界。

(2)几何组成：用于绘制具体的建模信息；

①线段：将两个坐标点进行连线，以代表坡面或支护结构的非封闭区域；

②在陡崖及其下部斜坡面边界插值加密线段以提高计算精度；

陡崖及其下部斜坡面与防护结构形成含N点的集合采用N-1条多段线表示，将表达多段 线的多项式以矩阵形式表达。

步骤S2具体为：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按 岩土体的不同采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵(不同岩土 材料具有自身固有的接触参数赋值区间法向恢复系数、切向恢复系数、滚动摩擦系数等)， 对包含预设防护结构的模型，同样采用在材料赋值区间内随机赋值接触参数的方法形成材料 参数矩阵。

包括以下步骤：

将实体落石试验得到的不同岩土体接触参数数据组成数据库，按地层岩性、地质构造对 岩土体的接触参数进行分类和整理，形成不同类别的材料矩阵；

引入接触参数数据组成的材料矩阵M_N；

V＝ω×PG×M_N

M表示材料矩阵；e_t为切向材料恢复调整系数；e_n为法向材料恢复调整系数；e_ω为转动 材料恢复调整系数；

②将陡崖及其下部斜坡轮廓按节点划分为多段线组连续曲线，对不同岩土体属性区域内 多段线按其与数据库对应的接触参数进行赋值(按岩土类别赋值)，实现对落石接触区域岩 土物理力学特性的准确表达。

在冲击碰撞时刻，落石速度满足以下方程：

V＝ω×r＝ω×PG

PG＝(X_G-X_P)·i+(Y_G-Y_P)·j+0·k

ω+0·i+0·j-ω_z·k

V＝ω×PG×T_N

V表示碰撞时刻速度；ω为落石转动角速度；r为冲击点处的矢量；PG为冲击点处的矢 量；X_G为冲击初始水平距离；X_P为冲击回弹后水平距离；Y_G为冲击初始垂直距离；Y_P为冲击 回弹后垂直距离。

步骤S3具体为：预测计算模型中的坐标节点及形态，对初始状态赋值(水平初始速度、 垂直初始速度、水平初始加速度、垂直初始加速度，落石质量、落石形态)，形成初始状态 矩阵。

包括以下步骤：

①通过数据生成落石位置初始矩阵(仅有坐标节点，其余初值为0)；

封闭多线段：将坐标点视为几何形心，代表边界轮廓即表征落石外边界轮廓；

②根据需求确定落石形态(提供两种形态：1.质点，2.球形；根据不同选择装载形态部 分参数)；

③块体尺寸分布存在冲突时，默认优先级顺序为：1)落石位置；2)陡崖及其下部斜坡 面边界线段；3)防护结构边界；

对落石初始状态进行赋值(水平初始速度、垂直初始速度、水平初始加速度、垂直初始 加速度)；将落石坐标节点信息、形态信息、拟迭代计算次数，形成初始状态矩阵。

步骤S4具体为：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；

包括以下步骤：

①根据落石初始条件采用ODE轨迹算法对落石在第一次碰撞地面前的轨迹进行求解(ODE 轨迹算法指采用变步长4阶与5阶Runge-Kutta算法的非刚性常微分方程)；

公式中，y_k为k点处的解，t_k为求解时间段，k_i为对应时间段t_k的斜率；a_i为对应系数，h为k+1段的y值增量；

②碰撞时刻，采用改进GJK算法将拟接触地面多段线视作单纯形，保证明科夫斯基解的 存在。在碰撞发生后，落石初始入射状态受接触材料参数影响并以不同出射状态弹出，弹出 状态作为第二次落石轨迹计算的初始状态进行下一步计算。

假设两凸体A和B(落石分析中将视作多段线)，两凸体之间的距离为d(A,B)，如图7所示，则d(A,B)可用式(1)来表示：

d(A,B)＝min{||x-y||；x∈A,y∈B} (1)

GJK算法还可以返回2个物体之间距离最近的2点a和b，满足：

||a-b||＝d(A,B),a∈A,b∈B (2)

如果定义ν(C)为凸体集C中离原点最近的一个点，满足：

ν(C)∈C且||v(C)||＝min{||x||:x∈C} (3)

则A和B之间的距离可用Minkowski差表示为：

d(A,B)∈ν(C) (4)

式中：C＝A-B即为A与B的交集；A、B分别代表2个凸体边界，x为凸体边界A 的界限点坐标，y为凸体边界B的界限点坐标。

③基于不同坡面岩土体接触参数随机区间赋值，对落石轨迹按迭代次数进行计算；将每 次迭代计算得到的落石轨迹按笛卡尔坐标系点集形式储存于矩阵。

步骤S5具体为：根据剖面范围内ODE轨迹算法和GJK改进算法得到的落石轨迹矩阵为 基础，引入包含“横移比”的落石三维预测算法对落石在三维空间内的分布规律预测。

包括以下步骤：

①据ODE轨迹算法和GJK改进算法得到的落石轨迹矩阵在剖面上的轨迹分布图，按高斯 过程对落石在剖面上落点分布特征进行分析，如图8所示；

②引入包含“横移比”的落石三维预测算法对落石在三维空间内的分布规律预测。

落石运动路径计算是落石威胁区域确定比较理想的方法，路径计算的结果为二维运动路 径，可以反映落石在计算断面水平方向可以达到的最大距离，实际落石为三维运动，其威胁 区域为一个空间上的地理区域。如果能够结合前述试验得到的“横移比”指标，即可确定落 石威胁区域，如图9所示。图9中S为代表性落石计算得到的最大运动斜长，η为落石运动 偏移比，依据现场试验得到的分布规律矩阵，对于不同风险等级的工程，选取η值。

步骤S6具体为：生成满足特定空间规律的防护结构设置区域(包含防护结构的初始模 型则根据区域更新原模型中防护结构位置)；

包括以下步骤：

①按S5步骤得到的落石在三维空间内的分布规律过滤出满足特定空间规律的防护结构 设置区域；

②将上述区域采用高亮显示，推荐设置被动防护结构的位置和类型。

最后说明的是，以上实施例仅用以说明本发明的技术方案而非限制，尽管参照较佳实施 例对本发明进行了详细说明，本领域的普通技术人员应当理解，可以对本发明的技术方案进 行修改或者等同替换，而不脱离本发明技术方案的宗旨和范围，其均应涵盖在本发明的权利 要求范围当中。

Claims (7)

1.一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，该方法包括以下步骤：

S1：将拟预测区域典型陡崖坡面几何形态转换为坐标数据，将坐标点排序后采用多项式对相邻点直线段进行描述并以张量表达；

S2：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按岩土体的不同采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵；

S3：将落石初始位置及状态参数装载为初始矩阵；

S4：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；

S5：在剖面预测概率模型基础上引入包含横移比的落石三维预测算法对落石在平面内的分布规律预测；

S6：生成满足特定空间规律的防护结构设置区域。

2.根据权利要求1所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S1具体为：将拟预测区域典型陡崖坡面几何形态转换为坐标数据对落石运动分析问题所需的地形及防护结构外形边界轮廓进行确定；

其中在建模中，应满足以下建模准则：

(1)模型边界：用于区分轨迹预测分析工作区与不参与信息拾取的非工作区；

数据的最低点高程的0.9倍高程作为底界；以图形界面绘制或导入数据的最高点高程的1.1倍作为顶界；数据的最小平距的0.9倍作为左界；以图形界面绘制或导入数据的最大平距的1.1倍作为右界；

(2)几何组成：用于绘制具体的建模信息；

①线段：将两个坐标点进行连线，以代表坡面或支护结构的非封闭区域；

②在陡崖及其下部斜坡面边界插值加密线段以提高计算精度；

陡崖及其下部斜坡面与防护结构形成含N点的集合采用N-1条多段线表示，将表达多段线的多项式以矩阵形式表达。

3.根据权利要求2所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S2具体为：根据实体落石试验得到不同岩土体接触参数，对构成坡面的多段线按岩土体的不同采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵，对包含预设防护结构的模型，同样采用在材料赋值区间内随机赋值接触参数的方法形成材料参数矩阵；

包括以下步骤：

①将实体落石试验得到的不同岩土体接触参数数据组成数据库，按地层岩性、地质构造对岩土体的接触参数进行分类和整理，形成不同类别的材料矩阵；

②将陡崖及其下部斜坡轮廓按节点划分为多段线组连续曲线，对不同岩土体属性区域内多段线按其与数据库对应的接触参数进行赋值，实现对落石接触区域岩土物理力学特性的准确表达。

4.根据权利要求3所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S3具体为：预测计算模型中的坐标节点及形态，对初始状态赋值，初始状态包括水平初始速度、垂直初始速度、水平初始加速度、垂直初始加速度、落石质量和落石形态，提供落石拟迭代计算次数，形成初始状态矩阵；

包括以下步骤：

①通过数据生成落石位置初始矩阵，仅有坐标节点，其余初值为0；

封闭多线段：将坐标点视为几何形心，代表边界轮廓即表征落石外边界轮廓；

②根据需求确定落石形态：

提供两种形态：一种是质点，另一种是球形；根据不同选择装载形态部分参数；

③块体尺寸分布存在冲突时，默认优先级顺序为：1)落石位置；2)陡崖及其下部斜坡面边界线段；3)防护结构边界；

④对落石初始状态进行赋值；将落石坐标节点信息、形态信息、拟迭代计算次数，形成初始状态矩阵。

5.根据权利要求1所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S4具体为：基于改进GJK算法与ODE轨迹算法对落石运动轨迹迭代计算；

包括以下步骤：

①根据落石初始条件采用ODE轨迹算法对落石在第一次碰撞地面前的轨迹进行求解，ODE轨迹算法指采用变步长4阶与5阶Runge-Kutta算法的非刚性常微分方程；

②碰撞时刻，采用改进GJK算法将拟接触地面多段线视作单纯形，保证明科夫斯基解的存在，避免普通GJK算法可能出现的无解情况；在碰撞发生后，落石初始入射状态受接触材料参数影响并以不同出射状态弹出，弹出状态作为第二次落石轨迹计算的初始状态进行下一步计算；

③基于不同坡面岩土体接触参数随机区间赋值，对落石轨迹按迭代次数进行计算；将每次迭代计算得到的落石轨迹按笛卡尔坐标系点集形式储存于预测矩阵。

6.根据权利要求1所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S5具体为：根据剖面范围内ODE轨迹算法和GJK改进算法得到的落石轨迹矩阵为基础，引入包含横移比的落石三维预测算法对落石在三维空间内的分布规律预测；

包括以下步骤：

①根据ODE轨迹算法和GJK改进算法得到的落石轨迹矩阵在剖面上的轨迹分布图，按高斯过程对落石在剖面上落点分布特征进行分析；

②引入包含横移比的落石三维预测算法对落石在三维空间内的分布规律预测。

7.根据权利要求1所述的一种结合岩土材料特性的落石轨迹三维预测分析方法，其特征在于，所述步骤S6具体为：生成满足特定空间规律的防护结构设置区域；

包括以下步骤：

①按S5步骤得到的落石在三维空间内的分布规律过滤出满足特定空间规律的防护结构设置区域；

②将上述区域采用高亮显示，设置被动防护结构的位置和类型。

Images