通过对滑坡变形破坏现场的地质分析和相似材料模型模拟试验研究,对滑坡变形破坏及滑坡对管道的作用机理有了一定的认识。为了更深入地研究滑坡与管道相互作用的机理,验证前述认识,采用有限元计算方法对管道沿线典型滑坡进行了数值模拟研究。
研究方法
对于典型滑坡的滑体、滑动面、滑床、滑坡周界等可以用有限元方法进行离散,着重考虑滑坡潜在滑动面的摩擦接触情况,采用摩擦接触力学方法,对接触单元进行改进,建立接触力学模型,对滑坡进行研究分析。
在考虑可能滑坡面的法向抗张强度和切向抗剪强度下,通过多步迭代计算的方法,获得非线性接触问题的收敛解。通过接触力学计算分析,模拟滑坡的孕育;从滑动面的接触状态的改变(粘着、滑动、未接触状态),获得滑坡的滑动情况,滑坡的动态过程(启动、滑动、停止),同时还进行坡体和管道的应力、变形,以及边坡滑面接触关系分析。在这一过程中重点研究管道在滑坡发生过程中受力破坏的力学机理。
数值模拟
改进接触单元
数值模拟研究中,接触单元的法向与切向模量不一样的,法向模量EN可由两接触物间不可能相互挤入这一事实来确定。在三维情况下,如取N表示接触单元的法向、T1表示切向方向1、T2表示切向方向2,则可推导出:
(1)法向模量EN可由接触物间不可相互挤入的条件确定;
(2)而模量 、剪切模量 ;
(3)在切向平面内的剪切模量 可由接触状态(粘着、滑动等)确定。于是可确定接触单元的弹性矩阵[D],进而可获得三维接触单元的刚度矩阵为:
山区油气管道地质灾害防治研究
式中:[B]是接触单元的几何矩阵。
对包含潜在滑面的整个边坡结构进行有限元离散,可得到边坡结构的刚度矩阵:
山区油气管道地质灾害防治研究
潜在滑面的强度
对于未出现任何破坏的边坡,其潜在滑面是具有一定强度的(实际上滑面还未形成。一般情况下是具有法向强度 和切向强度τ,的,一般情况下 。及 。法向强度的数值可取为边坡潜在滑面处的材料抗张强度,而对于切向强度通常在岩土力学中是取Mohr-Coulomb模型,认为切向强度τs与法向压应力之间σN是一个线性关系:
山区油气管道地质灾害防治研究
式中:c是边坡潜在滑面处的内聚力;ψ是内摩擦角。由上式可确定边坡潜在滑面处的切向强度。而在滑坡发生后,潜在滑面已经变为真实滑面,此时真实滑面位置退化为纯摩擦接触问题。
斜坡潜在滑面的接触状态
在滑坡发生前,潜在滑面处于完好的结合状态,在斜坡自重及外部条件的作用下,接触状态会改变。
非接触状态;当接触单元的法向应力是张应力,且该张应力大于潜在滑面处材料的抗张强度,即 ,此时,滑面将被拉开。接触单元的法向间隙会大于零,即Zg>0。接触单元的模量需取为一个很小的值EN=10-6~10-4Ee,Ee是接触单元的初始模量
山区油气管道地质灾害防治研究
当接触单元的法向应力是张应力,且当接触单元的切向应力τ大于潜在滑面材料的剪切强度τs时,潜在滑面将会破坏,此时Zg>0且EN=10-6~10-4Ee。
山区油气管道地质灾害防治研究
粘着状态:当接触单元的法向应力为压应力,且接触单元的剪应力小于潜在滑面处的抗剪强度时,此时Zg≤0。这意味着接触单元的法向隙Zg是一个很小的近似于零的值
山区油气管道地质灾害防治研究
当接触单元的法向应力为张应力,且接触单元的剪应力小于潜在滑面处的抗剪强度τs,此时
山区油气管道地质灾害防治研究
滑动状态:当接触单元的法向应力为压应力,且接触单元的剪应力大于潜在滑面处的抗剪强度时,此时潜在滑面将破坏,变为真实滑面,而Zg≤0。
山区油气管道地质灾害防治研究
由于剪应力的作用,潜在滑面变为真实滑面,滑坡发生了。接触单元的剪切模量 由下式决定:
山区油气管道地质灾害防治研究
非线性迭代过程的收敛计算
接触状态的最终确定,需满足以下条件:
非接触状态应满足Zg>0且EN=10-6~10-4Ee,而接触状态(粘着、滑动)应满足Zg=0(在指定精度ε下)。
每两次迭代循环计算时,接触单元模量的改变应满足
山区油气管道地质灾害防治研究
其中ΔE、ΔG是相对的模量改变量,ΔEp、ΔGp是许可的相对模量改变量,
山区油气管道地质灾害防治研究
对于滑动状态的接触单元,当下式满足时,滑动才处于稳定状态:
山区油气管道地质灾害防治研究
在以上3个条件同时满足时,该接触问题的计算才是收敛的。
滑坡过程的逐步计算
在边坡滑动过程的数值模拟中,使用了逐步计算的方法。用一系列时间步t1,t2,…,tK,…,tn来分解边坡的滑动过程。而对其中的某一个时间步tK,根据边坡下滑的加速度,对边坡的滑动部分施加相应的惯性力,所以,边坡将保持平衡状态。
在某一个时间步tK时,滑坡系统的势能改变ΔΠ可写为:
山区油气管道地质灾害防治研究
其中:
山区油气管道地质灾害防治研究
在tk时间步,上式中ΔWK是体积力bi和边界力Fs(不包含接触边界—滑面)所做的功,重力势能的改变量为ΔwE,ui为边坡系统的位移量,ΔUK是应变能改变量,σij和εij是应力和应变,Ω为边坡系统的总体积,ΔAK为滑动面上接触应力所做的功,p。和pt为接触面上的法向接触压力和切向力,r为接触面 和 分别为接触面的上、下边界, 和 分别为接触面上的上、下边界的切向位移,ΔKK为动能改变量。由于问题的复杂性,研究中未考虑材料破坏所消耗的能量。
在时间步从t1=0到tn计算完成后,可获得整个滑坡动态过程的应力、应变、接触状态改变情况、滑动距离等重要参数。