CN115935742A - 基于有限差分的非均质边坡地震位移概率分析方法及系统 - Google Patents

Abstract

本发明公开了基于有限差分的非均质边坡地震位移概率分析方法及系统，通过生成边坡滑动位移危险性曲线来评价重要工程边坡的抗震性能。本方法及系统首先确定边坡土体参数的统计量并选取同边坡场地条件相匹配的地震波记录；然后产生土体参数随机场样本，循环执行随机有限差分数值模拟计算边坡累积位移；接着利用地震动参数和位移数据构建边坡位移预测模型，根据模型标准差进行地震动参数的优选；最后在概率地震危险性分析框架内发展边坡地震滑动位移危险性曲线，估计目标年超越概率对应的位移值。本发明能有效表征土体参数空间变异性和地震荷载不确定性，并捕捉地震作用下边坡体的应力变形机制，为高重要性程度的工程边坡抗震设计提供一种有效途径。

Description

基于有限差分的非均质边坡地震位移概率分析方法及系统

技术领域

本发明属于灾害预测技术领域，涉及地震作用下边坡永久位移预测分析方法及系统，尤其涉及一种基于有限差分的非均质边坡地震位移概率分析方法及系统。

背景技术

地震诱发的滑坡灾害具有波及范围广、规模大及破坏性强的特点，对人民群众生命财产安全构成严重威胁。边坡永久位移是一种常用的边坡地震稳定性评价指标，所以边坡地震位移的有效预测具有重要的工程意义。相比于拟静力边坡稳定分析和Newmark滑块位移分析，考虑应力变形机制的数值模拟能够提供更准确的边坡性能估计，尤其适用于高重要性程度的边坡工程。边坡地震滑动位移危险性概率分析方法可考虑多种不确定性因素的影响，估算指定抗震设防标准(或重现期)对应的位移值，当前在国内外备受关注。然而，边坡地震位移预测通常由基于Newmark滑块法的半经验模型来实现，其预测准确度一般具有较大偏差，往往难以被重要工程项目的设计方所能接受。

近年来，一些研究将基于数值模拟的位移预测模型纳入到概率分析中，为重要边坡工程的抗震性能评价提供了一种更为准确的途径。然而，现有的基于数值模拟的位移概率分析方法通常假设同一土层内的土体参数空间分布均匀，无法反映土体空间变异性和地形效应对边坡动力响应与滑移模式的耦合影响。因此，有必要在基于数值模拟的边坡地震滑移危险性概率分析中合理表征土体空间变异性，发展一种非均质边坡地震位移概率分析方法及系统。

发明内容

本发明的目的在于针对现有方法的不足，提供一种概念明确、操作简便、准确度较高的基于有限差分的非均质边坡地震位移概率分析方法及系统。

为了解决上述技术问题，本发明采用如下技术方案：

基于有限差分的非均质边坡地震位移概率分析方法，包括以下步骤：

步骤1：根据现场勘察资料估计边坡土体参数的统计量，并选取同边坡场地条件相匹配的m条地震波记录；

步骤2：建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，利用中点法生成N组土体参数随机场样本X_i，i＝1,2,…,N；

步骤3：将随机场样本X_i赋给对应的边坡模型单元，在模型底部边界施加第k条地震波作为荷载，执行有限差分计算边坡坡面的最大累积位移D_i,k，循环N×m次执行本步骤，即i＝1,2,…,N，k＝1,2,…,m；

步骤4：选取不同类型地震动参数IM，计算每条地震波对应的地震动参数值IM_k，k＝1,2,…,m，建立边坡地震位移预测模型D＝f(IM)；

步骤5：对所有备选地震动参数对应的模型标准差σ_lnD进行排序，将σ_lnD最小的IM作为目标地震动参数IM_T，执行概率地震危险性分析得到IM_T的年平均率密度MRD，即IM值与频次间的关系；

步骤6：将年平均率密度MRD与位移预测模型D＝f(IM_T)进行卷积，生成边坡滑动位移危险性曲线，估计目标年超越概率对应的位移值。

进一步地，步骤1中，利用随机场理论描述土体参数空间变异性；

以向量X＝[X₁,X₂,…,X_n]指代分布在空间n个位置的某种土体参数，向量X中的各变量均服从同一边缘分布，均值和变异系数分别为μ_X和COV_X，变量间的相关系数ρ(X₁,X₂)可由指数型自相关函数得到：

