碾压混凝土重力坝结构非线性有限元分析

侯振宇，崔 炜，刘 毅，朱银邦，邓检强

(1.中国水利水电第十一工程局有限公司，河南 郑州 450001；2.中国水利水电科学研究院，北京 100038)

近年来，中国企业在国外建设了大量水利水电工程，成就了一批有影响意义的业绩[1]。历史原因造成，很多项目应当地政府要求须依照欧洲和美国标准设计。对此，我国正在积极推进中国技术规范“走出去”[2]，以使我们的技术标准被国际广泛接受。本文所研究的海外某水电站由中国企业承建，按美国标准设计，目前处在前期设计阶段。通过三维有限元静力和地震动力分析方法，研究该水电站最高坝段的结构受力情况，复核大坝现有结构方案的安全性，为后期方案的调整和优化提供基础。有限元法是分析重力坝应力、变形及抗滑稳定分析的重要途径[3- 5]。

1 工程概况

某水电站总装机容量为2120MW，最大坝高为130m，参照美国陆军工程兵团发布的文件设计。根据美国标准，该项目亦属于大型项目。挡水建筑为碾压混凝土重力坝，最大坝高在溢流坝段，顶高程为190.00m，建基面高程为60m，坝顶宽度为10m，上游坝坡1∶0.1，下游坝坡1∶0.7。上游正常蓄水位为184.00m，万年一遇洪水位上游/下游为186.60/93.60m，可能最大洪水位上游/下游为189.10/103.50m。

大坝下游朝向为NE52°。基岩以泥岩、砂岩为主，层理面大多为近水平向或轻微倾斜，倾向N和NW，倾角0～10°。河床坝段基岩自上而下为微风化/新鲜的砂岩和泥岩。砂岩厚度约20m，其间夹若干薄层泥岩和粉砂岩。砂岩以下的泥岩厚度超过60m。

2 计算条件分析与确定

2.1 荷载组合分析

根据美国陆军工程兵团重力坝设计手册[6]，大坝应力及稳定需按照常遇(Usual)、不常遇(Unusual)和极端(Extreme)等3种荷载组合所包含的共7种工况计算。其中，常遇荷载组合计算1种工况，即正常运行；不常遇荷载组合计算3种工况：完建，万年一遇洪水位，正常运行遭遇运行基准地震(OBE)；极端荷载组合计算3种工况：施工期遭遇OBE地震，正常运行遭遇最大可信地震(MCE)，可能最大洪水(PMF)。

综合比较手册[6]与中国标准[7]，二者对重力坝结构计算荷载设定的基本原则是相近的。但是，诸多具体规定存在差异，例如：对抗倾覆的要求，对地震荷载的规定，对扬压力的规定和对冰压力的考虑等。美国标准对地震荷载的考虑更多，对坝基帷幕加排水孔的扬压力消减作用持保守态度，忽略了冰压力等。笔者认为中美两种标准无优劣之分，其差异属于学术思想的多元化。实践证明，依照两种标准均可成功建设大坝。

表1 大坝基岩力学参数

本文依据手册[6]的规定施加重力、水压力、扬压力、淤沙压力和地震荷载等。地震动水压力采用动水附加质量的方式施加，依据Westergaard公式[8- 9]。

2.2 材料参数确定与本构模型

根据文献[10]，在早期设计阶段无试验测试的情况下，可依据混凝土静态力学参数推算地震动力状况下的一部分参数。对于文献[10]未给出的经验值，参考中国标准[11]的规定选取。最终确定：在地震条件下混凝土弹模、抗拉强度比静态提高30%，混凝土抗压、抗剪强度不提高。混凝土采用损伤塑性模型[12]，模拟混凝土拉压性能的差异和屈服后的软化、硬化和刚度退化等行为。岩体采用弹塑性本构模型，屈服准则为考虑拉伸破坏的修正Mohr-Coulomb准则。

大坝动力分析的临界阻尼分数取0.03，通过计算大坝第一阶自振频率为3.21Hz。由此得出坝体瑞利刚度阻尼系数为0.00297s。

2.3 模型网格与程序

建立溢流坝段二维有限元计算模型如图1所示，详细模拟了坝体各材料分区及坝基各岩层。模型的x轴指向下游，y轴为坝轴线方向，z轴垂直向上。上、下游边界距离坝体超过140m，底部边界距离建基面超过1倍坝高，共剖分7.5万个节点和6.5万个单元。计算程序采用国际知名的工程力学仿真分析软件ABAQUS。

图1 溢流坝段有限元模型

3 坝体与地基联合作用的静力分析

3.1 坝体位移

大坝完建期，坝体位移总体为垂直向下并略倾向上游，坝顶向上游的位移为15mm，坝基压缩沉降量为12mm。正常蓄水位运行期，在水推力、扬压力等荷载综合作用下，坝体向下游和垂直向上发生位移，如图2所示，坝顶向下游的位移增量为51mm，坝底垂直向上的位移增量为5.8mm；若坝基排水失效，则位移在此值的基础上再增加，坝顶向下游的位移增大约4mm，坝底垂直向上的位移增大约2mm。可能最大洪水工况，坝体下游面受到消力池水推力作用，大坝上、下游水位增长，使坝底扬压力增大。可能最大洪水引起坝体向下游和垂直向上的位移略大于正常蓄水运行工况，两个位移增大约1mm，增长不显著。

图2 正常蓄水位工况坝体位移增量(单位：m)

3.2 坝体应力

