基于Monte-Carlo法的尾矿库坝坡稳定可靠度分析

张 宁，李玉雷，李哲辉，郭大平，胡良才，张 娟

(1.中核第四研究设计工程有限公司，河北 石家庄 050021;2.河北省水利规划设计研究院有限公司，河北 石家庄 050000)

铀尾矿库是铀矿山安全生产的重要设施之一，具有容量大、势能高、堆存时间长等特点。铀尾矿库一旦发生溃坝事故，不但会造成严重的环境污染，也会给下游居民造成巨大且难以挽回的损失[1]，因此尾矿库的坝体边坡稳定性问题一直是工程设计和尾矿库安全分析的重点。边坡稳定性计算结果主要受物理参数的不确定性、模型的不确定性和统计的不确定性等因素影响[2]。目前，国内外采用传统的定值分析法来分析边坡稳定问题，此方法通常把经验确定的安全系数K作为边坡稳定性度量指标[3-4]；然而由于考虑不确定性影响而设置的“裕量”(如K≥1.25)往往带有经验性质，因此即使K相同，也不能反映相同的稳定性水平[5]。

近年来，基于概率统计理论的可靠度分析方法已逐步建立起来，将各种不确定性因素纳入到边坡稳定性计算中，用可靠性指标来衡量边坡的安全程度[6]。常用的可靠度分析方法主要有点估计法、一次二阶矩法、蒙特卡罗法(Monte-Carlo法)等[7]。Monte-Carlo法适用于随机变量概率密度分布函数已知或符合假定的情况[8]，在目前高速计算机上该方法通用性强，计算相对精确。

笔者以南方某退役铀尾矿库为例，分析坝体岩土层物理力学参数对尾矿库坝体稳定可靠度的影响。在考虑岩土参数不确定性和滑动面空间变异性基础上，选取服从正态分布的尾矿库坝体及基岩物理力学参数为随机变量，采用Monte-Carlo法模拟计算不同工况下尾矿库坝体的稳定可靠度，并对不同工况下坝体滑坡风险进行分析，为尾矿库坝体安全性评价提供参考。

1 铀尾矿库概况

某退役尾矿库的初期坝为碾压式均质土坝，坝顶标高316.4 m，坝高20 m，下游坡比(1∶2.6)～(1∶3.0)，坝脚设有8 m高的排水棱体并有排水褥垫伸入坝体内。堆积坝顶标高355.9 m，堆积坝坡比(1∶3.2)～(1∶3.6)，在标高约328.1 m处设有宽约22.0 m的平台。在初期坝顶内侧及标高约328.1 m堆积坝平台的坡脚处，各有一道伸入坝内且平行坝轴的排渗体，排出坝体渗水以降低浸润线。

坝体座处在坡积、残积黏土和基岩上。根据《中国地震动参数区划图》(GB 18306—2015)，尾矿库所在区域抗震设防烈度为Ⅵ度，设计基本地震加速度为0.05g。尾矿坝体剖面及土体各层信息如图1所示。

图1 某尾矿库坝体剖面示意图(单位：m)Fig. 1 Schematic diagram of dam section of a tailings pond (unit: m)

2 坝体稳定可靠度表达式

2.1 定值分析法安全系数表达式

极限平衡法是尾矿库坝体稳定分析的基本手段，本次采用瑞典条分法对坝体进行稳定计算。假定边坡滑动面为平面应变问题的圆柱面，且认为滑动土体为刚体、划分条块之间不会变形，通过条块受力平衡方程来求解坝坡最小抗滑稳定系数。对于简单边坡形式，圆弧滑动面的安全系数计算公式为[9]

(1)

式中：Ks为安全系数；c为粘聚力，kPa；L为滑弧长度，m；φ为内摩擦角，(°)；γ为土体容重，kN/m3；b为土条宽度，m；n为土条数；hi为第i个土条的平均高度，m；θi为第i个土条滑弧中点法线与竖直线的交角，(°)。

2.2 Monte-Carlo法失稳风险率表达式

Monte-Carlo法的理论基础是伯努利大数定理，即当实验次数足够大时，事件出现的频率几乎接近于其发生的概率[10]。利用Monte-Carlo法对坝体稳定可靠性进行计算时，可直接从给定的不确定性参数统计特征中产生随机样本，并计算边坡的可靠度，其结果具有无偏性。

计算时可根据岩土参数正态化处理后，得到N组随机样本，N组随机样本对应得到安全系数K(Ni)(i=1，2，…，N)，则坝体边坡失稳风险率可表示为[11]

(2)

