特高拱坝非平衡演化的变形稳定控制理论及应用

程 立，刘耀儒，吕庆超，杨 强

(1.水电水利规划设计总院，北京 100120；2.清华大学水沙科学与水利水电工程国家重点试验室，北京 100084；3.中国电建集团西北勘测设计研究院有限公司，陕西 西安 710065)

0 引 言

特高拱坝是水电开发中重要坝型。近年来，随着我国锦屏一级、溪洛渡等7座坝高超过200 m级的特高拱坝先后建成并蓄水运行，为特高拱坝建设积累了丰富资料和宝贵经验。而在工程实践中，卸荷松弛等勘测设计期未充分重视且常规方法难以定量分析的变形破坏问题日益突出。自然条件下的岩体结构经过漫长的地质演化过程而处于天然的平衡状态，Fairhurst等称之为预先存在的平衡(Pre-existing equilibrium)[1]。开挖、蓄水等剧烈的工程扰动破坏了这种平衡，扰动岩体出现与时间有关的应力和变形重分布过程，该过程可称为非平衡演化，直接影响特高拱坝等岩体工程结构安全建设及长期稳定运行[2]。

近年来，在岩体工程实践中，非平衡演化及其对工程整体稳定性影响研究的紧迫性日益突出。如，贯穿锦屏一级水电站“开挖-浇筑-蓄水”全周期的左岸高边坡时效变形[3]；小湾、二滩、白鹤滩等水电站特高拱坝建基面开挖产生严重的卸荷松弛现象[4]；金川矿区水平巷道变形破坏[2]以及深挖隧洞洞壁大变形、岩爆等，均可认为是工程实践中与时间有关的变形破坏现象。洞室开挖的新奥法施工原理(NATM)是岩体结构非平衡演化的很好例证，陈宗基院士指出：“新奥法原理是一般流变结构演化规律甚至地球动力学演化规律在隧道工程背景下的缩影”[2]。

目前，特高拱坝工程实践中稳定性数值分析主要基于有限元法，但其核心思想来源于极限分析，即超载法和强度折减法。主要的判据可分为2类：①依赖于算法的稳定性判据，如残余力大小[5]等；②经验稳定性判据，如塑性区贯通[6]等。但这些稳定性判据缺乏严格的理论基础。根本原因在于极限分析法采用时间无关的观点研究结构稳定性，认为结构在给定荷载作用下只有稳定或失稳2种可能的最终状态，忽略了结构达到最终状态的时效演化过程，而严格的稳定性理论应基于动态演化观点。此外，动态演化观点对结构长期稳定性具有重要意义，对此极限分析方法很难适用。

针对有限元法缺乏严格的结构稳定性判据这一问题，杨强等[7-9]提出变形加固理论(Deformation Reinforcement Theory)。变形加固理论的发展经历了2个阶段：早期是以弹塑性理论框架和结构失稳定义[7-9]为基础；后期基于稳定演化为基础的时效变形加固理论[10]，最终形成基于Duvaut-Lions粘塑性模型[11]的非算法依赖的结构稳定性理论[12]。并且可以进一步证明，时效变形加固理论的材料本构关系可由Duvaut-Lions模型拓展到非线性、非平衡态热力学框架Rice正则结构，从而适应裂隙岩体的多尺度特点[12]。

1 岩体结构非平衡演化及控制理论

传统的变形加固理论在数值上采用超出屈服面的残余不平衡节点力表征结构的非平衡状态[13-14]，这种静态分析忽略了岩体结构变形的时间效应。为克服上述缺点，杨强、刘耀儒等基于Duvaut-Lions粘塑性模型建立了考虑时效的变形加固理论[10-12]。

1.1 Duvant-Lions粘塑性模型

Duvaut-Lious模型采用率形式本构关系，在考虑关联流动法则情况时，一个变形体在小变形情况时弹-粘塑性材料率形式的应力-应变关系如下

(1)

(2)

1.2 时效的最小塑性余能原理

考虑时效的变形加固理论中，定义E为结构的塑性余能范数，可用下式表示

(3)

(4)

式(4)表明，在恒定外部条件下，理想粘塑性结构始终朝向塑性余能范数减小的方向变形，直到塑性余能范数达到最小值，此后塑性余能范数保持恒定，此即为最小塑性余能原理。