式中，

为第1个与第2个随机场单元间的水平向间距；

为第1个与第2个随机场单元间的竖直向间距；δ_h和δ_v分别为水平向波动范围和竖直向波动范围，反映土体参数空间自相关性即非均匀性的程度，当现场勘察数据较为充足时，统计量μ_X、COV_X、δ_h和δ_v，可通过矩估计或极大似然估计直接得到，否则可参考类似工程的现场资料或文献建议来取值。

进一步地，步骤2中，根据现场勘察资料建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，计算随机场变量的自相关系数矩阵R：

进一步地，步骤2中，利用中点法生成N组土体参数随机场样本X_i，i＝1,2,…,N，通过乔列斯基分解，将相关系数矩阵R分解为下三角矩阵L，可表达为：

LL^T＝R

式中，T代表矩阵转置，因此，一组随机场样本X可模拟为：

X＝F^-1[Φ(LU)]

式中，U为维度为n×1的独立标准正态样本向量；Φ(·)为标准正态分布变量的累积分布函数；F^-1(·)为土体参数的边缘累积分布函数的反函数，对上式重复执行N次可得到所有随机场样本，用于后续的概率计算。

进一步地，步骤3中，利用黏性边界即静态边界模拟底部半无限的地基条件，将第k条地震波的加速度时程a_k(t)转换为应力时程τ_k(t)后再施加在黏性边界上，转换公式为：

式中，G_max和ρ分别为模型底部材料的初始剪切模量和密度；v_k(t)为地震波加速度积分后的速度时程，执行一次动力有限差分模拟，记录边坡坡面若干监测点的累积位移时程，将时程结束时最大的累积位移值作为D_i,k。

进一步地，步骤4的具体实现包括以下子步骤：

步骤4.1：计算得到m组地震动参数值IM_k，由于每组地震动参数值对应N组随机场样本，分别将每种IM_k值复制N次得到IM_i,k，k＝1,2,…,m；

步骤4.2：利用N×m组数据D_i,k和IM_i,k回归得到各类地震动参数对应的位移预测模型D＝f(IM)，选用的模型具体形式为：

式中，ln代表自然对数符号；

为边坡滑动位移的预测值，单位为cm；IM代表所选地震动参数中的任意一种；a₀和a₁为回归系数；ε为标准正态分布变量；σ_lnD为模型标准差，反映预测不确定性，σ_lnD越小则表明模型对边坡位移预测的有效性越强。

进一步地，步骤5中，MRD(z)可视为地震动参数年平均超越概率λ_IM(z)随参数值z的变化率，通过下列方式获得：

式中，λ₀为地震震源活动率；f_M(m)和f_R(r)分别为震级M和断层距R的概率密度函数；f_IM(z|m,r)为IM在M＝m和R＝r时的概率密度函数，可表示为：

式中，μ_lnIM和σ_lnIM由某种地震动参数预测方程GMPE得到，GMPE一般通过对大量实测地震动记录进行回归分析得到，其本质为μ_lnIM或σ_lnIM关于地震工况参数的函数。

进一步地，步骤6中，滑动位移危险性曲线描述了位移年平均超越率λ_D(d)随阈值d的变化关系，其中λ_D(d)可计算为：

λ_D(d)＝∫P[D>dIM_T＝z]MRD_IM(z)dz

式中，P[D>d|IM_T＝z]为IM_T＝z时位移D超过阈值d的概率，可由下式得到：

式中，Φ(·)为标准正态分布变量的累积分布函数，将不同d值代入上述公式计算对应的λ_D(d)，便可绘制边坡滑动位移危险性曲线，可用来估计任意目标超越概率水平(重现期)对应的边坡位移值。

本发明还提供一种基于有限差分的非均质边坡地震位移概率分析系统，其特征在于，包括以下模块：

模块1，估计边坡土体参数的统计量和选取同边坡场地条件相匹配的地震波记录；

模块2，构建边坡有限差分数值模拟模型并生成土体参数随机场样本；

模块3，自动化计算非均质边坡在一系列土体参数随机场样本和地震加速度时程情况下的地震滑动位移；

模块4，构建以地震动参数作为输入的边坡地震位移预测模型；

模块5，实现地震动参数比选和基于最优地震动参数的概率地震危险性分析；

模块6，生成边坡地震滑动位移危险性曲线以及估计目标年超越概率对应的位移值。

