有限元法介入工程的程度是以往从来没有过的

包承纲

（长江科学院，湖北 武汉 430010）

本文的标题是沈珠江院士生前在参加三峡二期围堰计算工作时常说的一句话。因为对于极具挑战性的三峡二期围堰工程，有限元分析曾起过非常重要的作用。这在上世纪80年代，有限元法刚引入我国不久，岩土研究人员尚在学习和初始应用的阶段，有限元法就大规模在这个重要的工程中实际应用，确实十分难能可贵。当时我们在国家“七五”、“八五”重点科技攻关项目下，组织了全国性的大协作，先后有对岩土数值分析有经验的国内 15家单位、几十位专家参加了这项攻关活动，时间持续十多年，成果相当丰硕，效果十分显著，效益也很可观。这样的大协作在八、九十年代我国的岩土工程界，确是一件难得的盛事！

1983—1985年间，一批改革开放后即出国学习的岩土学者陆续回国参加了三峡工程的研究工作，他们希望采用较先进的技术解决工程中的疑难课题，因此在二期围堰课题中引入了数值分析技术，主要采用有限元法（也用有限差分法等）对设计方案的安全性、可靠性进行论证，并进行方案的进一步优化。二期围堰设计方案的许多决策都是借助数值分析作出的，可以说数值分析在三峡二期围堰建设中的作用是不可或缺的。沈院士的话也是他参加这项工作的切身体会。同时，这项工作对我国岩土数值分析技术的臻进和推广也有实证的作用。

三峡二期围堰的形状比较奇特，这是由特殊的功能要求和当地的地质、施工等建设条件决定的。由于长江主河道流急、水深，为减少截流水深，龙口采用平抛垫底提高江底高程，然后两边端进的截流方案。截流完成后，加高围堰，堰体中部采用花岗岩风化砂填筑，以方便打设防渗墙，而上下游侧则用块石、石渣等粗粒料筑成稜体。围堰在60 m水深中填筑，当填筑到围堰出水（73 m高程）后，再采用干填碾压法填筑30 m到要求的高程。然后，为了防渗在围堰中打两道混凝土防渗墙（深槽部位150 m长），或一道混凝土防渗墙（两侧滩地部位），见图1～2。围堰顶部15 m采用复合土工膜斜墙作为防渗体。由于断面形状奇特，组成材料复杂，填料又不均匀，水上水下堰体的密度差别很大，各部分堰体的变形模量或弹性模量也很悬殊，加之堰体与墙体之间存在相互作用问题，一道墙与两道墙连接处的特殊结构的应力较复杂，复合土工膜这种新材料的工作机理不是很清楚，如此等等都无经验可循。若欲研究其全服役期在受力过程中的性状变化规律，即使只研究二维平面问题，采用传统老办法（如仅核算极限破坏状态的稳定安全系数）是根本不可能的，其安全判别和设计优化也无法进行。此外，有些问题还需借助实地量测才能厘清，如堰体中复合土工膜的变形和应力，不是仅进行实验室测试就能解决问题，因为在堰体中的土工膜会受到周围土体的约束作用，其性状与实验室空气中无约束状态下的情况很不相同。为此，首先要解决测试技术和设备，特别研制了80 mm长的应变片（一般应变片仅几毫米至十几毫米长），而且黏接应变片的胶水还需专门研制，这种胶水干固后的弹性模量宜与土工膜相一致，以使应变片与膜同步变形。

图1 三峡二期围堰深槽段双防渗墙方案Fig. 1 Two cut-off walls in deep canal section of the second stage cofferdam of the Three Gorges project

图2 三峡二期围堰滩地段单防渗墙方案Fig. 2 Single cut-off wall in bottom land of the second stage cofferdam of the Three Gorges project