对于任意时刻，塑性余能范数E越大，则表明过应力越大，结构离平衡态越远。换句话说，塑性余能度量了结构当前非平衡演化状态到平衡态的距离(E满足Euclid距离的所有条件)。同时，式(3)表明，E是半正定的，当且仅当结构达到平衡态时(E=0)，也就是结构最终停止变形，达到稳定。式(4)表明，在结构变形过程中，E随时间是不增加的，即结构始终不会背离平衡态演化。结构最终发挥最大的自承载力抵抗外荷载，使其自身尽量趋于平衡态。

塑性余能范数随时间演化见图1。图中，t为时间；Emin为塑性余能范数E的最小值；tls为塑性余能范数达到最小值时的时刻。从图1可知，结构演化分为2个阶段：第一阶段称为自平衡演化阶段；第二阶段称为恒定演化阶段。2个演化阶段的分界点定义为极限恒定演化状态，其对应时刻即为tls。

图1 塑性余能范数随时间演化示意

(5)

基于最小塑性余能原理，塑性余能范数E和不平衡力分布U分别实现了结构非平衡演化的稳定评价与失稳控制。该原理直观表述为：结构通过自身变形达到自承载力最大化、所需加固力最小化的状态，不平衡力即是结构的最优加固力。

(6)

结构失稳控制：阻止结构失稳的加固力为U，从塑性余能范数最小的角度，根据最优不平衡力Umin确定的加固力具有全局最优性，或者说具有最高加固效率。通过有限元分析获得精确到节点的全局最优加固力配置是变形加固理论的最大特点，也是其较快得到工程应用的首要原因。

变形加固理论是基于有限元法的结构稳定性理论，但它并未对有限元法本身进行实质性更改。塑性余能与不平衡力均为有限元法固有中间变量，特别是不平衡力，它就是有限元法中的残余力(residual force)，即增量节点平衡方程的右端项。从计算角度讲，变形加固理论与传统有限元结构分析的唯一区别在于额外计算并记录塑性余能与不平衡力这2项结构响应，它们与位移、应力、塑性区等结构响应的地位是并列的，只是塑性余能与不平衡力反映了结构失稳的根源。

综上，变形加固理论并非更改有限元法，而是深入挖掘有限元法易被忽略的内涵。

1.3 与刚体极限平衡法的对比

对于单滑块单自由度的抗滑稳定问题，变形加固理论与刚体极限平衡法的分析结果是一致的。单滑块稳定分析见图2。图2中，T为滑动力，R为抗滑力，N为法向力，Q为加固力。从图2可知，假定抗滑条件满足Mohr-Coulomb准则，在法向力N的作用下，最大可能的抗滑力R=fN+cA。其中，f和c分别为摩擦系数和粘聚力，A为滑面面积。当T=R，滑块处于临界稳定状态。当T>R时，滑块无法保持平衡，需要施加与滑动方向相反的加固力Q来保持稳定，Q=-(T-R)。

图2 单滑块稳定分析

从变形加固理论的观点分析，当抗滑条件(R≤fN+cA)和平衡条件(T=R)无法同时满足，滑块失去稳定。由于滑面真实抗滑力的最大值为R=fN+cA，因此上述2个条件中，平衡条件不能满足，即T>R。此时作用在滑块上的不平衡力为F，F=T-R=-Q。根据潘家铮“最大最小原理”，此时的滑动面是唯一的滑面，失稳时滑面将发挥最大的抗滑作用。根据牛顿第三定律，不平衡力量值正比于滑块加速度。这就表明，滑块失稳时不平衡力趋于最小，而Q就是最优的加固力。

1.4 不平衡力性质的讨论

不平衡力是弹塑性有限元计算中的一个中间变量，本质是一种虚拟力，并不会在实际工程中出现。如果有限元计算中结构出现不平衡力，则实际工程中会通过损伤、开裂等方式使不平衡力得到释放。根本原因是由于在常用的基于位移有限元分析中，结构连续性和变形协调条件被人为强制满足，无法反映结构局部开裂带来的不连续性和非协调性。理想弹塑性材料是无损材料，出现不平衡力意味着结构无法继续满足在连续、完整、无损的条件下收敛到一个稳定解，此时的不平衡力可从加固和破坏2个角度理解：

