某铀尾矿库地震动力稳定性分析

2019-08-13郭大平胡良才李玉雷李哲辉

铀矿冶 2019年3期

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

(中核第四研究设计工程有限公司，河北石家庄 050021)

地震是影响铀尾矿库安全稳定的重要因素之一。在役铀尾矿库内有大面积的水体，坝体浸润线较高，大部分坝体尾矿处于饱和状态，对振动荷载敏感，在地震作用下易发生液化和滑坡等破坏作用[1]。中国对尾矿库的研究主要集中在坝体稳定、环境影响、安全监测和运行维护等方面[2-3]。

中国铀尾矿库库容普遍较小，工程等级大多为三等库，其动力响应分析大多沿用传统岩土边坡、土石坝等分析方法，如静力法、拟静力安全系数法、Newmark类永久位移计算方法[4]等；计算精度高、拟合度好但计算过程复杂、计算耗时较长的有限元时程分析方法还没有广泛应用于铀尾矿库动力分析中。笔者采用有限元时程动力分析法对某铀尾矿库坝体在设计地震条件下进行动力模拟计算，并对尾矿坝体的动力稳定性进行分析。

1 工程概况

1.1 尾矿库基本情况

某铀尾矿库处于中-低山区，区内沟谷纵横、多呈“V”字型，植被发育。尾矿库属三等库，采用上游法筑坝，由初期坝和尾矿堆积坝组成，总库容约260万 m3。初期坝为碾压式均质黏土坝，坝顶标高125 m，坝高25 m，坝长145 m，外坡坡度1∶2～1∶2.5，内坡坡度1∶2.5，坝外坡利用块石砌筑了排水棱体。尾矿堆积坝顶标高150 m，坝坡坡度1∶3.5，马道标高141 m，在初期坝顶以上部分的尾矿坝坡利用块石修筑了压护平台。库内堆积尾矿分为尾中砂、尾细砂、尾黏土三层，尾矿滩面水面标高147 m。尾矿库平面布置如图1所示。

图1 某尾矿库平面布置

1.2 尾矿库地质地震情况

根据GB18306—2015《中国地震动参数区划图》，尾矿库区域属于小于Ⅵ度地震烈度区；根据当地地震动参数区划工作用图，工程所在区域未来50年内设计基本地震加速度值小于0.05g。场地土属中软土，基岩埋深15～80 m，设计地震分组为第一组，等效剪切波速100 m/s，地震动反应谱特征周期为0.35 s。考虑一定安全裕度，本模拟计算中取地震加速度峰值为0.05g。

2 计算原理

2.1 动力运动方程

根据Seed-Lee-Idriss动力分析方法[5]，地震动力响应的有限元运动方程可表述为

(1)

为提高计算效率，程序中质量矩阵采用集中质量阵的形式进行计算，计算公式为

(2)

式中：ρ—质量密度，kg/m3；[ψ]—质量分布因数的对角矩阵。

2.2 弹性本构关系

根据弹性力学理论，动力系统的应力应变关系描述为

{σ}=[C]{ε}

(3)

式中：[C]—有限单元本构矩阵，其表达式为

式中：E—杨氏弹性模量，MPa；ν—泊松比。

2.3 等效线性材料本构模型

计算中采用等效线性本构模型来定义岩土材料，该模型以弹性模型为基础，考虑了土体刚度随动力变化的折减。程序最初为土体指定一个特定的刚度值，然后计算本时程的有限单元峰值剪应变。利用该剪应变值和事先定义的刚度折减函数，程序会生成一个新的材料刚度值，开启下一时程的计算。整个过程不断重复直到刚度的改变量进入事先指定的范围。

3 有限元计算

3.1 单元剖分及模型边界条件

为便于建模和节省计算时间，对尾矿库地层进行适当简化处理。在GeoStudio软件中建立简化的尾矿库模型，并对模型进行单元剖分。模拟计算的边界条件：库底边界位移固定；初期坝与排渗体交界面为潜在渗出面；滩面水头为147.0 m。尾矿库计算模型如图2所示，高程系采用1954年北京坐标系，1985年国家高程。

图2 尾矿库计算模型

3.2 尾矿库材料物理特性及本构模型函数

尾矿库坝体材料的静力学物理参数取自尾矿库工程地质勘察报告，如表1所示。

尾矿库动力分析采用等效线性模型，所用到的动力学参数以函数关系的形式确定。由于没有进行相关室内试验确定材料的动力变形特性，因此参考国外类似工程案例和相关文献的试验结果及建议取值[7-11]，并结合本尾矿库坝体材料特性对函数关系进行修正，最终确定了尾矿库材料的刚度折减、孔隙水压比和剪应力比的函数关系。

表1 尾矿库材料的岩土物理力学参数

坝体材料孔隙水压比的变化关系如图3所示，图中横坐标的循环数目比(N/NL)是指等效动力循环次数(N)与液化动力循环次数(NL)的比值，其中N为等效循环应力(最大剪应力的0.65倍)作用下换算得到的动力循环次数，NL为在等效循环应力下土体发生液化时的循环次数。

3.3 初始应力计算

建立QUAKE/W初始静态分析，设置材料刚度折减函数、孔隙水压比函数、剪应力比函数等相关本构模型参数，得到坝体初始应力状态如图4所示。可以看出坝体初始静态最大总应力基本呈层状分布，由坡内向坡外逐渐减小。

图3 尾矿库材料的孔隙水压比函数关系