面对如此复杂的围堰，欲预测它的工作性状，只有另辟蹊径。20世纪七八十年代，电子计算机的应用推动了计算技术的飞速发展，尤其是有限元法等近似计算方法，使科技领域中积存的许多疑难或无解的课题，获得了良好的近似解答。这种计算技术在水利土木领域也受到极大关注。20世纪80年代，美国最高的奥罗维尔土坝第一次采用有限元法进行了坝体应力应变分析，是世界上用现代计算方法进行设计的首例。我国也在这个时期开始在工程中尝试应用，但因经验不足，在工程中应用尚欠“火候”，怀疑者也大有人在。尤其在大型工程中获得实效的先例尚少，某些权威专家一听“有限元”，就连连摇头说：“别跟我谈有限元！”，表示出不屑一顾的姿态。但面对三峡二期围堰的问题，也只有借助数值分析的有限元法，别无它途。于是我们决心闯一闯，不仅为了三峡，也为了我国土力学和土石坝工程的发展。

但我们自知功力不够，难以独当此重任，于是联合国内有志的同行及专家相互学习共同攻关。这项工作正式自 1985年开始，除长江科学院和长江设计院作为牵头单位外，参与二期围堰计算的高等院校、大型科研单位和工程单位先后达15家，参与的教授、副教授或博士以上的学者多达四、五十人，我国著名的土力学计算专家沈珠江院士和多位著名的教授、专家都热心奉献，沈院士还参与总结报告的撰写工作。

从 1985—2005年最后一篇研究反分析的文章发表，前后持续十多年，对分析和优化二期围堰的方案起了不可替代的作用，并发表了许多研究论文，培养了一批博士、硕士研究生。其规模、深度和效果在国际、国内土石坝行业中也是罕见的。通过这次实践，在学术方面也收获颇丰，土的本构模型、土坝有限元计算方法和数值计算前后处理技巧等也有许多经验和成果。这个研究成果曾多次在国内外学术会议或期刊上交流和发表，反映良好。

计算分为基本方案（规定动作）和对比（自选动作）方案两类，基本方案各家都进行，以便对不同单位的成果进行比较，并事先对数学模型、计算方法、参数选择、计算工况、加载减载条件等都作了周密考虑与统一安排，规定采用同一本构模型（Duncan-Zhang E-μ模型）作为基本模型，采用同一网格剖分、一套基本参数、相同的加载方式进行计算。此外，各家还可根据自身特长，自选各自的本构模型，不同的工作条件，考虑不同的因素等设定方案，或作敏感性分析。因此计算成果具可比性，又兼顾各种变化条件，找出最优工况。这样使计算达到了“准确定性，粗略定量”的数值分析目的，成为工程设计和施工的重要依据。

初期的计算成果表明，防渗墙墙体的最大水平位移在1.0～1.5 m量级，变形过大。并且发现，墙体下端部分区域的应力水平超过 1.0，出现了塑性区，且该塑性区上下游贯通，这样围堰的安全性就难以保证，方案不够可靠。为此，曾考虑过两个措施：（1）设法加密堰体，增加墙体后面的支撑。（2）增加墙体本身的刚度和强度，即“固身强体”，这样做本以为会有效果，但相应的有限元分析并未给出预想的答案。因为前者堰体加密只能在墙体前后局部有限的范围（也就5 m左右）进行，且深度最大不能超过30 m，对防渗墙的应力和变形影响很小，同时它要占用直线工期，影响施工安排，故非根本之策；而后者增加墙体的刚度和强度，虽能使应力有所改变，但于墙体变形却无大补益。例如将防渗墙混凝土弹模由18 000 MPa增至20 000 MPa，再增到22 000 MPa甚至更高，墙体变形基本未变。可见，想以“固体强身”来抵抗变形，成效甚微，而墙体上下游贯通的塑性破坏单元虽有所减少但却不能完全消除。这种“怪”现象引起了我们的思索，追溯计算过程（幸好有限元分析还可作这样的追溯），并细加斟酌才发现，原来墙体的变形主要是受庞大的堰体控制，堰体会带着单薄的墙体一起变形，因此不改变堰体密度就于事无补。这个发现真使我们又兴奋又作难，因为如此庞大的堰体是没有办法使其全部加密的。无奈之下只好改变思路，另辟蹊径，可行之路只有一条，就是从改变墙体的材料入手，大幅度降低墙体的刚度，以增加其适应变形的能力。由此，我们的根本思路从“增强”而转为“增柔”，即“柔”的设计思想成为主导思想了。二期围堰这种设计思想的根本转变是在走了弯路并经有限元分析的启发下翻悟出来的。需要注意的是，根据计算，在降低墙体材料刚度的同时还必须保证墙体的一定强度，即墙体的应力要求也需满足。于是，研制一种“低弹高强”的材料就成为当务之急，这又给我们提出了新的难题。

