汶川震区桃关沟泥石流浆砌块石坝损毁数值仿真研究

蔡红刚，裴向军，唐秋元，何 妍

（1.中煤科工集团重庆设计研究院有限公司，重庆400039；2.成都理工大学环境与土木工程学院，成都610059）

桃关泥石流沟是一条老沟，该沟于1962年6月15日发生泥石流，位于北纬31°11′50″，东经103°38′5″，海拔1 100～1 600m，属沟谷狭长地段。汶川“5.12”地震后导致沟内不稳定斜坡、滑坡、崩塌等极其发育，为泥石流的暴发提供了丰富物源。相关部门于2010年在沟道下游修建了浆砌块石拦挡坝，该坝被后期暴发泥石流损毁。通过对汶川桃关沟泥石流流域特征、物源以及拦挡坝的损毁特征调查，反演分析该浆砌块石坝的损毁机理。

1 汶川桃关沟泥石流特征

1.1 泥石流沟道特征

桃关沟长11.4km，汇水面积51km2，平均坡降128.6‰，两河口经上至桶棚梁子为汇水区，植被发育，覆盖率达60%以上，两河口至桃关工业园区为形成流通区，坡降为90‰，植被覆盖率约30%，两侧山体坡度一般30°～50°，松散堆积层较丰富，特别是地震后在距园区200m形成了大量的松散碎块石堆积体，一般粒径为0.3～1.2m，最大粒径可达8m，为泥石流的发生提供了物源。流通区沟道顺直，未见陡坎，基本无卡口，物质组成以碎石土为主，黏性小。工业园区至桃关村段为流通区，坡降为70‰，建桃关工业园区时已建排导槽。堆积区域为桃关沟下游至岷江河。

图1 桃关沟泥石流分区图（参考勘查报告）

1.2 水文及物源特征

据有关气象资料统计，汶川一带多年平均降雨量为524.03mm，极强降雨量为55.6mm／d。该泥石流沟常年流水，部分水流渗入地下，补给地下水。

桃关沟泥石流的物源主要为两方面：一是地震后在沟道及坡脚等形成大量的松散崩坡积层，一般粒径为0.3～1.2m，最大粒径可达8m；二是沟道内存在的冲洪积层，磨圆度较好，粒径大小一般在2～200mm之间，少量漂石粒径达400mm。

1.3 泥石流的平均流速及撞击力

通过访问和对沟内物质进行分析判断，并对泥石流物质进行现场配比试验，容重为1.4～1.5之间，取值1.5，比重按2.65计算。

泥石流平均流速以及冲击力的计算参照规范DZ／T 0220—2006《泥石流防治工程勘查规范》中相关公式，计算过程见文献［1－3］，计算结果为VC＝5.09m／s；泥石流冲击力计算按最危险情况考虑，该冲击力计算包括整体冲击力和单块最大冲击力计算，计算只考虑正面冲击。计算结果为泥石流流体的冲击力F＝13.1kN／m2，泥石流大块石的冲击力为FS＝57.0kN。

2 桃关沟泥石流浆砌块石坝损毁特征

汶川震区桃关沟泥石流浆砌块石坝的坝体材料为M10浆砌片石，片石单抗压强度不底于30MPa，且为不易风化的硬质岩石。以桃关沟2＃浆砌块石坝的损毁特征来分析该类型的坝体受泥石流块石冲击损毁的特征和机理。现场调查发现，坝体受大块石的撞击后，坝体撞击区表面出现纵向裂隙和横向裂隙，延伸长度0.5～2.5m，宽度约3mm，直接撞击区受损严重，坝体出现张裂、掉块，导致坝体的整体性受到严重破坏。未受大块石直接冲击的区域，局部出现小的裂缝，坝体整体性受轻微破坏。受块石撞击损毁的坝体在后期泥石流的持续冲击下，凹坑逐渐变大，坝体整体极易出现溃坝，如图2。

图2 汶川桃关沟泥石流2＃坝

3 浆砌块石坝受冲击破坏数值仿真

桃关沟泥石流大块石多为4～5m3，形成的块石冲击力较大，导致坝体受冲击的变形较大。ANSYS／LS－DYNA软件比较适合于大变形的冲击动力模拟，结合该软件建立1∶1的仿真模型，反演分析泥石流浆砌块石坝受块石冲击损毁的过程和破坏机理。ANSYS／LS－DYNA软件的相关原理内容，可见参考文献［4－5］。

3.1 模型特征

汶川地震后，桃关沟内松散物源丰富，并于2009年7月23日暴发泥石流，其下游的2＃浆砌块石坝受到严重损毁，其损毁原因主要是受大块石冲击损毁。本文选择拦挡坝受块石冲击损毁段（最危险段）为研究对象，建立如图3所示的三维模型。野外调查的桃关沟泥石流中大块石为4.2m3，该数值仿真的主要目的是模拟坝体受大块石冲击的受力状态和损毁特征，块石外形对模拟的结果影响较小，故本文采用球体作为大块石，考虑体积相近的原则，选用直径为2m的实心球体等效。

图3 块石冲击拦挡坝三维模型