进一步地，所述模块4选用的边坡地震位移预测模型具体形式为：

式中，ln代表自然对数符号；

为边坡滑动位移的预测值，单位为cm；IM代表所选地震动参数中的任意一种；a₀和a₁为回归系数；ε为标准正态分布变量；σ_lnD为模型标准差，反映预测不确定性，σ_lnD越小则表明模型对边坡位移预测的有效性越强。

与现有技术相比，本申请具有如下有益效果：

(1)本发明提供一种非均质边坡地震位移概率分析方法，适用于不同的地震工况与边坡条件；

(2)本发明的方法基于有限差分数值模拟，相较于传统的Newmark滑块分析能更准确地捕捉地震作用下边坡的应力变形机制；

(3)本发明的方法实现了基于随机场理论的土体参数空间变异性刻画，以及基于概率地震危险性分析的地震荷载不确定性表征。

附图说明

为了更清楚地说明本发明实施例中的技术方案，下面将对实施例描述中所需要使用的附图作简单地介绍，显而易见地，下面描述中的附图仅仅是本发明的一些实施例，对于本领域普通技术人员来讲，在不付出创造性劳动的前提下，还可以根据这些附图获得其他的附图。

图1为本发明实施例的方法流程示意图；

图2为本发明实施例中的所选地震波谱加速度分布图；

图3为本发明实施例中的边坡模型及确定性拟静力边坡稳定分析结果图；

图4为本发明实施例中的摩擦角随机场典型实现云图；

图5为本发明实施例中的边坡位移与地震动参数相关性分布图；

图6为本发明实施例中的位移模型标准差对比图；

图7为本发明实施例中基于谱强度SI的地震危险性曲线；

图8为本发明实施例中谱强度SI的年平均率密度曲线；

图9为本发明实施例中边坡地震滑动位移危险性曲线。

具体实施方式

以下结合附图和实施例详细说明本发明的技术方案。

为使本发明实施例的目的、技术方案和优点更加清楚，下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例是本发明的一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动的前提下所获得的所有其他实施例，都属于本发明保护的范围。

下面将结合附图，对本发明进行详细说明，基于本发明的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

实施例：请见图1，本发明提供的基于有限差分的非均质边坡地震位移概率分析方法，包括以下步骤：

步骤1：根据现场勘察资料估计边坡土体参数的统计量(如均值、变异系数、波动范围)，并选取同边坡场地条件相匹配的m条地震波记录；

土体受到各种地质、环境和物理化学过程的综合作用，原位土体参数往往沿着竖直和水平方向变化，呈现出一种空间变异性。经过近几十年的发展，随机场理论已被认为是描述土体参数空间变异性最有效的工具。以向量X＝[X₁,X₂,…,X_n]指代分布在空间n个位置的某种土体参数(如摩擦角φ)，向量X中的各变量均服从同一边缘分布(如对数正态分布)，均值和变异系数分别为μ_X和COV_X，变量间的相关系数ρ(X₁,X₂)可由指数型自相关函数得到：

式中，

为第1个与第2个随机场单元间的水平向间距；

为第1个与第2个随机场单元间的竖直向间距；δ_h和δ_v分别为水平向波动范围和竖直向波动范围，反映土体参数空间自相关性(非均匀性)的程度。当现场勘察数据(如静力触探测试CPT)较为充足时，统计量μ_X、COV_X、δ_h和δ_v，可通过矩估计或极大似然估计直接得到，否则可参考类似工程的现场资料或文献建议来取值。

本实施例中，从全球地震动数据库NGA-West2中挑选地震波，并选取同边坡场地条件(如土质或岩质场地)相匹配的m条地震波记录。选取原则为涵盖较广范围的地震震级和断层距，从而尽可能考虑地震动的不同特性(如幅度和频率量)，使后续得到的边坡位移预测模型可较好地应用于不同地震情景。

步骤2：建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，利用中点法生成N组土体参数随机场样本X_i(i＝1,2,…,N)；

步骤2.1：根据现场勘察资料建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，计算随机场变量的自相关系数矩阵R：

本实施例中，采用有限差分软件FLAC进行边坡动力分析，土体的循环应力应变特性由FLAC中结合滞后阻尼模型Sig4的摩尔-库伦(Mohr-Coulomb)塑性准则所描述，考虑的不确定性土体参数为黏聚力(c)、摩擦角(φ)和初始剪切模量(G_max)。该准则可解释动力响应弹性阶段的土体非线性行为，也可模拟潜在的不可逆应变累积，其中模型Sig4的参数由文献中的非线性土体动力特性模型校准。为了最大限度地减小模型边界处波的反射对动力响应分析的影响，在模型左右两侧施加自由场边界条件，并将侧向边界设置在离坡面足够远的位置。