应当指出，上面所有的讨论都是建立在有限元分析成果基础上的，否则不可能有这种新认识的领悟。那么什么样的“低弹高强”材料才能满足工程的要求呢？这也只有仰仗有限元分析才有可能。经过几家反复计算分析得出，若模量控制在1 000 MPa以下，强度R28=4～5 MPa（依墙体高度不同），则墙体的工作状况良好，墙体变形和下部应力均满足要求，围堰可保安全。我们提出以“模强比”（模量与强度的比值）作为墙体材料的控制指标，则其值应在 200～250为宜。这种材料非一般塑性混凝土能够满足，而且当地可用的骨料主要是风化砂，因此，尽快研制出一种以风化砂为主要原材料的柔性墙体材料就是新的重要任务（此点本文不展开）。

有限元的作用不仅体现在二期上游围堰，下游围堰同样亦受其惠。围堰设计中，曾有权威专家对二期下游围堰仅一道1.0 m厚的防渗墙是否过于单薄表示疑虑。当时考虑过两种措施：（1）紧贴墙背增加一道高喷墙，这种办法比较“刚”性，而且要占用3个月直线工期，造价也高。（2）继续走“柔”的思路，适当加厚“柔”性防渗墙。这个决策论证也只有靠有限元。经计算认定，适当加厚墙体厚度至1.20 m，安全度可满足要求，这不仅省工省时，而且还节省好几千万投资，真是多快好省之美事！

三峡二期围堰的分析计算工作，是一个大型工程应用现代计算技术的一次成功实践，增强了同行在土建工程中应用数值分析方法的信心。同时，在计算理论和方法上也摸索了一些经验，为数值计算在工程应用的成熟添了一砖一瓦。例如，当时岩土界比较常用的本构模型是非线性 Duncan-Zhang模型，不少学者认为它与土的性状不完全相符，建议采用弹塑性模型。而国内外发表的弹塑性本构模型虽多，但都缺乏实用经验。为此，参与计算的各家除进行 Duncan-Zhang模型计算外，还自选了南水模型、河海模型、K-G模型、多项式模型等8种模型进行对比。结果表明，弹塑性模型的成果略小于非线性Duncan-Zhang模型，但比较接近，都在合理的范围内，且规律性较好。因此，为弹塑性模型在工程中应用提供了例证。

在防渗墙土石坝的计算中，还有一些不解的疑问，比如墙与堰体的接触关系，因为造墙时有泥浆固璧，墙的周围必存在泥皮，若泥皮很完整地包住墙体，则墙与堰体的接触应属光滑接触，若泥皮仅零星存在，则墙体与风化砂仍属粗糙接触，而后者将大大增加墙体的垂直应力。这个问题影响计算结果，但谁也无法回答。为此，不同计算者只好假定光滑、粗糙和部分光滑等多种情况进行敏感性对比分析，沈院士也曾为此而屡伤脑筋。这个问题直到二期围堰拆除时，将围堰开膛破肚查看才找到答案，表明“泥皮是完整的，厚度可达2 cm左右”，这样就弄清了两者的接触关系，使以后类似计算有所依循。

二期围堰计算工作最有成效的例子，还是在1998年 8次大洪水突然袭来时围堰是否安全这个致命问题的判别上，为回答这个紧急的难题，真真实实是经历了一场真刀真枪的考试，此事以后将另文介绍。

总之，三峡二期围堰数值计算工作，其规模之大、持时之长、方案之多、成果之丰、效果之好，影响之远，在国内岩土工程中确实是空前的，在国外也不多见。这是一次岩土工程与现代科技结合的成功实践，正是：围堰性状难预测，弄清唯靠有限元！