一种水库诱发地震概率计算方法

摘要

本发明涉及一种水库诱发地震概率计算方法，属于水库诱发地震分析领域。本发明基于断层面的库仑应力，建立断层发震的功能函数，将功能函数中断层面凝聚力，断层面摩擦系数，发震部位孔隙水的实际水头，坐标面应力，坐标面正应力，断层面的方向余弦视为随机变量；再应用可靠度理论，计算水库诱发地震的概率。本发明将影响水库诱发地震的因素视为随机变量，分析其分布类型；再运用广泛采用的库仑应力理论分析断层面应力变化，并基于可靠度理论计算水库诱发地震的概率，为分析水库诱发地震提供新思路。

说明书

技术领域

本发明涉及一种水库诱发地震概率计算方法，属于水库诱发地震分析领域。

背景技术

水库诱发地震是由于人类拦河筑坝，在坝前壅高河水，形成水库引发的地震活动。世界 上首次有关水库诱发地震的资料报道是美国的米德湖(LakeMead，胡佛大坝的水库)。目前全 世界见诸报道的水库诱发地震震例为130余起，得到较普遍承认的约100起，仅占已建坝高 在15m以上大坝总数的2‰左右；中国是水库诱发地震较多的国家之一，迄今已报道的有34 例，得到广泛承认的为22例。按我国坝高大于15m的水坝约25800座计，发生诱发地震的 仅占1‰左右。世界上先后发生了4次震级大于6级的水库诱发地震，即中国的新丰江水库 诱发地震(6.1级，1962年3月)，赞比亚的卡里巴水库诱发地震(Kariba，6.1级，1963年9月)， 希腊的克瑞马斯塔水库诱发地震(Kremasta，6.3级，1966年)和印度的柯依纳诱发地震(Koyna， 6.5级，1967年)。尽管水库诱发地震的震级不高，发生概率不大，但是一旦发生会引起一些 次生地质灾害，给水利工程安全带来巨大的隐患。

水库诱发地震主要包括3种类型：构造型、喀斯特(岩溶)型、浅表微破裂型，其中又以 构造型的水库诱发地震强度较高，对水利工程的影响较大，也是世界各国研究最多的主要类 型。

目前水库诱发地震分析模型都是把有关参数视为确定性变量，但是由于影响水库诱发地 震的众多因素大多与地质构造条件有关，比如：深部岩体的力学指标，地震发生时地应力的 大小、方向，孔隙水的分布等，包含大量随机不确定性因素，导致这些因素很难准确测量。 即使采用先进的方法，对这些量进行准确测量，但是它们的值也是存在比较大的变化，给水 库诱发地震分析的带来困难。

发明内容

本发明提供了一种水库诱发地震概率计算方法，以用于分析水库诱发地震的概率。

本发明的技术方案是：一种水库诱发地震概率计算方法，基于断层面的库仑应力，建立 断层发震的功能函数，将功能函数中断层面凝聚力c，断层面摩擦系数μ，发震部位孔隙水的 实际水头h，坐标面应力σ_ij，坐标面正应力σ_ii，断层面的方向余弦n_i视为随机变量；再应用 可靠度理论，计算水库诱发地震的概率。

所述方法的具体步骤如下：

Step1、将断层面的剪切强度τ_crit与断层面剪应力τ相减，得到断层面库仑应力σ_f：

σ_f＝τ_crit-τ(1)

式中：断层面的剪切强τ_crit＝c-μ(σ_n+p)，断层面剪应力c为断层面 凝聚力(Mpa)，μ为断层面摩擦系数，p为断层面孔隙水压力(Mpa)，断层面正应力σ_n＝σ_ijn_in_j(Mpa)(i、j＝1,2,3代表x、y、z三个方向)，断层面应力(i、j＝1,2,3 代表x、y、z三个方向)；n_i、n_j为断层面的方向余弦，σ_ij为坐标面应力(Mpa)；其中，断层 面剪应力τ取断层滑动方向为正；断层面正应力σ_n取拉伸为正，压缩为负；当τ＞τ_crit时，断 层失稳；当τ＝τ_crit，断层处于临界状态；当τ＜τ_crit时，断层处于稳定；

Step2、断层面孔隙水压力p：

p＝γ_wh(2)

式中：γ_w为水的重度(KN/m³)，发震部位孔隙水的实际水头h＝ξ(d+H)；ξ为水头系数，取 值为0～1；d为发震部位埋深(m)；H为库水深度(m)；

Step3、根据断层面库仑应力，建立断层面滑动功能函数：