由于坝体受法向冲击是最易造成损毁的，故本文主要模拟坝体受块石的正向冲击损毁。该仿真模型以竖向为Z轴，顺沟方向为Y轴，垂直沟道方向为X轴，建立1∶1的等效模型。模型由3个实体组合而成，考虑到砌筑坝体的块石粒径约为0.20m，网格划分为0.20m方体，总的单元为70 280个；土体采用映射网格划为9 850个单元；块石单元为3 986个，整个模型共有89 580个节点，总共为84 014个单元。

由于本文模拟的是坝体中部段，坝体两侧按固定边界条件与实际相符，大块石的初始速度按计算速度5m／s取值（模拟时按Z向速度为－2.2m／s、Y向速度为4.5m／s）采用。大块石与坝体结构面的接触关系为CONTACT＿AUTOMATIC＿POINT＿TO＿SURFACE（点—面自动接触）。模拟时间的确定主要依据考虑到泥石流中大块石冲击坝体的过程是实际极短的动力过程，冲击作用为几毫秒，坝体受冲击作用后应力—应变存在调整的过程，整个冲击模拟时间选择为1.2s，其实际冲击时间为0.1ms左右，冲击角度采用正向冲击，取冲击角度为90°。

3.2 材料参数选取

对于大块石，冲击过程中不考虑本身的变形，需选用变形小，强度较高的材料按刚体考虑；对于坝体和土体，为弹塑性变形体，选择弹塑性本构模型的材料即可满足要求。依据Hertz接触理论［6－8］，假设材料 满足Vonmises屈服准则条件，材料参数见表1。

表1 材料参数表

3.3 数值仿真分析

3.3.1 浆砌块石坝体冲击损毁特征分析

数值仿真切片显示，坝体受泥石流大块石的冲击过程在极短时间内完成，其损毁特征主要表现在坝体的变形和损毁区严重程度2个方面。首先，坝体在受大块石冲击的时候，坝体表部发生较大的变形，其变形随冲击的继续逐渐增大，在坝体的变形量达到5～7cm，即超过该类坝体的抗变形极限值时，块石将会脱落损毁，在多次冲击作用下，坝体损毁严重；其次，仿真结果显示，坝体受大块石直接冲击的区域损毁严重，而未受大块石直接冲击的区域受损程度较轻，为轻微损毁区。仿真结果与实地调查的坝体损毁特征基本一致。

图4 坝体损毁特征

3.3.2 大块石冲击力与时间的关系

根据能量守恒原理，大块石冲击力与其体积、冲击速度及块石岩性（强度）等成正相关。块石冲击开始时，其速度最大，随着冲击过程的发生，动能转化为变性能，撞击接触力达到最大值，之后接触力逐渐减小。

泥石流大块石对坝体的冲击为法向冲击，其冲击力—时间关系如图5。仿真结果显示，法向接触力在时间0.599 98～0.627 994s之间从0达到最大值66.31kN，其冲击接触时间极短，该冲击过程为瞬时冲击完成的。理论计算泥石流大块石冲击力为FS＝5.7tf，即为57.0kN，模拟的最大冲击力较计算值大9.31kN，偏差范围为14.04%，其主要原因是数值仿真模拟计算是一种理想的状态，较相关经验公式计算更与实际吻合。数值仿真的结果在一定程度上可以指导该类坝体的设计及后期修复治理工作。

图5 法向撞击接触力时间关系

4 结语

通过上述分析，可以得出以下结论：

（1）现场调查结果表明，坝体受大块石的撞击后，直接撞击区受损严重；未受大块石直接冲击的区域，局部出现小的裂缝，整体性受轻微破坏。数值仿真结果与实地调查的坝体损毁特征基本一致。

（2）浆砌块石坝体受大块石的冲击破坏过程是在极短的时间内完成的，为瞬态过程；坝体受冲击破坏的演变过程为：冲击变形—形达到极限状态—损毁。

（3）仿真结果显示，法向接触力的最大值66.31kN，理论计算的冲击力为57.0kN，模拟的最大冲击力较计算值大9.31kN，偏差范围为14.04%，其主要原因是数值仿真模拟计算是一种理想的状态，较相关经验公式计算更与实际吻合。

（4）数值仿真的结果在一定程度上可以指导该类坝体的设计及后期修复治理工作。

［1］田连权，张信宝，吴积善.试论泥石流的形成过程［C］／／泥石流论文集（1）.重庆：科学技术文献出版社重庆分社，1981：54－57.

［2］吴积善，田连权.泥石流及其综合治理［M］.北京：科学出版社，1993：13－15，17－20，182－185.

［3］蔡红刚，裴向军.汶川震区泥石流防护工程损毁特征及破坏机制研究［D］.成都：成都理工大学，2012.

［4］郝好山.ANSYS12.0LS－DYNA非线性有限元分析从入门到精通［M］.北京：机械工业出版社，2010.

［5］白金泽.LS－DYNA3D理论基础与实例分析研究［M］.北京：科学出版社，2005.

［6］Johnson K L.Contact mechanics［M］.2nded.NewYork：Cambridge University Press，1985.

［7］Brizmer V，Kligerman Y，Etsion I.The effect of contact conditions andmaterial propert ies on the elast icityterm inus of a spherical contact［J］.International Journal of Solids and Structures，2006，（43）：5736－5749.

［8］何思明，李新坡，吴永.考虑弹塑性变形的泥石流大块石冲击力计算［J］.岩石力学与工程学报，2007，28（8）：1664－1669.