（工程力学专业论文）基于非线性有限元的高拱坝整体稳定分析.pdf

摘要 随着国民经济的飞速发展和西部大开发战略的实旌，我国的西南和西北正在修建和即 将修建一批3 0 0 级的特高拱坝，如小湾、溪洛渡、锦屏一级、白鹤滩等。这些： 程坝高库 大，环境条件复杂，又处于高烈度地震区，因此，高拱坝的破坏模式研究与安全评估极其 重要。本文基于非线性有限元法，探讨使用综合法研究高拱坝整体稳定安全度。论文主要 工作如下： ( 1 ) 介绍了水工结构工程有限元分析中非线性计算的基本理论与方法，以及在 a b a q u s 中的实现手段。对拱坝整体稳定性问题进行了研究，并分析了综合法确定整体安 全度的有限元有效加载模式。 ( 2 ) 对于一个含有软弱夹层的重力坝算例，进行了整体安全度分析。在初始折减系 数分别取值为1 2 、1 2 5 、1 3 、1 3 5 、1 4 的条件下，采用综合法分析各自状况对应下的整 体安全系数，并比较基于不同安全判据下的安全系数的区别。提出以塑性应变能突变判据 为主，其余判据为辅的建议。 ( 3 ) 利用a b a 0 u s 对锦屏一级高拱坝进行了非线性有限元分析，利用综合法分析了 该拱坝的整体稳定系数，得出该拱坝的整体安全系数为4 5 5 5 1 7 。与已有实验结果进行比 较，发现研究成果合理，可为以后类似工程设计提供一定的参考。 关键词：拱坝：非线性有限元；稳定性；综合法；整体安全度 a b s t r a c t a c c o r d m g 诵mt l l el l i g m yd e v d o p m c n to f 甘1 en a t i o n a le c o n o m y 触l dt l l e 、v e s t e m d e v e l o p m e n t ，t l l e r ea r em a n yh i 曲a r c hd 锄sv v i t h3 0 0 mh c i 曲ta r eb u i l tn o w o rl a t e r ，s u c ha s x i a o w 抽、x i l u o d u 、j i n p i n g 、b a i h 咖n n l e s ep r o j e c t sa r ea l w 甜sw i 血h i 曲d 锄、b i gr e s e “o 址 锄dc o m p l e xe l w i 咖e n tc o n o n s ，e s p c c i a l l yi i lt l l eh i g hi n t e n s i t ye a r t h q u 罟出ea r e a s o ( h e f 缸i u r em o d ea n ds a f e t ya s s e s s m e n to ft h el l i g ha r c hd 锄i sv e 巧i m p o r t a n t t h e9 1 0 b a ls 撕 d e g r o fl l j 曲a r c hd a n lw a sd i s c u s s e db ys y n t h 鲥c 鲥m e t h o do nn o n “n e a rf e m t h em a j o r c o n t e m so ft l l i sp a p e ra r ef o l l o w s ： ( 1 ) 7 1 1 1 eb a s i ct l l e o 巧a n dm ea p p l i c a t i o no nh y d 咖l i cs t n 蜘r ce n g i n e e r i n go fn 舳l i n e a r f i n i t ee l e m e n tm e t h o da r es u m m e du p ，甜l dt 1 1 em e 柚so fr e a l i z a t i o n 抽a b a q u sa r ej n 协) d u c e d i nd e t a i j t h er e s e a r c ho nt h ew h o j cs 协b j j i 哆o f 盯c hd 姗h 船b e e nd o n e t h ee l j f e c t i v el o a d p a t t e mt od e t e m l i n eg l o b a ls a f e t yd e g r e eb ys y n t l l e t i c a lm e t h o dw a s s t u d i e d ( 2 ) ag r a v 畸d 锄w i t h1 1 7 l 忙a k 硫e r l a y e r si ss n l d i e da sa i l “啪p l e ，a n di t sg l o b a ls a f 呻 d 叼雠i sa n a l y z e d n e9 1 0 b a ld e g r e e sm a t c h 。dt od i 疵r e mi n i t i a lr c d u c t i o nc o e 伍c i e n t sw h i c h 1 2 、1 2 5 、1 3 、1 3 5 、i 4 ，w e r cc a l c u l a t c db ys y n t h e t i c a lm e t h o d ；柚dd i 蔬r e n tr e s u l t su n d e r d i f r e r e ms a f 毫锣c r i t e r i o n sw a sc o m p a r e d ，as a 岔够c r i t e r i o n 、7 燃a d v i s e d ：c a t a s t r o p h eo c c u r r e n c e s o fp l a s 矗cs t m i ne n e r g y 州m 撕l y ，趾di n e a n ，h i l em 觚yk i n d so fs a f e t yc r i 州o nt a k e ni n t o a c c 0 嘶 ( 3 ) n o l l l i n e a rf e m 锄a l y s i so nj i n p i n ga r c hd 锄i sc a 而e do u tb ya b a q u s n er e s u l t t h a tt l l e 百o b a ls 如哆d e g r e ei s4 5 5 5 17w a sc a c u l a t e db ys y n t h e t i c a lm c 廿1 0 d ，a n da c c o r d i n g w i me x p 耐m e mr e s u l t sw e l l ，i e m es y n t h e t i c 甜n l e m o da s s e s s m gt h eg l o b a ls a f e t yd e g r e eo f h i 曲a r c hd 锄i sr e a s o n a b l e ，a n dc a nb eu s e di ns i m i l a rp m j c c t s k e y w o r d s ： a r c hd 枷；n 砌i n e a rf e m ；s t a b i l 峨驴t h e t i c a lm e m o d ；g l o b a ls a f c t yd e 舻e j l 学位论文独创性声明： 本人所呈交的学位论文是我个人在导师指导下进行的研究工作及取得的 研究成果。尽我所知，除了文中特别加以标注和致谢的地方外，论文中不包 含其他人已经发表或撰写过的研究成果。与我一同工作的同事对本研究所做 的任何贡献均己在论文中作了明确的说明并表示了谢意。如不实，本人负全 部责任。 论文作者c 签孙j 选册月和 学位论文使用授权说明 河海大学、中国科学技术信息研究所、国家图书馆、中国学术期刊( 光 盘版) 电子杂志社有权保留本人所送交学位论文的复印件或电子文档，可以 采用影印、缩印或其他复制手段保存论文。本人电子文档的内容和纸质论文 的内容相一致。除在保密期内的保密论文外，允许论文被查阅和借阅。论文 全部或部分内容的公布( 包括刊登) 授权河海大学研究生院办理。 论文作者( 签名) 粘啦囊谬年5 月呷日 第一章绪论 1 1 问题的提出及研究背景 第一章绪论 随着经济的发展和对环境要求的提高，同时为了缓解西部乃至全国电力的紧张状态， 作为绿色资源的水电资源受到了人们极大的重视。为了让这些水电资源得到较好的运用， 我国决定在西部建立一批混凝土拱坝来实现这个目的。拱坝的使用历史最早可以追溯到古 罗马人利用功的原理修建高十余米的圆形圬工结构。一般来说，拱坝具有以下优点【l 】：受 力条件好，拱坝同时起拱的作用和悬臂梁的作用，受力比较合理；坝的体积小，相对于同 样高度的重力坝，拱坝所需的混凝土节省1 3 2 3 ；超载能力强，安全度高，拱坝属周边嵌 固的高次超静定结构，能进行坝内应力重分配，在合适的地形地质条件下具有很强的超载 能力；抗震性能好且能大量泄流。由于以上各种优点，使得拱坝在造福人类的历史长河上 得以不断发展与完善。 虽然实践中证明拱坝是较安全的一种挡水建筑物，但由于拱坝结构及地质条件的复杂 性等因素的影响，拱坝建设中仍然存在一些技术上不确定因素与认识上的片面性。若不加 以重视，不仅会给国家造成巨大的经济损失，还会给下游居民的生命财产安全带来巨大损 失，造成不良的社会影响。1 9 5 9 年法国的m a l p a s s e t 拱坝因坝体连同基岩发生深层滑动而 溃坝，导致4 0 0 多人死亡及重大经济损失；1 9 6 3 年意大利瓦依昂( v 萄o n t ) 双曲薄拱坝库 内左岸滑坡体突然下塌，导致3 0 0 0 人死亡。据不完全统计，至1 9 8 0 年止，国外已建拱坝 有4 8 座出现问题或破坏【2 ，7 2 1 。在我国西部地区正在修建和即将要兴建一批世界最高级别拱 坝，高度都达3 0 0 米左右，为了让这批拱坝得以更好的造福人类，以及解决高拱坝建设中 的问题，进行拱坝的安全度研究是十分必要的，并且在我国的连续三个五年计划中，高拱 坝的安全研究也占着重要的一席之地【3 】。 拱坝是一个由坝体、地基和建基面组成的三维超静定系统，但是我们使用现有的方法 对拱坝进行安全度校核的时候，往往将坝肩岩体的稳定性与坝体强度的计算和校核分开进 行。实际上，在坝体横缝开张的时候，坝肩交接面是坝体薄弱部位，但该部位可能远在坝 肩岩体作整体滑移前，由于其局部开裂、错动引起的坝肩变形就已导致该面的严重开裂， 坝体已经由于拱座过大变形而破坏，因此坝肩岩体的抗滑稳定安全系数满足设计标准并不 一定能保证坝体的安全，反之亦然。为此，研究高拱坝的破坏与安全，将拱坝强度安全判 别与坝肩岩体稳定判别进行有机的结合，考虑地基内的不良地质构造，以及坝体与地基的 l 河海大学硕士学位论文 相互作用，将坝体和包括关键潜在滑移岩块的坝肩岩体作为综合体系，建立坝体薄弱部位 开裂的数学模型和求解方法，才能切实反映拱坝整体失稳状态。 本论文主要结合位于四川省凉山彝族自治州盐源县和木里县境内的雅砻江干流上的 重要梯级电站锦屏一级水电站开展研究。该水电站以发电为主，兼有防洪、拦沙等作 用；拱坝坝高3 0 5 米，是目前世界上最高的拱坝，水压力巨大，其设计和分析都是超出现 行规范的。而且，坝址区工程地质条件复杂。两岸边坡高度近千米，虽然基岩裸露，风化 卸荷不强，岩体坚硬完整，但左岸存在深部裂缝和f 5 、f 8 ，右岸存在f 1 3 、f 1 4 等断层， 给拱坝坝身的工作性态和坝基岩体的变形、渗透稳定带来不利的影响。坝址右岸为大理岩 顺向坡，其中等倾角的大理岩层面断续分布有绿片岩透镜体，形成层间弱面，对拱坝右坝 肩的稳定不利。这些都使锦屏高拱坝的安全性，包括坝身的强度和稳定、坝肩岩体的稳定 性，成为人们极为关注的问题，也是锦屏拱坝建设中最重要的关键技术问题之一。 1 2 拱坝整体安全度研究现状 1 2 1 拱坝的失稳形式 近1 0 多年来国内外学者逐步重视高拱坝破坏机理的研究，积极探索与整体安全度相 适应的研究方法。当坝体由于各种原因处于极限平衡状态时，只要有一个小的扰动，它将 偏离原平衡状态而产生无法恢复性的破坏。根据运动稳定性定义：任何状态或事物的变化 都是一种运动，都存在是否稳定的问题。鉴于高拱坝整体破坏时状态从静止向可动状态转 变，一些学者认为高拱坝整体破坏问题就是失稳问题【4 l 。由实验可知，拱坝的破坏，如图 1 1 所示，基本可以分为四种【5 】：( a ) 坝体结构开裂，坝肩仍保持稳定；( b ) 坝肩局部开裂， 在渗水的挤压中出现崩裂，导致坝体结构破坏；( c ) 坝肩出现滑移体，引起整个坝体结构破 坏；( d ) 坝基深部软弱带变形过大，坝体拱端出现大变形，从而导致拱坝结构破坏。后三种 破坏均为稳定破坏，在研究拱坝的稳定性的时候，单独考虑任何一种破坏方式均有一定的 片面性，因此我们需要综合考虑各种可能破坏形态，且因为稳定理论可以解决高拱坝的整 体破坏问题，高拱坝稳定性可用整体安全度衡量【6 】。此外，沿建基面的滑移有时也被考虑 为拱坝的一种失稳方式，但是，对于双曲高拱坝，沿建基面的上滑失稳的可能性较低【7 1 ， 且该破坏方式一般伴随着坝体的破坏同时发生1 6 】，因此该失稳方式一般不单独作为一种拱 坝失稳方式加以研究。 2 第一章绪论 ( a ) ( b )( c ) 图1 1 拱坝的几种破坏形式 1 2 2 整体稳定性分析方法综述 从整体稳定的角度考虑拱坝的安全度的一般分析方法主要分为以下两种： ( 1 ) 地质力学模型试验法 地质力学模型试验法就是从力学的观点出发，采用试验的手段，研究地质构造条件对 工程的影响。地质力学模型试验又称为地力学模型试验或岩石力学模型试验【8 t9 】。该方法 可以较全面的模拟不连续岩体的构造、软弱结构面和断层破碎带等自然条件，以及岩体自 重、强度、变形模量、抗剪指标等岩石力学特性。通过试验可了解工程从加载开始直至破 坏的整个过程和破坏机理、拱坝的超载能力和变形特性、裂缝的分布规律和需要加固的部 位以及地基处理的效果等uo 】。它着重研究结构超出弹性范围以外的性能，故仅考虑一次加 荷效应；但它不仅限于已知荷载条件下的某一状态，更重要的是研究在渐增荷载作用下直 至破坏的整个变化过程，这就是地质力学模型与其它模型的区别所在。而且，随着试验量 测技术的发展，地质力学模型试验中的研究方向，已由最初的简单模拟实现的研究向新的 模型材料的研发、新的技术和加载方法甚至全新的实验理论研究方面发展【l l 】。 早在2 0 世纪6 0 年代，f 啪a g a l l ie 等专家在意大利贝加莫结构模型试验所( i s m e s ) 已 经成功地进行了多项地质力学模型试验【1 2 】。自2 0 世纪7 0 年代以来，该方法得到了广泛的 应用，扩大了结构模型试验研究的领域，世界一些国家和地区建立起了结构模型实验室， 对各自国家的重点工程也分别开展模型试验的研究，苏联的英古里、契尔克，日木的黑部 川、殿山等拱坝都进行过脆性材料模型破坏试验【1 3 】。我国在2 0 世纪7 0 年代中后期，也开 展了这方面的研究工作【1 4 】，主要代表试验为：张林等【1 5 t1 6 1 、姜小兰等【1 7 1 、李朝国等【1 8 】分 别通过试验方法对锦屏、溪洛渡、二滩等拱坝做了地质力学试验研究。这些试验结果都证 明地质力学模型试验是一种分析研究整体稳定安全度的有效方法。 由于拱坝坝肩周围岩体本身的复杂性，使其尺寸、材料类型及性质均变化较大，使得 室内模型试验的各种应力比尺、几何尺寸、荷载比尺和力学性质比尺等难以统一、相容； 河海大学硕士学位论文 并且模型试验成本昂贵，地震、扬压力、温度等荷载的模拟不可能与现实完全相符，其测 得的应力、变形有一定的误差，难以进行多种工况的研究。现场试验一般耗费巨资，且受 地形、地质、施工条件限制，测试工作往往只能在一个很小的范围内进行，并且小范围的 测试成果很难用来代表复杂的大范围的工程岩体，因而实验结果有一定的局限性l l 9 。 ( 2 ) 数值计算法 数值计算法自6 0 年代引入稳定性分析中，它主要包括边界元法【2 0 】、离散元法、刚体 弹簧元法【2 、有限元法及块体元法【2 2 1 等。 有限元法具有广泛的适用范围，可以对任何形状的几何体进行线性或者非线性的分 析。该方法通过对坝体基岩的应力、变形分析，进而分析工作岩体及坝体结构的稳定性。 有限元法从5 0 年代中期开始提出，在6 0 年代初期就己广泛应用于岩体力学中。克劳 夫、威尔逊、辛克维茨等人首先采用线弹性连续体有限元分析方法来求解岩体力学问题。 1 9 6 8 年古德曼等人提出岩体不连续面单元，这对模拟岩体中的断层、裂隙、节理等结构面 的存在有很重要的贡献。1 9 6 8 年，厄加塔迪斯等又采用有限元方法对拱坝及基础作了三维 分析。继1 9 6 3 年威尔逊首先提出应用非线性有限元进行结构分析计算，以后阿尔杰瑞斯、 德塞、马可尔等人进一步发展，邓肯、古德曼等人又相继研究了有节理岩石的边坡问题1 5 j 。 随着时代的发展，有限元法在各地科研人们的努力下得到了不断的扩充和完善，已日趋完 美和成熟。陈胜宏等【2 3 1 在开裂单元中引入不连续形函数，然后基于虚功原理推导出开裂单 元的广义刚度矩阵和广义荷载的具体表达式，以及有限元平衡方程，从而建立了一种三维 裂缝扩展的不变网格有限元分析方法。杨强等【2 4 】采用工程结构极限分析方法，利用已有的 有限元计算程序并结合数学规划方法来推求拱坝的整体安全系数，提出了适用于三维混凝 土高拱坝结构极限分析的有限元数学规划算法。基于t j r h u 曲e s 【2 5 1 在1 9 8 3 年提出有限 元e b e ( e m e n t - b v e l e m e n t ) 思想，刘耀儒等【2 6 1 将其运用于拱坝的三维有限元分析中。章 青等【2 7 1 从摩擦理论和矢量几何出发，导出了三维可能滑体的抗滑稳定安全系数计算公式， 并通过等参逆变换的方法，将非线性有限元法的计算成果纳入到所建立的公式中，形成了 “基于非线性有限元的高拱坝坝肩岩体抗滑稳定分析方法”。由于计算机水平的发展，自 1 9 7 2 年国内第一个拱坝计算机应用程序的发展以来，现在最新的拱坝c a d c 系统已经可 以为大型、综合性拱坝计算机辅助施工设计、施工放样、施工模拟、计算机图形处理系统， 适应各种型式拱坝设计、计算、施工和绘图的需要。同时包括尾水管与蜗壳的c a d c 与 c a d m 【2 8 l 。 按照研究方法的不同，有限元进行稳定分析的方法可以分为直接法和间接法两种【2 9 】： 4 第一章绪论 直接法是根据滑移面上的应力计算沿滑移面的稳定安全系数。具体做法为：在正常情 况( 正常荷载、正常参数) 下，计算系统的位移场和应力场，假定某个面是破坏面，由以 下极限平衡公式直接计算安全度： k ：蜂：竺! 竺( 1 - 1 ) n c o s 肚 式中：吒和f ，分别为滑移面内某点的法向正应力和切向剪应力，厂、c 是该点的摩擦系数、 凝聚力，口和分别为该点运动方向与整体滑移方向之间的夹角。 间接法是模拟某种破坏模式，使拱坝与地基整个系统达到实际破坏前的临界状态，通 过它与正常状态之间的比较确定其安全度。 目前常用的破坏分析是通过超载或强度折减来实现的，也就是通过超载或者改变材料 抗剪强度的方法使系统达到极限平衡状态，根据结构达到极限平衡状态时的荷载改变倍数 或者抗剪强度改变倍数来确定稳定安全度，前者为超载安全度，后者是强度储备安全度。 直接法只适用于滑移面较为简单的问题；间接法建立在仿真分析应力应变的基础 上，得到的位移、应力场满足基本方程和边界条件，使稳定性分析具有正确的基础。 虽然有限元法有着各种显而易见的优势，但其自身始终存在一些不足：有限元法是 一种对于材料参数极其敏感的计算方法，不同的材料参数的选择对结果的影响非常之大。 而计算中所用到的岩土体及混凝土均是复杂的、各向异性的、非均质和非线性的弹塑性材 料，因此用有限元法分析计算时，必须考虑相应的材料非线性问题。现在的数值分析法 都不研究拱坝的实际破坏模式，所以没有一种拱坝整体安全度的数值分析方法。目前有限 元分析中常用的超载法和强度储备法，两种方法即使选用相同的失稳判据所得到的结果相 差都比较大，所以分析的时候，选取合适的研究方法对最后的结果也有着很大的影响。 在有限元分析整体稳定的时候，一直没有一种确定的判据可以用来实施对最后达到平衡状 态情况的判断。一般情况下，不同的判据得出的结果也有着不同程度的差别，因此在极限 分析问题的时候，我们需要慎重对待判据选取的问题。 1 3 本文的主要工作 针对以上各种问题，本文的主要工作如下： ( 1 ) 在查阅大量资料的基础上，对拱坝的失稳方式和整体稳定分析的研究方法进行 总结和综述；系统介绍非线性有限元基本理论、基本方程及常用分析方法。 5 河海大学硕士学位论文 ( 2 ) 对大型商用软件a b a q u s 进行详细介绍。对于几何非线性和材料非线性在该软 件中的实现进行各自具体分析，介绍a b a q u s 中非线性分析方法、初始地应力的施加问 题及该软件自带的结果后处理技术。 ( 3 ) 介绍综合法分析整体安全度的基本理论和基本思路，通过分析，提出几种综合 法的具体实施方案。以坝的整体安全作为研究主体，讨论安全判据的类型。通过具体算例 分析比较不同加载方法和安全判据对整体稳定安全度研究的具体影响，并确定最实用方 案。 ( 4 ) 用大型商用软件a b a q u s 对锦屏一级高拱坝进行三维非线性有限元模拟计算， 对不同荷载情况下坝体、坝肩和建基面的位移场、应力场和屈服区进行分析，了解坝基中 软弱结构面破坏情况，以及整个拱坝系统的破坏发展过程，确定锦屏一级高拱坝的整体安 全度。 6 第二章非线性有限元基本理论及软件应用 第二章非线性有限元法基本原理及软件应用 有限元法被用到岩土和水工问题分析中起源于1 9 6 6 年，美国的克拉夫( c l o u 曲) 和伍 德沃德( w r o o d 唧d ) 将有限元法引入到土坝分析中，自此，该方法在岩土和水利工程领域 发展迅速，目前使用有限单元法探讨复杂的水工课题己经相当普遍，其中包括应用有限元 法研究拱坝的稳定分析问题。 2 1 有限元法基本原理 有限元的概念可以简述为：采用有限自由度的离散单元组合体模型去描述实际具有无 限自由度的考察体，是一种在力学模型上进行近似的数值计算方法。 一般说来，线弹性问题的有限元分析的主要步骤如下【3 0 川】： ( 1 ) 离散化 离散化就是将连续体变换为离散结构，用有限多个在节点处相连起来的单元体结构代 替原来的实体，从而使研究对象转为可以使用计算机计算的数学模型。结构被离散化后， 可以确定一些进行计算的必要参数，如节点总数和单元总数等控制参数、单元和节点间的 对应关系，单元材料性质参数等。划分单元的大小、形状和数目应该根据具体问题来确定。 此时，还需要对位移边界进行约束处理，将连续的位移边界条件离散为节点的位移边界条 件。 ( 2 ) 单元分析 在每个单元内，将单元内各点位移用一个假定函数来表示，此函数用节点位移的值来 表示，这个函数即为位移函数。通过位移函数，就可以用节点位移值唯一的表示单元内任 一点的应变，然后再利用广义虎克定律，用节点位移唯一地表示单元内任一点的应力。利 用能量原理，可建立单元内部应力的等效节点力与节点位移的关系。对每一单元所承受的 荷载按静力等效的原则移置到节点上。 ( 3 ) 整体分析 对于每一个节点，建立与未知节点位移相应的平衡方程，集合后可以得到一个以节点 位移为未知量的线性方程组，即为整个结构的平衡方程组： 【k p = 仁) ( 2 1 ) 其中：【k 】为整体劲度矩阵，由单元劲度矩阵整合得到；p 为整体节点位移列阵，即为所 7 河海大学硕士学位论文 需求解部分；仁 为整体等效节点荷载列阵，由单元等效节点荷载整合得到。 2 2 非线性有限元法的基本原理 有限元的线性分析已广泛作为设计工具被采用，但绝大多数实际问题是非线性问题。 在分析非线性问题的时候，一般需要首先再清楚问题属于何种非线性，是单一非线性问题 还是复合非线性问题，然后是要针对不同的非线性问题选取不同的解决方法。而且在非线 性问题的解决中，非线性方程的迭代方法的选择对最后的结果也有着很大的影响，下面将 依次对这一系列问题做出详细说明。 2 2 1 非线性有限元的分类 根据非线性产生的原因，非线性问题主要分为三种：一是由于结构大变形所引起的非 线性，叫几何非线性；二是由于材料特性引起的非线性，叫材料非线性；三是由边界条件 的可变性和不可逆性产生的非线性问题，称为接触非线性。 ( 1 ) 几何非线性 几何非线性属于大变形问题，此时位移的大小会影响到结构的响应情况，发生的原因 可能是：大挠度或者大转动；突然翻转；初应力或载荷刚性化，此时反映应变和位 移关系的几何方程是非线性方程。例如，正应变g ，可表示为： s ，= 罢+ 丢i ( 罢) 2 + ( 罢) 2 + ( 尝) 2l + ( 2 2 ， 剪应变以，表示为： ，。，：堡+ 丝+ 塑丝+ 塑堡+ 堂堂+ ( 2 3 ) y = 一+ 一+ 一一+ 一一+ 一一+ t ，j 7 叫 缸。却。锄却叙却。融却 如果应力和应变之间的关系也是非线性的，就变成了更复杂的双重非线性问题。不过，在 几何非线性问题中一般都认为应力在弹性范围内，仃与占之间呈线性关系。工程中的实体 结构和板壳结构都存在几何非线性问题，例如弹性薄壳的大挠度分析，压杆或板壳在弹性 屈服后的稳定性问题。 ( 2 ) 材料非线性 材料非线性有时也称为物理非线性，就是指材料发生屈服后，材料的响应成了非线性 8 第二章非线性有限元基本理论及软件应用 和不可恢复的，其特点是应力盯与应变占之间为非线性关系，通常来讲就是材料加载和卸 载不是同一路径( 如图2 1 ) 。 o o 仉 o 图2 1 非线性材料加载和卸载路径图 因而其物理方程可写成： 仃= d 占 ( 2 4 ) 式中的弹性矩阵d 是应变s 的函数。但材料非线性问题仍属于小变形问题，位移和应变是 微量，其几何方程是线性的。 在采用有限元方法分析非线性问题时，以上两种都表现为结构的整体劲度矩阵k 不再 是常量矩阵，而是结点位移万的函数。 本论文不考虑接触非线性问题，因此在此不介绍接触非线性的基本原理。 2 2 2 非线性有限元的基本方程 在非线性力学中，有多种类型的非线性问题，但无论是哪一类非线性问题，经过有限 元离散后，它们都归结为求解一个非线性方程组： ( 4 如瓯) = 0 y 2 ( 4 疋瓯) = 0 ( 2 5 ) 眠( 西岛瓯) = 0 其中磊，疋，瓯是未知量，y 。，：，是4 ，疋，瓯的非线性函数，引用矢量记号： 万= 4 疋瓯】r ( 2 6 ) y = 【5 c ，ly 2 。】r ( 2 7 ) 9 河海大学硕士学位论文 方程组( 2 5 ) 可表示为 该方程组还可以改写为以下形式： 沙( 万) = 0 ( 2 8 ) 少( 万) 三f ( 万) 一r 三k ( 艿) 万一尺= 0 ( 2 9 ) k ( 万) 是一个”刀的矩阵，其元素是未知节点位移艿的函数；r 为等效节点荷载， 是已知矢量；f ( 艿) 是等效节点力；方程沙( 万) = 0 表示节点的平衡方程。 在线弹性有限元中，方程组( 2 9 ) 为线性代数方程组，可以毫无困难的求得最后结果。 但在非线性方程中，计算结果的过程就比较复杂了。求解这类方程，我们通常采用线性化 方法，即把非线性问题转化成一系列线性问题，为了使这一系列线性解收敛且逼近非线性 解，人们用过且提出过许多具体解法，但每种解法都有各自的局限性，因此选用正确的解 法分析问题也成为非线性有限元分析的一个重要问题。 2 2 3 非线性有限元方程的解法 一般情况下，非线性有限元方程分析的方法分为以下几种【3 2 ，3 3 ，3 4 ，3 5 ，7 3 】： 非线性有限元求解方法 迭代法( 总荷载法) 直接迭代法 k w f o ”一j 己印加d ”法 修正的肫w 幻玎一r 印加o ”法 拟肫w 幻力一脚凰d ，l 法 增量法c 逐步法) 善意兰船砂法 混合法c 逐步迭慝二菡羔舢 其中，迭代法是对总荷载进行线性化处理，从而求解未知量。而将荷载增量进行线性 化处理，则称增量法，它的思想是将荷载分成许多小的荷载部分( 增量) ，每次施加一个 荷载增量。此时，假定方程是线性的，劲度矩阵k 为常矩阵。对不同的级别的荷载增量， k 是变化的。这样，对每次增量求出位移增量万，对它累加，就可以得到总位移，实际 上就是一系列线性问题代替了非线性问题。对同一非线性方程组混合使用增量法和迭代法 就是综合法，一般在总体上采用e u l e r 增量法，而在同一级荷载增量内，采用迭代法。 各种方法都有各自的利弊，但是优缺点一般都是相对的，优点往往是以某种缺点为代 1 0 第二章非线性有限元基本理论及软件应用 价而获得的，各种解法的适应性一般都有一定的限度，对某种问题适合，对另一种问题则 不一定适合，尤其在某些情况下，得不到唯一解，或者即使得到某一解答，但不是所期望 的解答，这是与线性问题不同的一个重要特点。因此，迭代方法的选取应视具体情况而定。 2 3 本构方程 本构关系或本构方程一般是给出作用于材料质点上的力和变形及温度历史的关系方 程。建立本构方程时应遵循以下基本原理【3 4 j ： ( 1 ) 确定性原理，物体中的物质点在时刻t 的应力状态由物体中各物质点的运动历史 唯一确定。 ( 2 ) 局部作用原理，物体中的物质点的应力状态与离开该物质点有限距离的其他物 质点的运动无关。 ( 3 ) 坐标不变性原理，本构关系与坐标系选择无关。张量形式的本构关系就自动满 足坐标不变性原理。 ( 4 ) 空问不变性原理，本构关系对于做刚体运动的参考系具有不变性，即本构方程 的函数形式不因参考系做刚体运动而变化。空间不变性原理又称客观性原理。 ( 5 ) 时间评议的无关性，在时间t 的平移变换f = ，一f ( f 为任意常数) 之下，本构 关系的形式不变。 2 3 1 非线性弹性本构方程 材料非线性本构方程也有弹性本构方程和塑性本构方程之分。非线性弹性本构方程可 以看成是线性弹性本构方程的推广，也可以看成是弹塑性本构方程的特例( 形变理论的本 构方程) 。非线性弹性本构方程的表示方式可以分成全量形式和增量形式。 全量形式的非线性弹性本构方程的表达式可写成以下形式【3 5 】： p = 晒卜b p 仁艇 = k 仁) l 叠 ( 2 1 0 ) 式中：p 和【d 】与线性弹性问题相同，【d 】为线性弹性矩阵，而【d p 仁) 】是修正弹性矩阵， 具体为： 河海大学硕士学位论文 捌_ ( 志一訾 2 12 1 一l 00 o 0 o 0 p ) = 【g l l9 2 2 s 3 32 占1 22 s 2 32 占3 l 】r 阪仁m 是割线非线性弹性矩阵，表达式为： k 仁地= 志 l + 2 口 1 一口1 + 2 口对称 1 一口1 一口1 + 2 口 00o 0 o 0 o 0 0 口：三! ! 二型盟 3 e s 3 一口 2 3 u 一口 2 oo三口 2 ( 2 1 1 ) ( 2 1 2 ) ( 2 1 3 ) ( 2 1 4 ) 把( 2 1 0 ) 式中的各个函数的改为增量形式，可得到非线性弹性本构方程的增量形式， 即为【3 5 】： 其中： 洳 = 【( ；) 】r 协) p 盯) = 【d 盯ll d 仃2 2c ，盯3 3d 仃1 2d 仃2 3 ( ，万3 l 】r k 仁) l = 志 1 + 2 口 l 一口l + 2 口7对称 1 一口1 一口71 + 2 口 oo三口， 2 ooo三口， 2 0ooo 1 2 三口， 2 ( 2 1 5 ) ( 2 1 6 ) ( 2 1 7 ) 3 2 称 3 2 o 对 3 2 0 0 2 o 0 o 第二章非线性有限元基本理论及软件应用 ，一厄 d f = 3 口，：三! ! 二型巡旦 3 ed 协 = 【如l l 出2 2 蜗32 如1 22 如2 32 如3 l 】r 2 3 2 非线性弹塑性本构方程 ( 2 1 8 ) ( 2 1 9 ) ( 2 2 0 ) 在弹塑性问题中，应力全量和应变全量之间的关系通常与加载历史有关，并不像非线 性弹性那样互为单值函数，而与非线性弹性本构方程相同的基于形变理论的弹塑性本构方 程只适用于简单加载情况，没有普遍意义。 由文献【3 4 】可知：理想刚塑性材料的本构方程可写为： f ， c r y ，如 ，= o 卜他硼苦 q 2 。 式中，c r y 为屈服应力，厂是应力偏量不变量：和以的函数，或者是主应力差的函数，一 般使用时，总是取= f 。 将总应变增量如 ，分解成弹性应变增量如；和塑性应变增量如f 之和，则理想弹塑性材 料的本构方程可以分成弹性部分和塑性部分的综合，即： 弹性部分： 蟛：鲁一素熹蛾m 磊 ( 2 2 2 a ) ” 2 g2 gl + v ”9 f o y f = o y 式( 2 2 2 ) 就是理想弹塑性材料的本构关系。 ( 2 2 2 b ) 视 。旦魄 =如 分部性塑 河海大学硕士学位论文 2 4a b a q u s 软件介绍 2 4 1 简介 a b a q u s 是一套功能强大的基于有限元方法的工程模拟软件，它可以解决从相对简单 的线性分析到极富挑战的非线性模拟等各种问题。 a b a q u s 具备十分丰富的单元库，可以模拟任意 实际形状；a b a q u s 也具有相当丰富的材料模型 库，可以模拟大多数典型工程材料的性能，包括 金属、橡胶、复合材料、聚合物、可压缩弹性泡 沫、钢筋混凝土以及类似土壤和岩石等地质材料。 作为一种通用的模拟工具，应用a b a q u s 不仅能 解决结构分析问题，而且能模拟和研究包括热传 导、质量扩散、电子元器件的热控制、声学、土 壤力学等广阔领域中的问题。有限元软件主要由 前处理、计算和后处理三部分组成，a b a o u s 有 限元软件也不例外。这三个步骤之间通过文件建 立的联系如图2 2 所示： 2 4 2 几何非线性 在几何非线性问题中，单元形状的变化对应 力平衡和外荷载条件产生一定的影响。几何非线 l前处理 fa b a q u 8 c a e 或其他软件 i后处理 ia b a q u s c a e 或其他软件 图2 2a b a q u s 中的文件联系 性问题是不采用小位移假设，从几何上严格分析单元体的尺寸、形状变化，得到非线性的 几何运动方程，由此造成基本控制方程的非线性问题。下面通过一个典型的悬臂梁的例子 来分析几何非线性对计算结果的影响( 悬臂梁结构如图2 3 所示) ： 梁长8 0 m ，宽1 m ，高l m ，在自由端施加7 0 0 l n 的集中荷载，梁的杨氏模量为7 0 g p a ， 由材料力学3 q 相关理论可知，在集中力，作用下，梁的自由端静态挠度近似为： w ：堡( 2 2 3 ) 矿= 一 l z z jj 3 e l 其中：，= 砌3 1 2 ，6 和厅分别是梁的宽度和高度；，是长度。 1 4 第二章非线性有限元基本理论及软件应用 x。j c 一万 f 1 -r1r a 1p 。 一一、爹 x 7 ，j 。 矽1_-7 r y 图2 3悬臂梁受力图 当f = 7 0 0 剧时，自由端的挠度根据公式( 2 2 3 ) 计算出来的理论值为2 0 4 8 m ，相比 较原始长度8 0 m 来说，此时该问题已经是大挠度问题，因此，基于小变形理论提出的公式 ( 2 2 3 ) 是明显不适用于此种情况的，该公式计算出来的挠度结果只可作为结果分析时的 参考。若想得到真实位移，则需要考虑几何非线性时的悬臂梁的真实位移函数。 综合参考文献【3 7 】和文献【3 8 1 ，则水平向位移甜可写成【3 7 】： 删一厣i ( 2 2 4 ) 式( 2 2 4 ) 中，吼参考小挠度方程口= 警 竖向位移w 可写成【3 8 】： w = 一去警+ ；菇( 棚；一丢蠢( 铭帼；+ 素铭 眩2 5 ， 引2 苈肿，铭= c 器) 2 + 善一器 在使用a b a q u s 计算过程中，我们先不考虑几何非线性情况计算出自由端的挠度为 2 0 4 8 0 m ，与基于小变形理论的梁的挠度计算结果符合，但整个梁构件沿杆轴线方向位移 为零，此结果明显不符合实际情况；在考虑了几何非线性影响后，自由端的挠度减小为 1 9 2 5 0 m ，自由端沿杆轴线方向位移达到2 8 3 7 m 。通过公式计算出理论结果分别为1 9 2 7 3 m 和2 8 4 5 m 。虽然与计算值不能完全符合，但该误差友卜于可接受范围之内。 河海大学硕士学位论文 2 4 3 材料非线性 2 4 3 1 a b a q u s 中的本构模型 a b a q u s 提供了多种本构模型，极大地方便了混凝土、岩体、土体等水利工程材料的 模拟，常用的有m o h r - c o u l o m b 塑性模型、常规的d r u c k * p r a g e r 模型、扩展的d m c k e r - p r a g e r 模型、混凝土弥散开裂模型、混凝土损伤塑性模型等。 其中，m o h r - c o u l o m b 塑性模型用于模拟服从经典m o h 卜c o u l o m b 屈服准则的材料，在 岩土工程领域，可用来模拟单调荷载作用下材料的力学性状。扩展的d m c k e r - p r a g e r 模型 可用来模拟土、岩石等摩擦材料，以及压缩屈服强度远大于拉伸屈服强度的材料等相关材 料。混凝土损伤塑性模型能够用来模拟各种结构类型( 梁、桁架、壳和实体) 中的混凝土 和其他准脆性材料。 下面主要介绍了a b a q u s 中几种常用于模拟岩土、软弱夹层及混凝土的屈服准则。 ( 1 ) m o h 卜c o u l o m b 准则 c o u l o m b 于1 7 7 6 年首次提出了关于土的强度理论线性破坏准则，即： f = c t a i l 口 ( 2 2 6 ) m o l l r 于1 9 0 0 年将上式进行了推广，推广后所得表达式通称莫尔一库伦( m o l l r c o u l o m b ) 强度准则，该准则自提出之后即成为岩土、混凝土类材料常用的屈服准则。该 准则假定：作用在某一点的剪应力等于该点的抗剪强度时，该点发生破坏，剪切强度与作 用在该点的正应力成线性关系。 如图2 4 所示，m o h r - c o u l o m b 强度准则形式为： f = c 一仃t a i l 矽 ( 2 2 7 ) 式中，r 一剪切强度；仃一正应力；c 一材料的凝聚力；矽一材料的内摩擦角。 将f 和万用主应力带入，则准则可写成： 丢( 盯。一+ 丢( q + c r 3 ) s i n 矽= 咖s ( 2 2 8 ) 1 6 第二章非线性有限元基本理论及软件应用 j “。 qq吧巧矿 = 半i m o i 小c o u l o m b 准则假定材料的破坏与中间主应力无关，但典型的岩土材料的破坏通 常会受中间主应力的影响。通常情况下，这种影响较小，所以该准则在大部分情况下具有 足够的精度。 a b a q u s 采用的本构模型是经典m o h r - c o u l o m b 屈服准则的扩展【4 9 1 ，采用m o l m c o u l o m b 屈服准则，包括粘聚力的各向同性的硬化和软化，但该函数的流动势函数在子午 面上的形状为双曲线，在万平面上没有尖角，因此势函数完全光滑，确保了塑性流动方向 的唯一性。 ( 2 ) d r u c k e r p r a g e r 准则 m o h r - c o u l o m b 准则在万平面上的屈服轨迹为六角形，它在主应力空间的屈服面有一个 奇异的顶点。为消除角点，1 9 5 2 年d 1 1 l c k e r 和p r a g e r 对其提出修正，他们建议用一个正圆 锥面来代替上述的不规则六角锥面，如图2 5 所示。 o n 出 图2 5各屈服准则在万平面上的曲线 1 7 河海大学硕士学位论文 在万平面上的屈服轨迹为一圆。其屈服函数表示为： 八，2 ) = 咀+ 以一k = o ( 2 2 9 ) 其中，口与r 为材料常数。它们与内摩擦角够和凝聚力c 的关系为： 口：下型型l ，r ：毒堕塑l ( 2 3 0 ) 口= 芦l 一，r = 芦_ 二一 l 厶j u ， 3 ( 3 s i n 伊)3 ( 3 s i n 伊) 式中，“+ ”号对应于d m c k e r - p r a g e r 圆锥面与m o h r - c o u l o m b 体的内角点相接，“一 号则 对应于d m c k e r p r a g e r 圆锥面与m o h r - c o u l o m b 锥体外接。 a b a q u s 对经典的d m c k e r p r a g e r 模型进行了扩展，扩展的d m c k e r p r a g e r 模型的 屈服面在刀平面上不是圆形，屈服面在子午面上包括线性模型、双曲模型和指数模型，本 文计划采用线性模型。 该模型由三个应力不变量表示，同时提供了偏平面上的相关联的非线性流动、单独的 剪胀角和摩擦角。 j p 图2 6子午面上的屈服轨迹 子午面上的屈服轨迹如图2 6 所示，表达式可写成： f = ，一p t a n 一d = 0 ( 2 3 1 ) 式中，f = 丢g 【l + 去一( 1 一去) ( 争3 】为偏应力参数，不同的f 对应石平面上拉伸和压缩的不同 应力值，所以增加了拟合实验数据的灵活性；( 秒，) 一线性屈服轨迹在p f 应力平面上 的倾角，通常指材料的摩擦角；d 一材料的凝聚力，与硬化参数仃，有关。当硬化参数由单 轴压缩实验参数仃。定义时，d ：( 1 一昙t a n ) ；当硬化参数由单轴