步骤2.2：利用中点法生成N组土体参数随机场样本X_i(i＝1,2,…,N)。通过乔列斯基分解，将相关矩阵R分解为下三角矩阵L，可表达为：

LL^T＝R

式中，T代表矩阵转置。因此，一组随机场样本X可模拟为：

X＝F^-1[Φ(LU)]

式中，U为维度为n×1的独立标准正态样本向量；Φ(·)为标准正态分布变量的累积分布函数；F^-1(·)为土体参数的边缘累积分布函数的反函数。对上式重复执行N次可得到所有随机场样本，用于后续的概率计算。

步骤3：将随机场样本X_i赋给对应的边坡模型单元，在模型底部边界施加第k条地震波作为荷载，执行有限差分计算边坡坡面的最大累积位移D_i,k。

利用黏性边界(静态边界)模拟底部半无限的地基条件，将第k条地震波的加速度时程a_k(t)转换为应力时程τ_k(t)后再施加在黏性边界上，转换公式为：

式中，G_max和ρ分别为模型底部材料的初始剪切模量和密度；v_k(t)为地震波加速度积分后的速度时程。执行一次动力有限差分模拟，记录边坡坡面若干监测点的累积位移时程，将时程结束时最大的累积位移值作为D_i,k。循环N×m次执行步骤3，即i＝1,2,…,N，k＝1,2,…,m，共获得N×m个边坡位移值。

步骤4：选取不同类型地震动参数IM，计算每条地震波对应的IM_k值(k＝1,2,…,m)，建立边坡地震位移预测模型D＝f(IM)。

步骤4.1：以考虑三种地震动参数为例，包括峰值地面加速度(PGA)、峰值地面速度(PGV)和阿里亚斯强度(I_A)。需要自编程序或使用商业软件(如SeismoSignal)计算得到m组地震动参数值IM_k，即PGA_k、PGV_k和I_A,k(k＝1,2,…,m)。由于每组地震动参数值对应N组随机场样本，分别将每种IM_k值复制N次得到IM_i,k，即PGA_i,k、PGV_i,k和I_A,i,k(k＝1,2,…,m，i＝1,2,…,N)。

步骤4.2：利用N×m组数据D_i,k和IM_i,k回归得到各类地震动参数对应的位移预测模型D＝f(IM)。为了保证工程实践中的简易性，选用的模型具体形式为：

式中，ln代表自然对数符号；

为边坡滑动位移的预测值，单位为cm；IM代表所选地震动参数中的任意一种；a₀和a₁为回归系数；ε为标准正态分布变量(均值为0、标准差为1)；σ_lnD为模型标准差，反映预测不确定性，σ_lnD越小则表明模型对边坡位移预测的有效性越强。

步骤5：对所有备选地震动参数对应的σ_lnD进行排序，将σ_lnD最小的IM作为目标地震动参数IM_T，执行概率地震危险性分析(PSHA)得到IM_T的年平均率密度(MRD)，即IM值与频次间的关系。

PSHA可概括为4个基本步骤：划分潜在震源区、建立潜在震源区地震复发模型、描述地震动参数的概率分布、推求场址的地震危险性(如绘制地震动参数的年超越概率曲线)。经过数十年的理论研究和实践，PSHA已成为目前大多数抗震设计规范中地震动参数区划图编制的基础，如我国的《建筑抗震设计规范》(GB50011-2010)。MRD(z)可视为地震动参数年平均超越概率λ_IM(z)随参数值z的变化率，通过下列方式获得：

式中，λ₀为地震震源活动率；f_M(m)和f_R(r)分别为震级(M)和断层距(R)的概率密度函数；f_IM(z|m,r)为IM在M＝m和R＝r时的概率密度函数，可表示为：

式中，μ_lnIM和σ_lnIM由某种地震动参数预测方程(GMPE)得到。GMPE一般通过对大量实测地震动记录进行回归分析得到，其本质为μ_lnIM或σ_lnIM关于地震工况参数(如M和R)的函数。