式中：Pf为边坡失稳风险率；Nf为K(Ni)<1失效样本个数；I(·)为指示函数，可根据《核工业铀水冶厂尾矿库、尾渣库安全设计规范》(GB 50520—2009)[12]规定的边坡稳定安全系数Ks确定。若某组样本K(Ni)≥Ks，则I(·)=0，认为边坡处于稳定状态；否则I(·)=1，认为边坡存在失稳风险。

3 尾矿坝坡稳定可靠度分析

3.1 定值分析法分析坝坡静力稳定性

根据尾矿库工程勘察结果，尾矿库各土层物理力学参数见表1。坝坡稳定分析采用摩尔-库伦模型，在不影响坝坡稳定分析结果前提下对坝体土层结构进行适当简化，简化后坝坡稳定性计算结果如图2所示。

图2 尾矿库坝坡稳定计算结果Fig. 2 Calculation results of dam slope stability of tailings pond

表1 某尾矿库岩土物理力学参数Table 1 Geotechnical physical and mechanical parameters of a tailings pond

由静力分析结果可知，在正常工作条件和非常工作条件(考虑地震影响)下坝坡抗滑安全系数分别为1.506和1.256，均大于《核工业铀水冶厂尾矿库、尾渣库安全设计规范》规定的限值，尾矿坝坡稳定性满足要求。但在计算过程中，极限平衡分析的参数选取往往是采用土工试验结果的平均值，没有考虑参数的不确定性影响[13-14]，这显然与实际有所出入。因此，仅从安全系数难以对尾矿坝坡稳定程度进行判定，无法定量分析坝坡的失稳风险大小。

由图2尾矿库坝坡滑弧范围可知，滑弧主要经过尾粉细砂和黏土初期坝部位。为确定滑弧经过部位岩土力学参数对安全系数的影响，对滑弧范围内岩土参数进行敏感性分析，结果如图3所示。结果表明坝体土内摩擦角对安全系数影响较大，这与实际相符，因为滑动面主要穿过黏土初期坝。

图3 坝坡滑弧范围敏感性分析结果Fig. 3 Sensitivity analysis results of dam slope slip arc range

3.2 Monte-Carlo法分析坝坡稳定可靠度

岩土力学参数可能服从正态分布、对数正态分布、极值分布等有限的几种分布类型[15]。影响坝体边坡稳定性的岩土参数主要为尾矿细砂和初期坝体土的物理力学指标，通过对随机岩土力学参数样本拟合检验分析，黏聚力c、内摩擦角φ均服从正态分布。物理力学指标特征值见表2，在可靠度分析时引入参数特征值，将岩土物理力学参数的不确定性考虑进来。

表2 岩土随机物理力学参数指标特征值Table 2 Characteristic values of random physical and mechanical parameters of rock and soil

对坝坡滑动面随机抽样10 000次，进行模拟计算，得到尾矿库不同工况下坝坡抗滑稳定安全系数的概率密度分布函数和概率分布函数，如图4～5所示。可以看出，坝坡抗滑稳定安全系数呈较好的正态分布特征，且由概率密度正态分布峰值可得到坝坡安全系数值，由概率分布横坐标中找到规范规定的安全系数值，其所对应纵坐标即为失稳风险率。

图4 N=5 000时正常工况下安全系数概率密度和概率分布Fig. 4 Probability density and probability distribution of safety factor under normal working conditions when N=5 000

图5 N=5 000时非常工况下安全系数概率密度和概率分布Fig. 5 Probability density and probability distribution of safety factor under abnormal working conditions when N=5 000

尾矿库坝坡稳定可靠度计算结果见表3。可以看出，在不同工况下坝坡安全系数均值满足规范要求，坝坡处于稳定状态。依据规范规定的尾矿库安全系数限值Ks确定统计功能函数K(Ni)≤Ks的次数Nf，由Pf=Nf/N可得到尾矿库坝坡在正常工况下的失稳风险率为0.20%，坝坡风险率很小；安全系数标准差为0.102，极差为0.804，安全系数离散程度较大。尾矿库坝坡在非常工况下的坝坡失稳风险率为13.10%，失稳风险率远大于正常工况；安全系数标准差为0.085，极差为0.673，安全系数离散程度较正常工况小。

表3 尾矿库坝坡稳定可靠度计算结果Table 3 Calculation results of dam slope stability reliability of tailings pond

4 结论

尾矿库坝坡稳定性主要影响因素为坝体土的内摩擦角。采用Monte-Carlo法与定值分析法计算得出的安全系数一致，坝坡在正常和非常工况下均处于稳定状态。相较于用安全系数判定坝坡稳定性，可靠度分析从安全系数和风险率双维度进行判定，可以更全面反映尾矿库坝坡的长期稳定性，并且能够充分利用岩土物理力学参数，更客观可靠地评价坝坡稳定性。