《拱坝-地基破坏的数值模型与溃坝仿真》由会员分享，可在线阅读，更多相关《拱坝-地基破坏的数值模型与溃坝仿真（6页珍藏版）》请在金锄头文库上搜索。
1、1拱坝-地基破坏的数值模型与溃坝仿真摘要：本文使用刚体弹簧元模型来模拟拱坝-地基这类非连续介质和连续介质组成的系统，对一个悬臂梁的位移场和应力场的计算并与有限元方法计算结果的对比表明，这种方法具有良好的精度.应用这一模型对玛尔帕塞拱坝的溃坏进行了仿真，重现了拱坝的破坏过程，探讨了拱坝坝肩滑动失稳的原因，计算结果表明坝肩岩体裂隙和断层的强度不足是导致整个拱坝溃毁的关键因素. 关键词：三维刚体弹簧元 拱坝地基破坏 拱坝 仿真 近年来，用非连续介质力学模型研究高坝地基、高边坡与地下工程的变形稳定有长足进展，离散元14、DDA5、刚体弹簧元6等的发展提供了求解非连续介质力学问题的有效工具.关于高坝-地
2、基系统的破坏机理与仿真分析，由于问题的极端复杂性，迄今仍未见诸报导，本文利用改进的刚体弹簧元研究玛尔帕塞拱坝溃坏的机理与过程，是试图用连续-非连续介质统一的数值模型研究高坝-地基系统破坏过程的尝试. 法国玛尔帕塞拱坝 1959 年 12 月 2 日的溃坝事故是坝工史上的重大事件7,8，引起了坝工界的极大重视.世界各国专家和学者从未停止过对玛尔帕塞拱坝失事原因的探究，并提出了各自见解9，归纳起来有以下 3 种看法：(1)坝基变形，认为左岸坝基修建在弱风化顶层上，岩体弹性模量低，仅 1000MPa 左右，地基的过大变形使拱的推力向上向下转移，其中向上转移的推力使左岸重力墩超载而产生位移，导致拱坝失
3、去支撑而破坏.(2)坝体上滑失稳，认为岩体浅层的裂隙在不利的条件2下不能提供足够的抗滑力，致使坝体沿建基面向上和向下游滑动，拱圈拱弦被拉长，使坝沿一岸或两岸拱端转动，导致拱冠断裂，坝体崩毁.(3)坝肩岩体滑动，认为左岸坝基存在一个潜在的由下游断层、发育的节理面，以及下游临空面组成的滑动楔形体.由于地质条件在楔形岩体上游侧形成很高的静水扬压力使岩体滑移失稳.本文利用改进的三维刚体弹簧元对玛尔帕塞拱坝的失稳破坏过程进行仿真，对该坝的失稳机理作进一步的探讨.1 拱坝破坏仿真模型1.1 三维刚体弹簧元模型 假定三维刚体弹簧元的单元块体为不可变形的刚性体，在引入了刚度、强度等效准则后，将块体变形与强度指
4、标统一集中于接触弹簧上.块体的运动符合刚体的平动和转动定律，块体与块体依靠相互之间连接的切向和法向弹簧和阻尼器来传递相互之间的作用，如图 1.弹簧的作用力由块体间的相对位移确定.如图 2，块体 P 对块体 Q 的相对位移为 u，则他们之间的弹簧力变化量为式中：Ks、Kn 分别为切向和法向弹簧的刚度；s、n 分别为相对位移增量的切向和法向分量；s、n 分别为切向和法向弹簧力的变化量.图 1 接触模型图 2 相对位移增量示意 1.2 使用三维刚体弹簧元模拟非连续介质 三维刚体弹簧元模型基于离散元方法提出，是非连续介质模型的 1 种，由于块体之间可以错动、分离，因此可以方便地模拟非连续体如裂隙发育岩
5、体的大变形、大位移问题.离散元方法中接触关系一般可分为角面、角边、角角和边边四类，关于各类接触关系的检索和其它离散元的详3细内容已在 Cundall 和鲁军等人的文章中有详细叙述，参见文献14.1.3 使用三维刚体弹簧元模拟连续介质图 3 刚体弹簧元模拟连续介质原理示意图 4 中心受拉杆件1.3.1 基本原理 如图 3，通过将连续介质的变形和应力凝聚到单元之间连接的法向和切向弹簧上，可以用刚体弹簧元模型来模拟连续介质.即由弹簧的变形刚度(Kn、Ks)来等效连续介质的变形模量(E、G)，由弹簧力系统(Fn、Fs)来推求与之静力等效的连续介质的应力(、)分布，从而将刚体弹簧元系统和连续介质联系起来
6、.(1)刚度等效.如图 4 所示，一连续的 ll1l2 方杆，弹性模量为 E，泊松比为 .将杆沿长度方向等分成 n 个 l3l1l2 的块体，其中 l3=ln，每个块体间的面面接触被视为四个角点处的点面接触，当 n 足够大，取弹簧的法向刚度 Kn=l1l2E/4l3 时，使用三维刚体弹簧元计算它在中心作用力 P 下的变形与使用材料力学方法的结果是一致的.弹簧的切向刚度借用弹性模量和剪切模量的关系可取 Ks=Kn/2(1+).单元应尽可能采用长方体，对于不规则的块体可以将其均化成长方体使用上述估算式，或另外推求相应的计算式.(2)强度等效.如果刚体弹簧元网格划分得比较密，则可以取单元面上每个弹簧
7、力合力除以该面的面积求得其平均应力作为该面上的应力值.因为单元网格较密，可以认为当面上的平均应力达到材料的强度时该面上的弹簧全部断开，则单元在该面发生破坏.由于刚体弹簧元的计算工作量比较大，当为了节省计算工作量，网格划分得较稀疏时，可以按照材料力学方法假定单元面上的应力分布方式，由弹簧力根据静力等效原则求出该4面上的应力分布.(3)荷载等效.对于集中荷载，将荷载作用在相应的块体形心上，并将荷载对中心的力矩加入到该块体的合力矩中即与原荷载的作用等效.对于分布荷载也要等效成集中荷载加入到相应的块体中.图 5 简单算例1.3.2 简单算例 如图 5，一长 40m，高 4m，宽 1m 的悬臂梁受均匀荷
8、载 P=105N/m2作用，支座 A 固定，不考虑自重，坐标见图中所示.悬臂梁弹性模量 E=1010N/m2，泊松比 =03，密度 =2.4103kg/m3，离散网格如图 5(在厚度方向为一层，厚 1m).使用上述的刚体弹簧元方法计算悬臂梁的变形和应力并与有限元的结果进行对比，图 6 为悬臂梁沿 y 轴的位移对比结果，图 7 和图 8 分别为 x 方向的正应力 x 和剪应力 xy 的对比结果.从图中可以看出两种方法计算的位移和应力结果除有限元法在固定端能反应局部应力集中外，整体应力具有很好的近似性.图 6 y 方向位移等值线(单位：m) 图 7 x 方向正应力 x 等值线(单位：MPa) 图
9、8 剪应力 xy 等值线(单位：MPa)刚体弹簧元通过引入单元表面应力分布假设和变形近似等效来模拟连续介质，利用这些假设，刚体弹簧元通过上述的刚度等效、强度等效和荷载等效，能够把连续介质和非连续介质统一在一个模型之中，为研究坝肩失稳导致拱坝断裂破坏需要仿真连续与非连续耦合介质的大变形与破坏断裂的复杂过程提供了有力的工具.2 玛尔帕塞拱坝失稳机理的研究52.1 模型的建立2.1.1 坝肩和坝体的离散和计算条件的设定 法国玛尔帕塞拱坝坝高 66m，坝顶高程 102.55m，坝顶长 222.7m，长高比 3.3，坝顶厚 1.5m，底厚 6.78m，厚高比为 0.1.坝址岩体由带状片麻岩组成，含千节理
10、倾向下游右岸，倾角 3050.该坝于 1954年建成，蓄水历时 5 年至 1959 年 12 月蓄至最高水位约 100m 时，坝在短瞬间溃决，左岸与中部坝体全部冲走，右岸底部坝体残存.根据玛尔帕塞拱坝的地形和地质资料以及拱坝的体形和横缝设置将坝址处的岸坡简化成如图 9 所示，坐标系设置和拱坝离散网格如图 10.在分析中边界条件设定为除“二面沟”内的岩体外，其他岩体均被固定.坝体、右岸坝肩和中间岩体依照连续介质进行模拟，坝体弹性模量取 E=1010Pa，坝体材料的抗拉强度取 2.0MPa，抗压强度为 20MPa，摩擦系数f=1.0；考虑到右岸坝基事后调查仍有坝体残存，表明该区岩体强度较高，故将坝
11、体右岸及中间部分与基础的交界面强度较坝体提高 1 倍.左岸坝肩依照非连续体计算，重力墩和地基之间的摩擦系数设定为 1.0，坝体和左岸坝肩之间以及“二面沟”内的岩体两组节理面与水平层面取相同的摩擦系数 f1，根据不同的工况选取，上、下游断层摩擦系数为 f2，亦根据不同的工况进行选择.计算中未考虑粘结强度 c 的作用.主要荷载为库水作用，稳定分析中库水水位均取 100m，考虑左坝肩扬压力时，按静水压力作用施加在沿上游面坝基线的裂隙处，作用面见图 9 中虚线段所示.为了简化计算，所有荷载均一次施加.计算中阻尼参照 Cundall 提供的阻尼模式选取刚度与质量成比例的 Rayleigh 阻尼型式，阻尼
12、系数 、 由试算确定，试算表明：如不出现过阻尼或阻尼过小情况计算结果均收敛于一稳态解的近似值，详见参考文献3.2.1.2 位移分析 计算拱坝位移时左坝肩暂按连续介质考虑，并且不考虑扬压力作用，得到拱坝在 1955 年 9 月(水位 79.75m)蓄水至 1959 年 7 月(水位 94.19m)拱坝6产生的位移，并与相应实测位移值进行比较，图 11 列出了顶拱和拱冠梁处的比较结果.可以看出二者最大值基本接近，但在拱冠梁坝基附近，实测的位移值较计算大 1cm 左右，顶拱计算结果较实测结果略偏大.这一比较结果表明选择的弹簧刚度参数 Kn、Ks 以及单元等效模型基本符合实际.2.2 玛尔帕塞拱坝稳定
13、性分析及失稳过程的仿真 通过改变坝体和左岸坝肩之间以及左岸坝肩“二面沟”内岩体的摩擦系数 f1 和上下游断层面的摩擦系数 f2 可以模拟和分析拱坝在各类参数条件下失稳破坏的情形.图 9 岸坡离散图 10 刚体弹簧元离散网格图 11 拱坝位移计算结果与实测结果比较2.2.1 工况 1不考虑扬压力的浅层滑动 不考虑左岸坝肩处的扬压力，将坝体和左岸坝肩之间以及左岸坝肩“二面沟”岩体的摩擦系数 f1 和上、下游断层面的摩擦系数 f2 设为相同，从 1.0 开始折减，逐步计算，直到拱坝发生破坏.计算结果表明当摩擦系数 f1 和 f2 折减到 0.6 时坝肩失稳，拱坝破坏.图 12 列出了左坝肩附近的块体
14、 1(见图 10)的位移过程.f=0.7 时，块体位移量陡增，接近 30cm，此时拱坝处于失稳的临界状态，f=0.6 时，块体位移不再收敛，开始沿着坝肩向下游、向上滑动，整个拱坝随之破坏.图 12 块体 1 的位移过程(工况 1) 图 13 f=0.6 时左坝肩重力墩位移过程(工况 1)7分析图 13图 15，拱坝的整个破坏过程可以分为两个阶段，第一阶段是大约 3s以前拱坝整体位移，第二阶段是 3s 以后由于过大的变位使拱坝的工作方式发生改变，应力恶化，坝体溃毁.描述如下：在第一阶段，拱坝利用拱的作用将水压力转化为水平作用在坝肩上的推力，沿坝肩向上的分量使拱圈产生向上滑动的趋势，由于左坝肩岩体
15、的摩擦系数较小，不能提供足够的抗力，巨大的拱推力转移到了左岸的重力墩上，使重力墩沿着墩轴线偏向下游移动了近 2m(见图 13，3s 前重力墩主要以平动为主)，同时整个坝体沿建基面向下游、向上滑动，带动坝体以右岸为轴向下游转动(由图 14 中 y 方向的位移分布，坝体向下游位移由右到左依次增大，显示以右坝肩为轴向下游转动)和向上旋转(见图 14 中 z 方向的位移分布，坝体上抬量由右到左依次增大，显示以右坝肩为轴向上转动).此时，整个拱坝基本上保持整体发生变位，同时有沿右岸转动趋势.(a)y 方向位移等值线(b)z 方向位移等值线 图 14 f=0.6，3s 时拱坝坝体位移分布(工况 1，单位：
16、m)拱坝在经过了第一阶段后，拱圈拱弦拉长，逐渐被压扁，承载能力极度降低，应力发生恶化，坝体破裂，拱坝的破坏进入到第二阶段破坏阶段.由于中、下部块体受水压力较大，拱坝坝体从该部位偏左岸发生断裂，向下游弯折，首先被冲出，随后其周围坝块跟着折向下游，上部坝块受水压力较小，发生翻转后随在中下部坝块之后冲出，拱坝彻底破坏.与重力墩连接的坝段在被向下游冲出的过程中又带动重力墩发生了向下游的位移和旋转(见图 13 中 3s 以后的位移，向下游移动 1m左右，在平面上旋转约 3).图 15 给出了上述过程8一些典型时刻拱坝的状态，拱坝破坏后最后残留状态见最后一幅.图 15 f=0.6 时拱坝破坏过程(工况 1)上述拱坝破坏过程的计算结果和拱坝失事后的调查结果基本相符：(1)f 值的大小，“实有摩擦角为 3035”，本文计算表明 f 值取 0.7 时拱坝处于临界状态，对应的摩擦角为 35，当拱坝左坝肩的实有平均摩擦角小于该值时，拱坝会发生上滑失稳，与拱坝破坏事实相符；(2)坝体的残留部分，拱坝