步骤6：将步骤5得到的MRD与步骤4得到的位移预测模型D＝f(IM_T)进行卷积，生成边坡滑动位移危险性曲线，估计目标年超越概率对应的位移值。滑动位移危险性曲线描述了位移年平均超越率λ_D(d)随阈值d的变化关系，其中λ_D(d)可计算为：

λ_D(d)＝∫P[D>dIM_T＝z]MRD_IM(z)dz

式中，P[D>d|IM_T＝z]为IM_T＝z时位移D超过阈值d的概率，可由下式得到：

式中，Φ(·)为标准正态分布变量的累积分布函数。将不同d值代入上述公式计算对应的λ_D(d)，便可绘制边坡滑动位移危险性曲线，可用来估计任意目标超越概率水平(重现期)对应的边坡位移值。

接下来结合具体实施例对本发明做进一步的阐述。

一、边坡和震源模型简介

基于本发明所提方法，本实例将对一个简化边坡模型的地震滑动位移进行概率分析。边坡坡高和坡角分别为20m和30°，黏土层和基岩层厚度分别为40m和10m。黏土的黏聚力c、摩擦角φ、初始剪切模量G_max、重度γ的基准值分别取为55kPa、14°、180MPa和19.62kN/m³，基岩的G_max和γ的基准值分别取为1472MPa和22.56kN/m³。以一个点震源模型作为潜在震源，震源的活动率取为λ₀＝0.1，共考虑12个震级工况，范围为5～8(增量步为0.25)。假设边坡位于距震源10公里的岩石类场地。

二、具体实施流程

步骤1，将黏土的c、φ和G_max视作不确定性参数，根据勘察数据或工程资料估计边坡土体参数的统计量(如均值、变异系数、波动范围)，如表1所示。另外，从地震动数据库NGA-West2中选取83条地震波，对应的震级(M)范围为5.6～7.4，断层距(R)范围为0.9～83km，近地表30m的等效剪切波速(V_s30)大于600m/s。图2绘制了所选地震波的反应谱，其中峰值地面加速度(PGA)的范围为0.034g至1.406g。由图可知，不同地震动的特性(如幅度、周期)差异较大，可充分考虑地震动变异性。

表1边坡岩土体参数取值及统计量

步骤2，建立边坡有限差分数值模拟模型，将模型剖分成12066个材料单元，提取各单元的中心点坐标，利用Darendeli(2001)提出的半经验黏土模型标定Sig4的参数，从而表征土体剪切模量与阻尼比随剪应变的变化趋势。此外，为土体材料额外设置少量与刚度成比例的瑞利阻尼(0.2％)以去除高频噪声，为弹性基岩材料设置与质量和刚度成比例的瑞利阻尼(0.5％)。图3展示了边坡模型几何形状以及拟静力强度折减法所产生的的剪切带(最大剪切应变增量云图)，可知地震作用下边坡易呈深层滑动模式，滑动体体积较大。随后，利用中点法生成100组c、φ和G_max的随机场样本X_i(i＝1,2,…,100)。

步骤3，将随机场样本X_i赋给对应的边坡模型材料单元，在模型底部边界施加第k条地震波作为动力荷载。图4给出了边坡模型中φ的一组随机场实现，可以看到随机场模型能有效地描述土体参数的空间非均匀特性。沿边坡坡面布置一系列位移监测点，执行有限差分计算坡面的最大累积位移D_i,k(i＝1,2,…,100；k＝1,2,…,83)。

步骤4，选取8种常见的地震动参数IM，包括峰值地面加速度(PGA)、峰值地面速度(PGV)、阿里亚斯强度(I_A)、累积绝对速度(CAV)、平均周期(T_m)、谱强度(SI)、加速度谱强度(ASI)、70％能量持时(Ds_5-75)。计算每条地震波对应的IM_k值(k＝1,2,…,83)，利用位移和IM_k数据回归得到各类地震动参数对应的位移预测模型D＝f(IM)及其标准差σ_lnD，各模型的回归系数a₀和a₁如表2所示。地震动参数和位移的分布以及位移模型的拟合情况见图5。由图可知位移数据离散度较大，除了T_m和Ds_5-75以外，基于其它地震动参数的位移模型均可合理地预测边坡位移，但参数SI显然产生最佳的拟合效果。

表2位移预测模型的回归系数