$g\left(X\right)=c-\mu ({\sigma }_{ij}{n}_i{n}_j+{\gamma }_{w}h)-\sqrt{{\left({\sigma }_{ij}{n}_j\right)}^2-{\left({\sigma }_{ij}{n}_i{n}_j\right)}^2},(i,j=1,2,3)---\left(3\right)$

若功能函数g(X)>0断层稳定；g(X)＝0断层处于临界状态；g(X)<0断层滑动；

Step4、诱发地震的概率计算：

将功能函数式(3)中断层面凝聚力c，断层面摩擦系数μ，发震部位孔隙水的实际水头 h，坐标面应力σ_ij，坐标面正应力σ_ii，断层面的方向余弦n_i视为随机变量；

功能函数g(X)对各变量的偏导数为：

$\frac{\partial g}{\partial c}=1---\left(4\right)$

$\frac{\partial g}{\partial \mu }=-({\sigma }_{ij}{n}_i{n}_j+{\gamma }_{w}h)---\left(5\right)$

$\frac{\partial g}{\partial h}=-{\mathrm{\mu \gamma }}_w---\left(6\right)$

$\frac{\partial g}{\partial {\sigma }_{ii}}=-{\mathrm{\mu n}}_i^2-{({\left({\sigma }_{kl}{n}_l\right)}^2-{\left({\sigma }_{kl}{n}_k{n}_l\right)}^2)}^{-\frac{1}{2}}\left(\right({\sigma }_{ij}{n}_j)n_i-({\sigma }_{kl}{n}_k{n}_l\left)n_i^2\right),(i,j,k,l=1,2,3)---\left(7\right)$

$\begin{array}{c}\frac{\partial g}{\partial n_i}=-2{\mathrm{\mu \sigma }}_{ij}{n}_j-{({\left({\sigma }_{kl}{n}_l\right)}^2-{\left({\sigma }_{kl}{n}_k{n}_l\right)}^2)}^{-\frac{1}{2}}\left(\right({\sigma }_{kl}{n}_l){\sigma }_{ik}-2({\sigma }_{kl}{n}_k{n}_l\left)\right({\sigma }_{ik}{n}_k\left)\right)\\ (i,j,k,l=1,2,3)\end{array}---\left(9\right)$

得到功能函数的梯度后，变量的方向余弦α_X为：

${\alpha }_X=\frac{-{\rho }_X{s}_X▿g_X}{{(▿g_X)}^T{C}_X▿g_X}---\left(10\right)$

式中：ρ_X为X的相关系数矩阵；s_X为X的标准差矩阵；C_X为X的协方差矩阵；

得到功能函数的方向余弦后α_X，通过迭代求解计算可靠指标β，迭代求解过程如下：

①首先选择初始验算点x*；

②计算功能函数梯度▽g_X；

③对于非正态变量则需通过JC法或映射法，转化为当量正态变量，然后通过下式计算；

$\beta =\frac{g_X\left(x^{*}\right)+{(▿g_X(x^{*}\left)\right)}^T({\mu }_X-x^{*})}{{({(▿g_X(x^{*}\left)\right)}^T{C}_X▿g_X(x^{*}\left)\right)}^{1/2}}---\left(11\right)$

式中：μ_X为随机变量的均值；

④利用式(12)计算新的验算点x*，然后返回第②步再计算，直至前后两次之差小于允 许误差；

x^*＝μ_X+βs_Xα_X(12)

得到β后，根据标准正态分布函数Ф求出诱发地震的概率p_f：p_f＝Φ(-β)。

本发明的有益效果是：本发明将影响水库诱发地震的因素视为随机变量，分析其分布类 型；再运用广泛采用的库仑应力理论分析断层面应力变化，并基于可靠度理论计算水库诱发 地震的概率，为分析水库诱发地震提供新思路。

附图说明

图1为本发明断层稳定的可靠指标与发震概率随孔隙水压力的变化。

具体实施方式

实施例1：如图1所示，一种水库诱发地震概率计算方法，基于断层面的库仑应力，建 立断层发震的功能函数，将功能函数中断层面凝聚力c，断层面摩擦系数μ，发震部位孔隙水 的实际水头h，坐标面应力σ_ij，坐标面正应力σ_ii，断层面的方向余弦n_i视为随机变量；再应 用可靠度理论，计算水库诱发地震的概率。

所述方法的具体步骤如下：

Step1、将断层面的剪切强度τ_crit与断层面剪应力τ相减，得到断层面库仑应力σ_f：

σ_f＝τ_crit-τ(1)