图4 初始静态尾矿库最大总应力分布

3.4 地震波的输入

选用地震波时考虑的主要因素有地震动强度、地震波的频谱特性和地震波的持续时间。1971年2月9日在美国南加州圣费南度地区发生的圣费南度(San Fernando)地震，里氏震级5.4级，地震动峰值加速度0.09g，地震持时35 s，地震区场地与目标尾矿库所在场地同为中软性状土。圣费南度地震波加速度峰值和本次计算设计峰值接近，地震波的频谱特性相近(场地相近)，且地震持时较长，因此可选择该波作为本次计算使用的地震波。

地震荷载对坝体的破坏主要是由水平地震应力引起，根据GB50863—2013《尾矿设施设计规范》，竖向地震峰值按水平地震峰值的2/3进行折减。对该波进行滤波调幅，调整后地震水平加速度峰值为0.05g，竖向加速度峰值为0.34g。截取地震持续时间为10 s，输入水平地震波(如图5)。

图5 水平地震动加速度时程曲线

3.5 动力反应计算结果

3.5.1坝体应力

地震过程中的坝体最大有效应力和超孔隙水压力分布情况如图6、7所示。

图6 坝体最大有效应力分布

图7 坝体超孔隙水压力分布

由计算结果可知，坝体最大总应力分布规律与初始静态下的分布情况基本一致，且在整个过程中均为正值，即土体单元始终为受压状态。最大有效应力从下至上呈先降后升的变化规律，这是由于该部位的尾矿厚度最大，对地震动有较强的削弱作用；同时该区域的孔隙水压最低，抑制了有效应力的降低。在动荷载作用下，尾矿库内出现超孔隙水压力，自坝体中部呈放射状分布。坝体最大剪应力分布规律与最大总应力的分布较接近，但最大剪应力的应力集中现象较最大总应力更为明显。

3.5.2地震下坝体加速度反应

为反映尾矿坝加速度随时间的变化趋势，在初期坝顶(125 m)、堆石压重平台(131 m)、马道顶(141 m)和堆积坝顶(150 m)各设置一个监测点，监测点分布位置如图8所示。

图8 坝体监测点示意

通过计算得到监测点位置的加速度放大倍数如表2所示。初期坝顶和堆积坝顶位置的水平方向加速度时程曲线如图9、图10所示。可以看出监测点位加速度峰值出现的时间与地震波时程曲线基本相符，即均出现在强震阶段；初期坝坡与排水棱体顶面附近的加速度放大效应最为明显。随着坝体高度的增加，坝体底部土体单元的自振周期逐渐加大，孔隙水压力也在增加。由于本次模拟所采用的地震波强度较低，地震剪应力难以透过尾矿坝体向上大距离传播；因此，坝高较高区域与模拟地震的卓越周期相差越来越远，加速度放大倍数随坝体高度的增加而略有降低。

表2 监测点加速度放大倍数

图9 初期坝顶加速度时程曲线

图10 堆积坝顶加速度时程曲线

3.5.3坝体变形

根据程序计算结果，在整个地震过程中，坝体X方向的最大峰值位移为3.1 cm，出现在堆积坝顶附近；Y方向的最大峰值位移为1.6 cm，出现在靠近堆积坝顶的滩面附近。坝体X方向峰值位移分布云如图11所示。

图11 坝体X方向峰值位移分布

坝体最大残余变形集中分布在坝体中部的尾中砂和尾细砂两层上，地震发生后坝体中上部的尾中砂和尾细砂两层发生了较明显的指向下游方向的水平向变形。根据程序计算结果，X方向最大残余变形量为1.8 cm，Y方向最大残余变形量为0.8 cm。坝体残余变形矢量分布如图12所示。

图12 坝体残余变形矢量分布

3.5.4坝体液化

尾矿库坝体浸润线以上的部分是非饱和区域，可以直接判定为非液化区域；浸润线以下区域则通过引起残余孔压所需要的最小剪应变幅值决定。计算后得出的尾矿库液化区域分布情况如图13所示，库尾部浅色区域即为液化区域。

图13 尾矿坝体液化区域分布

液化区域主要集中分布在沉积滩面尾部的尾黏土浅层区域，且分布面积很小，不会对坝体稳定构成威胁。由于本次模拟地震强度较低，地震动造成的坝体孔隙水压力上升的范围和幅值有限，这就决定了坝体液化的范围不会很大。同时该区域尾黏土是颗粒最细的坝体材料，长期位于水面以下，浅层土体呈松散状，在地震荷载作用下容易产生液化，说明本次模拟计算的结果基本合理。

3.5.5坝坡动力稳定分析

建立SLOPE/W分析，按动力时程法对尾矿库坝坡的抗震稳定安全系数进行计算，得到时程曲线如图14所示，地震过程中的最小安全系数滑移面如图15所示。

图14 坝体稳定安全系数时程曲线

图15 坝体稳定最小安全系数Kmin对应的滑移面(t=5.79 s)

在设计地震发生的过程中(10 s)，坝体稳定安全系数大部分时间都处于较高的范围，大于1.3的累积时间达到6.43 s，其中最大安全系数Kmax=1.43，发生在5.93 s；最小安全系数Kmin=1.197，发生在5.79 s。根据GB50520—2009《核工业铀水冶厂尾矿库、尾渣库安全设计规范》，三等尾矿库在地震条件下的最小安全系数限值为1.15，可见该尾矿库在遭遇设计地震时是稳定的。