图3 预设裂纹的试件在逐步加载过程中开裂破坏

图4 预设裂纹的试件在逐步加载过程中不平衡力矢量

(1)通过施加大小相等、方向相反的加固力来消除不平衡力，进而结构连续、完整、无损性可得以维持。

(2)如不施加加固力，结构只能通过在不平衡力出现处产生损伤开裂(改变本构关系)来消除不平衡力，以此来弱化强制的变形协调条件。因此，不平衡力实质是一种材料损伤开裂程度的等效表述。

以不平衡力表征损伤类似于损伤力学中以有效应力表征损伤，但这种方法在求解过程中只需提供常见的刚度和强度参数，无需额外增加损伤演化相关参数。该方法具有一阶无条件收敛的特性，特别适用于大规模的复杂岩体工程计算。

图3为单轴压缩含预制裂纹的试件开裂破坏过程[15](单轴压力由0逐级增加到2 MPa)。单轴压力为1.6 MPa时有明显的开裂，开裂首先出现在预制裂纹的尖端，随着继续加载，新增裂缝近似垂直于预制裂纹；单轴压力为2.0 MPa时，相邻预制裂纹之间的裂缝贯穿，试件最终破坏。图4为对应的不平衡力矢量图[15]。从图4可知，单轴压力为1 MPa时，不平衡力首先出现预制裂纹的尖端；继续增加荷载，不平衡力扩展如图4的粗实线，与预制裂纹的走向近似垂直。综上，不平衡力的出现和扩展可用于近似模拟新裂缝的出现与扩展，不平衡力量值与裂缝开裂程度近似正相关。

2 特高拱坝建基面卸荷松弛研究

白鹤滩水电站特高拱坝坝高289 m，枢纽区地质条件十分复杂，地应力高，横河向剖面见图6。在枢纽区发育有较多的结构面。其中，LS331、F17、LS3319、LS3318等在拱间槽出露。坝基处于高地应力区，且侧压力系数较高。2014年12月，左岸建基面开挖至628 m高程，巡视发现左岸拱间槽松弛卸荷严重，部分断层出现错动现象，工程紧急停工，进行卸荷松弛研究。

为精细模拟白鹤滩水电站左岸建基面开挖卸荷松弛情况，采用大范围、小单元精细化网格计算。模拟了从950 m高程(天然边坡)连续分步开挖至538 m高程共19步开挖工况，模拟了C3等主要结构面。

图5 白鹤滩水电站左岸坝基开挖剖面

图6 典型剖面开挖不平衡力矢量

建基面开挖至630 m高程，不平衡力较大的结构面为LS331(558 990 kN)、LS3319(234 680 kN)，各断层与错动带的不平衡力和屈服区主要集中在结构面交汇处及邻近边坡或开挖面出露区域。断层F17与错动带C3-1在交汇区相互作用，引起拱间槽下游侧720 m高程附近临空面部分岩体屈服，这与现场观测中拱间槽下游侧坡沿F17产生裂缝相符。LS3319在下游侧545～570 m高程之间有比较大的塑性区，不平衡力也集中于此处。结构面卸荷松弛的现场监测结果与计算结果比较吻合。变形加固理论预测成果表明，从630 m高程继续向下开挖时应重点关注LS331、LS3319、F17的松弛破坏。

图6为从590 m高程开挖至538 m高程过程中，建基面某典型剖面的不平衡力矢量图。计算结果也表明，从590 m高程开挖至570 m高程，LS331的不平衡力由825 690 kN激增到1 179 710 kN，而这也是陡坎成形和LS331开始出露临空的过程。开挖到560 m高程以后，LS331的总不平衡力开始降低。因此，建议将左岸建基面570～590 m高程陡坎坡度放缓，并加强对LS331出露段附近勘测与监测，该建议被建设方采纳。

3 特高拱坝整体稳定性研究

变形加固理论中采用超载法评价拱坝稳定性，对每1个超载倍数K通过非线性迭代计算出相应的塑性余能范数E。以K—E曲线评价拱坝稳定性，可以近似说明岩体结构工程在局部失稳后的整体稳定性情况。引入工程类比法，将某拱坝计算的K—E曲线与之前类似拱坝的K—E曲线进行对比，可以说明拱坝整体稳定性的情况。图7为国内主要高拱坝基础和坝体超载倍数与塑性余能范数关系[16]。