步骤5，对所有备选地震动参数对应的σ_lnD进行排序，排序结果如图6所示，可见SI对应的σ_lnD最小，所以将SI作为目标地震动参数。考虑前述点震源模型，执行概率地震危险性分析PSHA得到SI对不同阈值的年平均超越概率曲线和年平均率密度MRD曲线，分别在图7和图8中给出。由图可知ln(SI)近似服从正态分布，与计算MRD中的假设一致。

步骤6，将MRD与位移预测模型lnD＝-5.47+2.05ln(SI)+0.65ε进行卷积，生成SI对应的边坡滑动位移危险性曲线。为了对比，同时推求了导致最差位移预测性能的地震动参数Ds_5-75对应的位移危险性曲线。图9展示了由最佳参数SI和最差参数Ds_5-75产生位移危险性曲线，可以看到选取Ds_5-75会显著高估边坡地震位移，进而导致过于保守和不经济的边坡抗震设计方案。若设计年限为475年，则可根据SI的边坡位移危险性曲线得到边坡在此年限内可能产生的最大位移为48cm，工程师需结合该边坡的具体抗震设计标准决定是否采取护坡加固等措施。

本实施例还提供一种基于有限差分的非均质边坡地震位移概率分析系统，其特征在于，包括以下模块：

模块1，估计边坡土体参数的统计量和选取同边坡场地条件相匹配的地震波记录；

模块2，构建边坡有限差分数值模拟模型并生成土体参数随机场样本；

模块3，自动化计算非均质边坡在一系列土体参数随机场样本和地震加速度时程情况下的地震滑动位移；

模块4，构建以地震动参数作为输入的边坡地震位移预测模型；

模块5，实现地震动参数比选和基于最优地震动参数的概率地震危险性分析；

模块6，生成边坡地震滑动位移危险性曲线以及估计目标年超越概率对应的位移值。

上述针对较佳实施例的描述较为详细，并不能因此而认为是对本发明专利保护范围的限制，本领域的普通技术人员在本发明的启示下，在不脱离本发明权利要求所保护的范围情况下，还可以做出替换或变形，均落入本发明的保护范围之内，本发明的请求保护范围应以所附权利要求为准。

Claims (10)

1.基于有限差分的非均质边坡地震位移概率分析方法，其特征在于，包括以下步骤：

步骤1：根据现场勘察资料估计边坡土体参数的统计量，并选取同边坡场地条件相匹配的m条地震波记录；

步骤2：建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，利用中点法生成N组土体参数随机场样本X_i，i＝1,2,…,N；

步骤3：将随机场样本X_i赋给对应的边坡模型单元，在模型底部边界施加第k条地震波作为荷载，执行有限差分计算边坡坡面的最大累积位移D_i,k，循环N×m次执行本步骤，即i＝1,2,…,N，k＝1,2,…,m；

步骤4：选取不同类型地震动参数IM，计算每条地震波对应的地震动参数值IM_k，k＝1,2,…,m，建立边坡地震位移预测模型D＝f(IM)；

步骤5：对所有备选地震动参数对应的模型标准差σ_lnD进行排序，将σ_lnD最小的IM作为目标地震动参数IM_T，执行概率地震危险性分析得到IM_T的年平均率密度MRD，即IM值与频次间的关系；

步骤6：将年平均率密度MRD与位移预测模型D＝f(IM_T)进行卷积，生成边坡滑动位移危险性曲线，估计目标年超越概率对应的位移值。

2.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于：步骤1中，利用随机场理论描述土体参数空间变异性；

以向量X＝[X₁,X₂,…,X_n]指代分布在空间n个位置的某种土体参数，向量X中的各变量均服从同一边缘分布，均值和变异系数分别为μ_X和COV_X，变量间的相关系数ρ(X₁,X₂)可由指数型自相关函数得到：

式中，

为第1个与第2个随机场单元间的水平向间距；

为第1个与第2个随机场单元间的竖直向间距；δ_h和δ_v分别为水平向波动范围和竖直向波动范围，反映土体参数空间自相关性即非均匀性的程度，当现场勘察数据较为充足时，统计量μ_X、COV_X、δ_h和δ_v，可通过矩估计或极大似然估计直接得到，否则可参考类似工程的现场资料或文献建议来取值。

3.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于：步骤2中，根据现场勘察资料建立边坡有限差分数值模拟模型，将模型剖分成n个材料单元，提取各单元的中心点坐标，计算随机场变量的自相关系数矩阵R：