大坝完建期，坝体最大压应力发生在坝踵，数值约4MPa；坝体出现拉应力的区域很少，主要位于廊道周围和大坝底板浅层的结构转折处。廊道周围的拉应力为结构尖角处易出现应力集中使然。底板浅层的拉应力值约0.3MPa。水库正常蓄水位运行期，坝体较大压应力发生在下游溢流面下部，该处位于坝体表面且弹模较高使然，最大压应力约5MPa。坝体较大拉应力出现在底板靠近坝踵的浅层(最大拉应力0.9MPa)、闸墩牛腿上游侧(最大拉应力0.6MPa)。坝体上游面的富浆混凝土区有少量拉应力，最大拉应力约0.05MPa。正常蓄水位工况坝体主拉应力情况如图3所示。坝基排水失效对坝体本身的应力影响较小，坝体应力与坝基排水正常时基本一致。可能最大洪水工况，由于闸门敞开，闸墩基本无拉应力，坝体其它部位的应力与正常蓄水位工况差异很小，坝体拉、压应力比正常蓄水运行工况增加0.01～0.05MPa。静力荷载工况的坝体拉压应力小于设计强度值，承载力满足要求。

图3 正常蓄水位工况坝体主拉应力(单位：Pa)

图4 正常蓄水位的坝基塑性应变

3.3 坝基强度状况

大坝完建期，基岩内无塑性区产生。水库正常蓄水位运行期，坝踵下部的基岩浅部有少量的塑性区，如图4(a)所示。这是由于该处岩体受到了坝体位移的影响，产生沿水平略倾斜向上的拉应力。坝基排水失效使扬压力增加，减少了坝体对坝基的压力，坝踵附近的岩体塑性变形略小于排水正常的情况。可能最大洪水工况的库水推力和扬压力同时增长，在综合作用下，坝踵附近的岩体塑性变形的与正常蓄水位运行情况基本相当，差异不明显。静力荷载工况的基岩满足承载力要求。

3.4 坝基深层抗滑稳定情况

有限元法分析坝基抗滑稳定有两种常用的方法，一种是对坝基假定滑面上的法向应力和切向应力进行积分或求和，然后代入抗滑稳定公式，求得安全系数[13- 14]；另一种是坝基强度折减法[15]，保持结构荷载不变，将基岩的强度指标除以折减系数K，逐渐增加K值使基岩强度降低，直至坝基失稳，失稳时的K即抗滑稳定安全系数。K与手册[6]定义的坝基抗滑稳定安全系数FS在力学原理上是一致的。本文采用坝基强度折减法，无需事先假定滑裂面，坝基于抗剪能力最弱的位置率先失稳，塑性变形较大的连通区域可视为滑裂面。

正常蓄水位工况，坝顶、踵、趾的水平位移随基岩强度折减系数K的变化情况见图5，当K大于3.8后位移骤增，且沿坝底以下基岩5～10m的深处延伸至坝趾河床表面出现明显的塑性贯通区，如图4(b)所示。因此，确定坝基深层抗滑稳定安全系数为3.8。

算得其它工况的坝基深层抗滑稳定安全系数如下：正常蓄水位运行遇排水失效工况为3.6，可能最大洪水工况为4.0，均满足手册[6]的要求。这三个工况坝基深层抗滑稳定安全系数的大小关系与针对本项目所做的刚体极限平衡分析一致。

图5 坝体水平位移随基岩强度折减系数K的变化

4 坝体结构地震动力分析

大坝工程区OBE和MCE两种地震的地表峰值加速度为0.034g和0.295g，研究后者。地震动力分析共分三步：第一步向坝体施加重力；第二步向坝体施加正常蓄水位条件下静荷载；第三步向坝基输入水平及垂直两个方向的地震加速度。

印度Koyna水库地震发生在1967年，导致重力坝严重受损[8]。该地震记录时间为10s，水平向地震加速度最大值为0.47g。本文等比例缩小了该地震实测的水平、垂直向加速度，水平加速度峰值缩减至0.295g，时程曲线如图6所示。

图6 水平向地震加速度时程

经计算，在地震过程中，坝顶相对于坝底的最大水平位移为25mm，指向下游，发生在5s左右。坝体经历地震后，坝体混凝土基本完好。仅在上层碾压混凝土分区靠近上游侧和坝底的坝踵浅层出现塑性区(如图7所示)，塑性应变基本小于1×10-4。坝体上游面和其它部位无塑性区。发生塑性应变的上层廊道的上游侧混凝土的拉、压应力变化过程见图8，最大拉应力出现在地震发生的5s左右，最大拉应力约0.9MPa，超过了该处混凝土的抗拉设计强度0.7MPa，最大压应力为2.7MPa。

图7 坝体塑性区

图8 坝体上层廊道的上游侧混凝土应力过程

5 结论

(1)高度最大的溢流坝段静力工况的坝顶最大水平位移约为50mm，坝基沉降量为数毫米。若坝基排水失效或遭遇可能最大洪水，坝体位移增加1～4mm，增长不显著。坝体最大压应力为5MPa，发生在下游面下部。底板靠近坝踵的浅层、闸墩牛腿上游侧拉应力相对较大。坝体承载力满足要求。

(2)坝踵下部的基岩浅部有少量的塑性区，与该处岩体水平向受拉有关。基岩满足承载力要求。强度折减法分析得出坝基深层抗滑稳定安全系数均超过了3.0，满足规范要求，滑动面为沿坝底之下5～10m的深处延伸至坝趾河床表面。

(3)大坝遭遇MCE地震后混凝土基本完好，仅在上层碾压混凝土分区靠近上游侧和坝底的坝踵浅层出现小面积的塑性区。