图7 主要高拱坝超载倍数与塑性余能范数关系

基于地质力学模型试验和非线性数值模拟，可以采用3K安全系数对高拱坝进行稳定性评价[16]，这和基于变形加固理论的稳定分析是对应的：

(1)K1为起裂安全系数。变形加固理论中指出，如果外荷载大于结构极限抗力，超出屈服面的应力无法完全转移，即不平衡力不能迭代到0，在不平衡力的作用下，岩体结构进行非平衡演化，出现局部起裂。岩体结构首先出现不平衡力的位置就是首先起裂破坏的位置；不平衡力越大，表明此处的开裂越严重。

(2)K2为非线性变形安全系数。当超载倍数超过一定值时，屈服区体积增加失去线性，认为此时结构整体进入非线性状态。在高拱坝整体进入非线性阶段之前，岩体结构的塑性屈服区体积相比于整个结构的体积很小，可近似认为屈服区的体积与超载倍数为线性关系。所以，可通过超载倍数与屈服体积曲线失去线性关系作为结构整体进入非线性阶段的判据。

(3)K3为极限承载力的安全系数。到达极限承载力后，拱坝全面进入失稳，不能继续超载，结构整体破坏；在非线性有限元计算中，表现为高拱坝的屈服区体积迅速增加，坝踵等薄弱区域发生较大的塑性流动。实践表明，地质力学模型试验中不能继续加载和拱坝塑性屈服区快速增加的2个极限承载力判据是一致的。

4 坝肩高边坡的稳定分析

大岗山水电站坝址区岩脉、挤压破碎带、断层和6组节理裂隙发育。在右岸坝肩区边坡开挖过程中，先后于2009年5月和8月份发生2次开裂，经过紧急支护处理后，边坡变形得到控制。然而2009年9月1日，在进行拱肩槽边坡开挖时，伴随着一声巨响，在右岸1 135 ～1 165 m高程出现裂缝，边坡已开始出现持续的宏观变形，导致工程停工近1年。由于大岗山水电站右岸边坡断层、岩脉密布，且初始地应力较高，开挖所引起的边坡变形破坏问题十分突出。

采用三维非线性有限元对开挖前后大岗山水电站右岸边坡的整体稳定性进行分析，以塑性余能范数来判定整体稳定性，以不平衡力来确定失稳部位及所需加固力。有限元计算网格及模拟的主要结构面见图8。天然边坡与开挖边坡在降强计算中，折减系数—塑性余能范数曲线见图9。开挖前后主要结构面的不平衡力见表1。

图8 计算模型及主要结构面示意

图9 塑性余能随折减系数变化

表1 开挖前后边坡主要断层及岩脉不平衡力对比10 kN

由图9和表1可知，开挖前后边坡的整体安全度均在1.5左右，发生整体失稳的可能性不大。对比开挖前后主要结构面的不平衡力，在整个边坡开挖及降强过程中，缓倾角断层f231不平衡力均有所增加，特别是拱肩槽边坡开挖时，需要在施工时重点关注。而β209不只在开挖过程中产生10 000 kN左右不平衡力，且在降强过程中不平衡力增长迅速，故可以初步判断，开挖过程中的可能破坏模式主要为β209的开裂引起的沿着缓倾角断层f231的局部错动。

5 结 语

处于临界平衡态的岩体结构在受到开挖卸荷、水荷载变化等工程扰动后，平衡态被打破，结构进入非平衡演化过程。在向非平衡演化过程中，会伴随损伤和大变形，直接影响到岩体结构工程的长期安全。时效变形加固理论基于非平衡演化规律的研究，建立了一个能用于结构稳定性评价的理论框架。时效塑性余能范数是评价结构整体偏离平衡态的程度，可用于结构稳定性演化的评价；时效不平衡力是偏离程度在节点上的体现，分布和大小表明了结构破坏的程度，可以用于指导岩体结构工程的加固设计。该理论目前在高拱坝、高边坡等工程的稳定性评价和加固设计中得到了广泛的应用。在特高拱坝建基面开挖卸荷松弛、整体稳定性、长期运行安全性等水工岩石力学的关键问题上，均得到了较好的验证。