式中：断层面的剪切强τ_crit＝c-μ(σ_n+p)，断层面剪应力c为断层面凝聚力， μ为断层面摩擦系数，p为断层面孔隙水压力，断层面正应力σ_n＝σ_ijn_in_j(i、j＝1,2,3代表x、 y、z三个方向)，断层面应力(i、j＝1,2,3代表x、y、z三个方向)；n_i、n_j为断 层面的方向余弦，σ_ij为坐标面应力；其中，断层面剪应力τ取断层滑动方向为正；断层面正 应力σ_n取拉伸为正，压缩为负；当τ＞τ_crit时，断层失稳；当τ＝τ_crit，断层处于临界状态； 当τ＜τ_crit时，断层处于稳定；

Step2、断层面孔隙水压力p：

p＝γ_wh(2)

式中：γ_w为水的重度，发震部位孔隙水的实际水头h＝ξ(d+H)；ξ为水头系数，取值为0～1； d为发震部位埋深；H为库水深度；

Step3、根据断层面库仑应力，建立断层面滑动功能函数：

$g\left(X\right)=c-\mu ({\sigma }_{ij}{n}_i{n}_j+{\gamma }_{w}h)-\sqrt{{\left({\sigma }_{ij}{n}_j\right)}^2-{\left({\sigma }_{ij}{n}_i{n}_j\right)}^2},(i,j=1,2,3)---\left(3\right)$

若功能函数g(X)>0断层稳定；g(X)＝0断层处于临界状态；g(X)<0断层滑动；

Step4、诱发地震的概率计算：

将功能函数式(3)中断层面凝聚力c，断层面摩擦系数μ，发震部位孔隙水的实际水头 h，坐标面应力σ_ij，坐标面正应力σ_ii，断层面的方向余弦n_i视为随机变量；

功能函数g(X)对各变量的偏导数为：

$\frac{\partial g}{\partial c}=1---\left(4\right)$

$\frac{\partial g}{\partial \mu }=-({\sigma }_{ij}{n}_i{n}_j+{\gamma }_{w}h)---\left(5\right)$

$\frac{\partial g}{\partial h}=-{\mathrm{\mu \gamma }}_w---\left(6\right)$

$\frac{\partial g}{\partial {\sigma }_{ii}}=-{\mathrm{\mu n}}_i^2-{({\left({\sigma }_{kl}{n}_l\right)}^2-{\left({\sigma }_{kl}{n}_k{n}_l\right)}^2)}^{-\frac{1}{2}}\left(\right({\sigma }_{ij}{n}_j)n_i-({\sigma }_{kl}{n}_k{n}_l\left)n_i^2\right),(i,j,k,l=1,2,3)---\left(7\right)$

$\begin{array}{c}\frac{\partial g}{\partial n_i}=-2{\mathrm{\mu \sigma }}_{ij}{n}_j-{({\left({\sigma }_{kl}{n}_l\right)}^2-{\left({\sigma }_{kl}{n}_k{n}_l\right)}^2)}^{-\frac{1}{2}}\left(\right({\sigma }_{kl}{n}_l){\sigma }_{ik}-2({\sigma }_{kl}{n}_k{n}_l\left)\right({\sigma }_{ik}{n}_k\left)\right)\\ (i,j,k,l=1,2,3)\end{array}---\left(9\right)$

得到功能函数的梯度后，变量的方向余弦α_X为：

${\alpha }_X=\frac{-{\rho }_X{s}_X▿g_X}{{(▿g_X)}^T{C}_X▿g_X}---\left(10\right)$

式中：ρ_X为X的相关系数矩阵；s_X为X的标准差矩阵；C_X为X的协方差矩阵；

得到功能函数的方向余弦后α_X，通过迭代求解计算可靠指标β，迭代求解过程如下：

①首先选择初始验算点x*；

②计算功能函数梯度▽g_X；

③对于非正态变量则需通过JC法或映射法，转化为当量正态变量，然后通过下式计算；

$\beta =\frac{g_X\left(x^{*}\right)+{(▿g_X(x^{*}\left)\right)}^T({\mu }_X-x^{*})}{{({(▿g_X(x^{*}\left)\right)}^T{C}_X▿g_X(x^{*}\left)\right)}^{1/2}}---\left(11\right)$

式中：μ_X为随机变量的均值；

④利用式(12)计算新的验算点x*，然后返回第②步再计算，直至前后两次之差小于允 许误差；

x^*＝μ_X+βs_Xα_X(12)

得到β后，根据标准正态分布函数Ф求出诱发地震的概率p_f：p_f＝Φ(-β)。

实施例2：如图1所示，一种水库诱发地震概率计算方法，基于断层面的库仑应力，建 立断层发震的功能函数，将功能函数中断层面凝聚力c，断层面摩擦系数μ，发震部位孔隙水 的实际水头h，坐标面应力σ_ij，坐标面正应力σ_ii，断层面的方向余弦n_i视为随机变量；再应 用可靠度理论，计算水库诱发地震的概率。