4.根据权利要求3所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于：步骤2中，利用中点法生成N组土体参数随机场样本X_i，i＝1,2,…,N，通过乔列斯基分解，将相关系数矩阵R分解为下三角矩阵L，可表达为：

LL^T＝R

式中，T代表矩阵转置，因此，一组随机场样本X可模拟为：

X＝F^-1[Φ(LU)]

式中，U为维度为n×1的独立标准正态样本向量；Φ(·)为标准正态分布变量的累积分布函数；F^-1(·)为土体参数的边缘累积分布函数的反函数，对上式重复执行N次可得到所有随机场样本，用于后续的概率计算。

5.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于，步骤3中，利用黏性边界即静态边界模拟底部半无限的地基条件，将第k条地震波的加速度时程a_k(t)转换为应力时程τ_k(t)后再施加在黏性边界上，转换公式为：

式中，G_max和ρ分别为模型底部材料的初始剪切模量和密度；v_k(t)为地震波加速度积分后的速度时程，执行一次动力有限差分模拟，记录边坡坡面若干监测点的累积位移时程，将时程结束时最大的累积位移值作为D_i,k。

6.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于，步骤4的具体实现包括以下子步骤：

步骤4.1：计算得到m组地震动参数值IM_k，由于每组地震动参数值对应N组随机场样本，分别将每种IM_k值复制N次得到IM_i,k，k＝1,2,…,m，i＝1,2,…,N；

步骤4.2：利用N×m组数据D_i,k和IM_i,k回归得到各类地震动参数对应的位移预测模型D＝f(IM)，选用的模型具体形式为：

式中，ln代表自然对数符号；

为边坡滑动位移的预测值，单位为cm；IM代表所选地震动参数中的任意一种；a₀和a₁为回归系数；ε为标准正态分布变量；σ_lnD为模型标准差，反映预测不确定性，σ_lnD越小则表明模型对边坡位移预测的有效性越强。

7.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于，步骤5中，MRD(z)可视为地震动参数年平均超越概率λ_IM(z)随参数值z的变化率，通过下列方式获得：

式中，λ₀为地震震源活动率；f_M(m)和f_R(r)分别为震级M和断层距R的概率密度函数；f_IM(z|m,r)为IM在M＝m和R＝r时的概率密度函数，可表示为：

式中，μ_lnIM和σ_lnIM由某种地震动参数预测方程GMPE得到，GMPE一般通过对大量实测地震动记录进行回归分析得到，其本质为μ_lnIM或σ_lnIM关于地震工况参数的函数。

8.根据权利要求1所述的基于有限差分的非均质边坡地震位移概率分析方法，其特征在于，步骤6中，滑动位移危险性曲线描述了位移年平均超越率λ_D(d)随阈值d的变化关系，其中λ_D(d)可计算为：

λ_D(d)＝∫P[D>d|IM_T＝z]MRD_IM(z)dz

式中，P[D>d|IM_T＝z]为IM_T＝z时位移D超过阈值d的概率，可由下式得到：

式中，Φ(·)为标准正态分布变量的累积分布函数，将不同d值代入上述公式计算对应的λ_D(d)，便可绘制边坡滑动位移危险性曲线，可用来估计任意目标超越概率水平(重现期)对应的边坡位移值。

9.基于有限差分的非均质边坡地震位移概率分析系统，其特征在于，包括以下模块：

模块1，估计边坡土体参数的统计量和选取同边坡场地条件相匹配的地震波记录；

模块2，构建边坡有限差分数值模拟模型并生成土体参数随机场样本；

模块3，自动化计算非均质边坡在一系列土体参数随机场样本和地震加速度时程情况下的地震滑动位移；

模块4，构建以地震动参数作为输入的边坡地震位移预测模型；

模块5，实现地震动参数比选和基于最优地震动参数的概率地震危险性分析；

模块6，生成边坡地震滑动位移危险性曲线以及估计目标年超越概率对应的位移值。

10.根据权利要求9所述的基于有限差分的非均质边坡地震位移概率分析系统，其特征在于，所述模块4选用的边坡地震位移预测模型具体形式为：

式中，ln代表自然对数符号；

为边坡滑动位移的预测值，单位为cm；IM代表所选地震动参数中的任意一种；a₀和a₁为回归系数；

ε为标准正态分布变量；σ_lnD为模型标准差，反映预测不确定性，σ_lnD越小则表明模型对边坡位移预测的有效性越强。

Images