罗昌泰,李栋伟,吴 霞
(1. 东华理工大学 土木与建筑工程学院,江西 南昌 330013; 2. 南昌工程学院 土木与建筑工程学院,江西 南昌 330099; 3. 江西应用技术职业学院 建筑工程学院,江西 赣州 341000)
根据防洪能力和坝体的稳定性,尾矿库的安全度分为危库、险库、病库和正常库.在环境、地质、人力破坏等因素的影响下,危库的整体安全性和稳定性最差,最容易产生失稳.历经常年堆积,尾矿危库的容量、堆填厚度一般较大,一旦产生失稳溃坝,容易造成较大范围的地质灾害,危及人民生命财产安全[1-2].因此,对尾矿危库坝体进行维护加固,并对边坡的稳定性进行预测和监控是非常必要的.
现场监控量测技术是一种有效的监测技术,但同时也存在耗时、耗力及现场布置复杂的缺点,不适合长期实施,目前多用于工后评价.数值模拟分析技术近年来被广泛应用于土建工程的各个领域,尤其适用于各类复杂岩土体位移、应力、强度和安全系数的预测评估,常与现场监控量测技术相互配合应用.物理计算参数的准确性在很大程度上决定了计算结果的可靠性[3].而为了获取岩土体的计算参数,采用相关规范、指南、准则中的推荐值进行模拟分析的做法符合计算要求,亦是目前广泛采用的取值方法.然而,自然界中的岩土体是非常复杂的[4],若按照规范进行赋值,那么将难以保证计算结果的准确性、合理性及可靠性.
为了解决该问题,对岩土体真实变形、受力特征实施监测,并对物理特征进行反演的思想相继被提出,即反分析.根据当前的反分析理论,反分析分为位移反分析、应力反分析和位移-应力反分析[5-6].由于位移在现场实测中容易获得,亦是较为可靠的变形特征参数,故基于现场实测位移对力学参数进行反演的位移反分析方法的应用与研究最为普遍.王少杰等[7]基于差分进化法,对在各向异性面上开挖隧洞时的位移反分析相关问题进行研究,建立了参数灵敏度、参数可惟一辨识性及测线布置方案.吴铭辉等[8]利用典型工程案例实测数据进行拟合反分析模型参数.为了模拟分析竖向变形轨迹及范围,ZHANG Z. Z.等[9]以隧道富水破碎带为对象,基于反分析程序对围岩的弹性模量进行反演计算,并利用离散元模拟分析程序对变形追踪算法进行了修正.为了提高收敛速度,GAO W.[10]基于全局优化方法构建了一种免疫连续蚁群优化算法,并据此对岩石参数反分析算法进行了优化修正.刘英棨等[11]为了解决监控数据难以直观反馈施工优化的技术难题,基于位移反分析建立了平面有限元模型,并利用应力释放率与荷载偏压参数,推导了围岩偏压分布公式.
细粒土尾矿坝具有高厚、大容量、大体积特征,采用现场监控量测难以实施长期监控.若依然采用常规的规范法进行参数赋值,即便使用数值模拟分析,由于坝体内部物理特征的巨大差异,也很难保证分析结果可靠.为此,笔者提出一种位移反分析方法,以江西省新干县冷坑冲细粒土尾矿库坝体边坡为实际依托,以现场位移实测值为依据,利用“试凑”式的BMP90位移反分析程序对岩体的等效弹性模量和侧压力系数进行反演,并将反演结果应用于坝体边坡稳定性、排渗加固措施的数值模拟分析,以期为细粒土尾矿库坝体边坡稳定性及加固提供一定的参考和借鉴.
冷坑冲尾矿库是江西省新干县的一个萤石矿尾矿库,萤石原矿年产量为6×104t,萤石精粉为3×104t,尾矿库总库容为7.17×105m3.经长年累月堆填,坝体最终堆积高程达190 m,相对堆高达64.9 m,为四等库.2019年5月,在运行过程中因库底排水斜槽发生断裂垮塌,造成尾砂和尾废水泄露.泄露后在库尾形成塌陷坑,经测算溃泄方量约5.00×104m3.
该矿所排放的萤石尾砂颗粒粒径极细,多为尾粉砂、尾粉土和尾粉质黏土.通过钻孔取得原尾矿砂样进行粒度成分分析,结果如表1所示.
表1 尾矿砂样粒度成分统计表 %
分析发现粒径大于0.075 mm的颗粒质量分数为31.83%,属于典型的细粒尾矿[12].细粒尾矿库的普遍特征是透水性极差,浸润线不易降低,导致振动液化,且软泥层中可能存在较高的超孔隙水压力,抗剪强度很低,致使坝体结构的静、动稳定性弱于一般尾矿,病害率偏高[13-14].
事故发生时,冷坑冲尾矿库没有新的尾砂浆液排放,周边亦无其他地质灾害发生.根据当地气象资料,全年丰水期(4— 6月)降水量占全年的59.5%,易形成洪涝灾害.事故发生前,当地连续降雨累计天数达到18 d,降雨使得库区积水无法及时排出库外,浸润线持续升高,导致尾矿坝体强度降低,加之如果洪水漫顶,坝体将被冲刷破坏而失稳.只有水位降落和上升历时越长,水力梯度的变化幅度越小,对渗流安全影响才越小[15].因此,浸润线升高是导致坝体失稳的主要风险因素,事故的发生导致该尾矿库处于危库状态.
为了对事故发生之前的坝体稳定性和渗流进行分析,采用Midas GTS NX有限元分析软件进行模拟推算.该软件操作界面简洁,模型库功能分析齐全,内置的FS模块具有优越的流-固耦合分析计算能力,运算时解的稳定性与收敛性较好[16].
尾矿库堆积坝的高度为64.9 m,初期坝至库尾的长度约为458.0 m.通过钻探手段发现,尾砂分层为尾中砂、尾细砂、尾粉土和尾粉质黏土,以下为花岗斑岩.采用二维模型进行模拟分析.网格划分单元总数为789个,节点数量为1 092个,如图1所示.图中1#-6#位置分别为干砌块石初期坝、尾中砂层、尾细砂层、尾粉土层、尾粉质黏土层和原花岗岩山体.
图1 模型网格划分
与粗粒土相比,正常固结条件下的细粒土层的结构更加致密,应力、位移传递更加有效[17].因此,在2#-5#位置的网格划分更加紧密,1#、6#位置的网格单元数量与2#-5#的数量比为1 ∶1.3.根据GB 50330—2013《建筑边坡工程技术规范》,认为在滑坡体周围0.5~1.5倍坡高范围外不受滑坡的影响,因此网格划分的横向长度取1.5倍的坡高.实际计算时认为网格边界上不存在应力传递与水头,选择摩尔-库伦应力准则,不设置固定的迭代次数,视应力分配达到平衡为基本原则,并由软件自行运算[18].考虑到滑坡垮塌事故已经发生,因此可以确定破裂面的位置,故选择应力极限平衡法对应力进行分配迭代.经运算后发现,安全系数为1.125 0.根据GB 50021—2001《岩土工程勘察规范》可知,重要工程边坡的安全系数取值应达到1.3~1.5,因此该尾矿库安全系数不符合要求,在现有的地质条件下极易发生垮塌灾害.
为分析不同浸润线高程所对应的安全系数,需要根据水头状态计算浸润线分布情况.在Midas GTS NX计算中找出尾矿库的浸润线(孔隙水压力p为0的线),根据坝顶与临时排洪设施的高程差,确定正常水头H=71.5 m,洪水期水头H=73.5 m.图2和3分别为尾矿库当前和洪水期浸润线计算云图,其中白色部分为孔隙水压力为0时的等值线区域,即浸润线.对比图2、3,可知浸润线整体稍微降低.
图2 尾矿库当前浸润线分布云图
图3 尾矿库洪水期浸润线分布云图
利用Midas GTS NX软件中“应力-渗流-稳定”计算模块,分别计算水头H=71.5、73.5 m时坝体的渗流安全系数,计算结果分别如图4、5所示,其中F为安全系数,d为位移.
图4 H=71.5 m时坝体渗流安全系数计算结果
图5 H=73.5 m时坝体渗流安全系数计算结果
对比图4、5可知:水头增高导致浸润线持续抬升;在渗流状态下,土体内部的孔隙水压力增大,有效应力降低,使得土体抗滑力减小;同时,土体内部含水率增大,土体容重也将增大,加大了边坡的下滑力,致使安全系数有所下降.经计算,水头H=71.5、73.5 m时的安全系数分别为1.012 5和0.938 2.可见,连续降雨导致浸润线持续抬升,弱化了坝体边坡的整体抗滑力,加大了滑塌风险.
基于前述分析,本次隐患处理不仅需着重增设排渗设施降低浸润线外,还需对坝体进行加固处理以增强坝体强度.表2、3分别为常见的排渗和坝体加固措施.
结合该尾矿库细粒尾砂渗透系数低的特性,对比表2、3可知:在排渗方面,水平排渗和辐射井排渗的效率更高,且施工便利;在固坝方案中,钻孔碎石桩、抗滑桩支挡法和辐射井加固法均适用.综合考虑排渗和固坝效果,同时兼顾施工便利性,选择辐射井为本项目的除险加固措施.图6为辐射井的三维模型示意图.图7为辐射井的平面布置示意图.
辐射井由竖井、水平排渗管和导水管组成,竖井为钢筋混凝土圆柱形结构,在现场配筋,并分节浇筑.成型后逐节下沉至设计高程,然后再封底就位.
表2 常见排渗措施
表3 常见坝体加固措施
图6 辐射井模型示意图
图7 辐射井平面布置示意图
根据坝体范围和浸润线高程,需布置3个辐射井,即1#、2#和3#辐射井,井身高度分别为23、25和27 m,内径均为3 m.导水管采用Φ114的钢管,长度分别为78、80和98 m,坡度为1%.水平排渗管设置3层,管长50 m,纵坡为5%,且每层垂直间距为4 m.
采用Midas GTS NX软件对排渗效果进行模拟分析,以掌握浸润线的变化状态.本次模拟分析主要考虑辐射井的支挡、加筋和固结作用.采用三维模型计算加固后坝体的稳定性.图8为建成的数值计算网格模型,其中单元数为137 352个,节点数为125 199个.
图8 数值计算网格模型
图9为辐射井施工后浸润线位置变化云图.图10为辐射井施工后渗流稳定性云图.
图9 施工后浸润线位置变化云图
图10 辐射井施工后渗流稳定性云图
由图9可知,辐射井施工后,坝体内部的浸润线快速下降,孔隙水压力消散较快.由图10可知,渗流安全系数的模拟分析结果为1.323 4,较排渗加固前有显著提高.
基于BMP90位移反分析程序来获取坝体边坡的力学参数,即以边坡位移实测值为反演依据,对边坡的等效弹性模量和侧压力系数进行反演,并将反演结果应用于有限元数值计算中.根据GB 50021—2001的有关规定,同时考虑到现场监控量测任务重、设备不足等因素,实际监测时仅对边坡的下滑位移进行监测,精度为0.1 mm.共设置3个监测点,分别位于1#、2#和3#辐射井的左侧,距离排渗井的水平距离均为3 m.BMP90反分析软件的基本原理是通过不断地调整参数取值,最终使得实测值与拟定值的差值满足某种容许范围,此时拟定值视为可用[19],即
式中:R为容许范围值;K1为位移拟定值,mm;K2为位移实测值,mm.
根据GB 50021—2001的规定,初步决定将观测时间定为30 d.由于实际施工中诸多因素的不利影响,最终实际观测得到24 d的实测值.表4为反分析结果.
表4 反分析结果
由于全体样本取值的波动幅度较大,因此从数理统计角度考虑[20],与求取全体样本的平均值相比,选择样本值波动趋于平稳阶段的样本平均值会更加合理.分析表4可知,从第21天开始侧压力系数波动趋于平稳,故选择第21-24天的平均值为0.57;从第20天开始等效弹性模量波动趋于平稳,故选择第20-24天的平均值为650.97 MPa.
分别在模型左右两侧、底部施加x、y方向的约束,地表为自由边界,假定边坡为均质的弹塑性体,并遵循摩尔-库伦准则[21-22].网格模型见图8.首先采用GB 50021—2001给出的V级边坡土体计算参数进行位移计算;然后采用由位移反分析得到的侧压力系数和等效弹性模量进行位移计算.将监控量测所得实测位移与位移的有限元计算值进行对比分析.计算结果如表5所示.由表5可知:通过反分析所得计算结果比根据规范给定参考值的计算结果更为接近实测值,这说明通过反分析所得计算结果更加真实有效;1#、2#和3#监测点的位移值逐渐下降,这是因为1#监测点距离凹陷坑最近,也就是距离自由变形边界最近,因此其位移值最大;模拟分析结果小于实测结果,这是因为在模拟分析时,在模型左右两侧、底部施加了x、y方向的约束.
表5 位移计算值与实测值对比 mm
综上,坝体边坡经加设排渗井后,浸润线明显下降,数值模拟分析和实测结果也表明坝体边坡整体趋于稳定,表明加固措施有效可行.
1) 细粒土尾矿库透水性差,浸润线不易降低,导致振动液化,在超孔隙水压力的作用下,坝体结构的静、动稳定性弱于一般尾矿.
2) 辐射井排渗系统能有效降低细粒土尾矿库内浸润线,同时发挥加固、排渗作用,能够较快稳定坝体边坡,防止垮塌.
3) 位移反分析是以实测位移值为依据,对土体的计算参数进行反演,避免了规范法赋值的不足,因而所得模拟计算结果与实测结果更加吻合.