实施例3：如图1所示，一种水库诱发地震概率计算方法，设某河道型水库，坝址区河 道蓄满水后水面宽5km，最大水深200m，坝址区岩性主要为花岗闪长岩，库底存在一断层， 断层产状为80°∠60°，岩体密度取2700kg/m³，弹性模量取30Gpa，泊松比取0.26，断层 面摩擦系数取1.2，凝聚力取600kpa，发震部位位于库底下5km，现分析水库蓄水诱发地震 的概率。

1、发震部位应力计算

(1)垂向应力计算

竖向自重应力近似采用上覆岩层重量σ₃₃＝σ_z＝γd＝2700×9.8×5×10³＝130.95Mpa；γ为 上覆岩层重量，d为发震部位埋深；

(2)水平向应力计算

水平向应力分布十分复杂，目前也提出很多计算方法，有采用广义虎克定理计算、海姆 (Heim)假说计算、采用实测资料的统计分析等。本文采用有关文献提出的统计规律计算

σ_x＝0.0216d+6.7808

σ_y＝0.0182d+2.2328

式中：d为发震部位埋深(m)；σ₁₁＝σ_x、σ₂₂＝σ_y分别为发震部位水平向大、小主应力(Mpa)。 取d＝5000m，得σ_x＝114.78Mpa，σ_y＝93.23Mpa。

(3)库水产生的附加应力

大量文献分析，库水对深部岩体产生的附加应力十分微弱。因此，这里不考虑库水产生 的附加应力。

2、随机变量分布类型和相关性

水库诱发地震主要涉及9个变量，由于缺少深部岩体这些参数的统计分析，因此本文假 设这些变量都服从正态分布，均值、标准差、变异系数如表1所示。断层面产状与其它变量 可视为相互独立，各变量间的相关性如表2所示。

表1各变量参数取值

变量 分布类型 均值 标准差 变异系数 c 正态 0.60Mpa 0.06 0.10 μ 正态 1.2 0.12 0.10 σ_x正态 114.78Mpa 22.956 0.20 σ_y正态 93.23Mpa 18.646 0.20 σ_z正态 130.95Mpa 26.19 0.20 h 正态 1300m 650 0.50 n₁正态 -0.1504 -0.015 0.10 n₂正态 0.8529 0.0853 0.10 n₃正态 0.5 0.05 0.10

表2各变量间的相关系数表

变量 c μ σ_xσ_yσ_zh n₁n₂n₃c 1.0 0.2 0 0 0 0 0 0 0 μ 0.2 1.0 0 0 0 0 0 0 0 σ_x0 0 1.0 0.8 0.8 0 0 0 0 σ_y0 0 0.8 1.0 0.8 0 0 0 0 σ_z0 0 0.8 0.8 1.0 0 0 0 0 h 0 0 0 0 0 1.0 0 0 0 n₁0 0 0 0 0 0 1.0 0 0 n₂0 0 0 0 0 0 0 1.0 0 n₃0 0 0 0 0 0 0 0 1.0

3水库诱发地震概率计算

选择不同孔隙水压力p(如：断层孔隙水压力均值分别取162.5m、325m、650m、1300m、 1625m、1950m、2275m、2600m、3400m、4300m、5200m)、得到水库诱发地震概率及可靠 指标如图1所示。

由图1可知，随着断层面孔隙水压力的增加，发震概率大幅增加，断层稳定的可靠指标 明显减小。当孔隙水压力均值低于990m时，断层发震概率在10^-4、10^-5数量级，远低于我国 水库诱发地震平均发震概率10^-3，因此可以认为诱发地震的可能性十分微小；当孔隙水压力 均值达到2600m时，发震概率增加到2.72％；当孔隙水压力均值达到5200m(总静水头)时， 发震概率增加到20.8％，较孔隙水压力较低时增加了10³～10⁴倍，诱发地震的可能性已很高了。 断层面孔隙水压力的大小受库水的下渗影响，如果库区岩性和断裂面渗透性差，库水渗流很 缓慢，下渗水流水头损失也很大，相应孔隙水压力就比较小，其诱发地震的概率也就小；相 反，如果库区岩性和断裂面渗透性好，比如：岩溶发育的地区，岩溶管道连通性好，下渗水 流水头损失小，形成的孔隙水压力大，极易诱发地震的产生，发震时间在水库蓄水后的滞后 时间也短。

上面结合附图对本发明的具体实施方式作了详细说明，但是本发明并不限于上述实施方 式，在本领域普通技术人员所具备的知识范围内，还可以在不脱离本发明宗旨的前提下做出 各种变化。