（水力学及河流动力学专业论文）基于ansys的土石坝稳定渗流场的数值模拟.pdf

大连理工大学硕士学位论文 摘要 建国以来，我国修建了大量的水利水电工程，土石坝因具有造价低、结构简单、对 自然条件的适应性强、抗震性能好、工作可靠、寿命长、施工管理简便等优点而被广泛 采用。在土石坝中，因渗流引起的滑坡、渗透破坏尤为严重，所以合理、正确地对土石 坝进行渗流分析是必不可少的。 本文在总结前人研究成果的基础上，对利用有限单元方法计算土石坝稳定渗流场进 行了研究和分析，并阐述利用体积加权法计算稳定渗流场的原理和步骤。然后以大型通 用有限元软件a n s y s 为平台，应用其自带的二次开发工具a p d l 语言，编写了渗 流计算命令流程序( s e e p a g e ) ，并利用u i d l 语言将之界面化，形成渗流分析模块。通过 算例验证了该程序的可行性和正确性，接着对均质和非均质，各向同性和各向异性的土 石坝以及烟集粘土心墙土石坝进行分析计算，并和其它渗流软件相互比较，结果较为相 近，证明了a n s y s 计算渗流场的实用性。 在渗流分析的基础上，采用瑞典圆弧滑动方法，利用已有的f o r t r a n 程序，确定烟 集粘土心墙土石坝的最危险滑弧条分及安全系数，计算得出的结果满足规范的要求，表 明坝体是稳定的。在文章的最后部分，对a n s y s 在三维渗流场中的应用以及如何计算 渗透体积力进行探讨和分析，并编写相应的命令流程序，得到的结果较为准确。本文的 研究成果不仅可以应用于计算土石坝的稳定渗流场，也可以应用到其它坝比如混凝土坝 的稳定渗流场的分析，为解决渗流问题提供参考。 关键词：渗流；有限元；a n s y s ；渗透体积力 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 n u m e r i c a ls i m u l a t i o nw i t ha n s y so ns t a b l es e e p a g ef i e l do f e a r t hr o c k f i l ld a r n a b s t r a c t s i n c et h en a t i o nw a sf o u n d e d ，a1 a r g en u m b e ro fh y d r a u l i ca n de l e c t r i ce n g i n e e r i n g sh a v e b e e nb u i l ti no u rc o u n t r y b e c a u s eo fm a n yv i r t u e ss u c ha sl o wc o s t ，s i m p l es t r u c t u r e ，g o o d a d a p t a b i l i t y ，g o o dc a p a c i t yo fr e s i s t i n ge a r t h q u a k ea n dl o n gl i f e ，e a r t hr o c k - f i l ld a m sa l e a d o p t e da b r o a d l y f o re a r t hr o c k f i l id a r n s ，p r o b l e m sa b o u ts e e p a g ea r ev e r ys e r i o u s ，s oi ti s n e c e s s a r yt oa n a l y z es e e p a g ei ne a r t hr o c k f 1 1 ld a m sf e a s i b l ya n dc o r r e c d y i nt h i st h e s i s ，i ti sa n a l y s e dh o wt oc a l c u l a t e dt h es t a b l es e e p a g ef i e l do fe a r t hr o c k 一删b y f i n i t ee l e m e n tm e t h o d ，t h e nt h ep r i n c i p l ea n dt h es t e po f t h ev o l u m em e t h o di sa l s oe x p a t i a t e d o nt h eb a s i so f f o r m e rr e s e a r c ha c h i e v m e n t s t h e nt h ea r t i c l et a k e sl a r g e s c a l eg e n e r a lf i n i t e e l e m e n ts o f t w a r ea n s y sa st h ep l a t f o r m ，w h i c hh a si t so w r ld e v e l o p m e n tt o o l a p d l l a n g u a g e ，a n dc o m p i l e st h es e e p a g ep r o c e d u r e ( s e e p a g e ) t h a ti sr e p r e s e n t e da si n t e r f a c eb y u s i n gt h eu i d ll a n g u a g ea n d i sf o r m e dt h es e e p a g ea n a l y s i sm o d u l e t 2 l i sp r o c e d u r eh a sb e e n c o n f i r m e df e a s i b l ea n da c c u r a t eb ye x a m p l e s ，t h e na n a l y s e sa n dc a l c u l a t e st h eh o m o g e n e o u s a n du n h o m o g e n e o u se a r t hr o c k f i l ld a m , i s o t r o p i ca n da n i s o t r o p i ce a r t hr o c k f i l ld a m ，y a n j i e a r t hr o c k f i l ld a mw i t hc l u n c hc o r e t h er e s u l t sb ya n s y si sc o m p a r e dw i t ht h er e s u l t sb y o t h e rs e e p a g es o f t ，w h i c hi sq u i t es i m i l a ra n dp r o v ea n s y sp r a c t i c a b i l i t yf o rc a l c u l a t i n g s e e p a g ef i e l d o nt h eb a s i so fs e e p a g ea n a l y s i s ，t h em o s td a n g e r o u ss l i p p e r ya r cs l i c e sa n ds a f e t y c o e f f i c i e n t so fy a n j ie a r t hr o c k f i l ld a mw i t hc h i n c hc o r ei sd e f i n e db yr u n n i n gt h ef o r t r a n p r o c e d u r ew h i c hu s et h es w e d e nc i r c u l a rs l i d i n gm e t h o d n 坨r e s u l t sw h i c hs a t i s f i e t h e r e q u e s to ft h es t a n d a r di n d i c a t ee a r t hr o c k - f i l ld a ms t a b l e i nt h ee n do ft h i st h e s i s ，i ti s a n a l y s e dh o wt ou s ea n s y s i nt h et h r e e d i m e n s i o n a ls e e p a g ef i e l da n dh o wt oc a l c u l a t et h e s e e p a g ef o r c e t h e nc o m p i l e sc o r r e s p o n d i n gp r o c e d u r er e s u l t so fw h i c ha r eq u i t ec o r r e c t t h e p r o c e d u r ei sn o to n l yu s e di nt h es t a b l es e e d a g ef i e l do fe a r t hr o c k f i l ld a r nb u ta l s ou s e di n t h es t a b l es e e p a g ef i e l do fo t h e rd a m ss u c ha sc o n c r e t ed a mi no r d e rt op r o v i d et h er e f e r e n c e f o rs o l v e i n gt h es e e p a g ep r o b l e m k e yw o r d s ：s e e p a g e ；f i n i t ee l e m e n t ；a n s y s ；s e e p a g ef o r c e 独创性说明 作者郑重声明：本硕士学位论文是我个人在导师指导下进行的研究工 作及取得研究成果。尽我所知，除了文中特别加以标注和致谢的地方外， 论文中不包含其他人已经发表或撰写的研究成果，也不包含为获得大连理 工大学或者其他单位的学位或证书所使用过的材料。与我一同工作的同志 对本研究所做的贡献均已在论文中做了明确的说明并表示了谢意。 作者签名：透啦日期：趁：! 生：如 大连理工大学硕士研究生学位论文 大连理工大学学位论文版权使用授权书 本学位论文作者及指导教师完全了解“大连理工大学硕士、博士学位论文版权使用 规定”，同意大连理工大学保留并向国家有关部门或机构送交学位论文的复印件和电子 版，允许论文被查阅和借阅。本人授权大连理工大学可以将本学位论文的全部或部分内 容编入有关数据库进行检索，也可采用影印、缩印或扫描等复制手段保存和汇编学位论 文。 作者签名 导师签名 驾虚丝 到皿 ( ；函年_ l 月且日 大连理工大学硕士学位论文 1 绪论 1 1 论文研究背景 从2 0 世纪初开始，渗流对工程的影响以及工程对环境的影响已为工程界广泛重视， 许多工程技术人员及学者从工程实践和理论两方面进行了大量研究，并取得了许多有价 值的成果，既解决了工程中的实际问题，又丰富和发展了渗流理论，随着相关学科的不 断发展和完善，各类工程实践提出的渗流问题日益广泛和复杂，伴随着实验方法的逐渐 改进，以及计算机在渗流计算中的普遍应用，渗流逐步形成具有自己的理论、研究方法 和应用范围的独立学科。同时随着人类工程活动类型的增多，活动空间不断扩大，出现 的工程渗流问题和工程诱发的环境渗流问题越来越多，越来越复杂，所造成的各方面损 失也越来越严重。主要表现有【1 】： ( 1 ) 在水利工程中，闸坝渗漏损失引起下游浸没，工程效益下降：土坝中渗流作用 引起的冲蚀和滑坡等破坏：渗透破坏和渗透变形导致垮坝等。 ( 2 ) 在建筑工程中，渗流引起的基坑边坡失稳：基坑支护结构的破环：基坑降水：基 坑涌水造成的淹没：地基处理：渗透变形和破环以及引起的周围土体移动，引发的地面 建筑和地下设施破环。 ( 3 ) 在采矿工程中，矿坑突水造成的矿井淹没：矿山排水造成的塌陷等。 ( 4 ) 在农业工程中，灌渠渗漏造成的工程效益下降：灌渠渗漏和灌渠渗入引起的土 地湿化、盐碱化：农田排水工程设计和盐碱地改良。 ( 5 ) 渗流有关的环境问题，地震引起的砂土液化：水库蓄水诱发的地震；大量抽取 地下水或开采液体矿采引起的地面沉降。 从以上各类工程和环境中的大量渗流问题中，可以看出渗流是工程设计、施工以及 安全使用的重要因素，也是评价工程效益、经济效益和环境效益的重要内容。 尤其是对于土石坝，渗流分析更显得重要。目前在国内外，土石坝所占比重达6 2 以上，高于2 3 0 m 的1 6 座大坝中，土石坝就有9 座【2 】。在我国，已建成的八万多座大坝 中，土石坝也占据很大比例。在导致土坝垮坝的各要素中，渗流起着主导的作用，据国 内外土坝失事的1 4 0 例统计分析，其中4 5 是由于渗流问题引起的【2 】，可见，渗流问题 是土石坝安全的关键。因此，如何更方便、更准确地计算土石坝的渗流是大家普遍关心 的问题。 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 1 2 渗流分析的发展状况 1 2 1 国外发展状况 1 8 5 6 年，法国工程师达西( h e n r id a r c y ) 通过试验提出了线性渗透定律，为渗流 理论的发展奠定了基础。1 8 8 9 年，h e 茹可夫斯基首先推导了渗流的微分方程。此后， 许多数学家和地下水动力学科学工作者对渗流数学模型及其解析解法进行了广泛和深 入的研究，并取得了一系列研究成果。但解析解毕竟仅适用于均质渗透介质和简单边界 条件，在实用上受到很大限制 5 0 i 。 1 9 2 2 年，h h 巴普洛夫斯基正式提出了求解渗流场电拟法，为解决比较复杂的渗流 问题提供了一个有效的工具。过去所沿用的电网络法都是基于差分原理，而最近研究了 基于变分原理的电网络法，使该方法得到进一步改进。随着电子计算机的迅速发展，数 值方法，即有限差分法、有限单元法和边界元法在渗流分析中得到了愈来愈广泛的应用。 有限差分法是在1 9 1 0 年由里查森首先提出的。经过长时期的研究和广泛应用，目前该 方法已具有较完善的理论基础和实用经验。有限单元法的基本思想是在1 9 4 3 年由柯朗 提出的，1 9 6 0 年克劳夫最先采用“有限单元法”这个名称，以与有限差分法相区别。1 9 6 5 年，滓克维茨提出有限单元法适用于所有可按便分形式进行计算的场问题，为该方法在 渗流分析中的应用提供了理论基础。目前该方法在渗流分析中的应用已十分广泛而有 效。边界元法是建立在经典力学理论基础上的。贝蒂互换定理及弗雷德霍姆积分方程早 在上世纪末及本世纪初就已提出，建立在这些理论基础上的边界元法初见于6 0 年代后 期，当时被称为边界积分方程法( b o u n d a r yi n t e g r a le q u a t i o nm e t h o d ) 。直到1 9 7 8 年，边界元法( b o u n d a r ye l e m e n tm e t h o d ) 这个名称才被确立并得到公认。目前，国 内外对边界元法的研究和应用都给了足够重视，在渗流分析中也得到一定应用。这里还 应着重指出，今年许多学者针对所研究问题的不同特点，研究和提出了能集合上述各数 值方法优点的杂交元法( h y b r i de l e m e n tm e t h o d ) ，使实际工程渗流问题得到更加合 理的解决。 1 9 6 2 年，m i l l e r 提出非饱和介质的渗透系数是含水量或压力水头的函数，这就为 达西定律应用于非饱和区提供了理论基础。1 9 7 3 年，n e u m a n 首先提出了用有限元法求 解土坝饱和一非饱和渗流的数值方法：后来，日本的赤井浩一采用了n e u m a n 的数值模 型和有限元法进行了试验和数值计算。1 9 8 6 年，r a n k w e r n e r 首先将自适应理论引入 到渗流的分析中，并将线性误差估计方法推广到求解具有自由面的二维非线性渗流问题 中。1 9 8 7 年c h u n g & k i k u c h i 讨论了二维非均匀多孔介质渗流问题的网格自适应问题。 通过迭代计算，提出了优化网格和确定自由面位置的方法。1 9 9 1 年b u r k l e y b r u c h 大连理工大学硕士学位论文 利用局部误差和整体误差的概念，发展了一种利用三角形单元进行无压渗流自适应分析 系统，在此基础上，1 9 9 6 年c h e n 发展了用四边行网格进行自适应渗流分析的方法，大 大提高了渗流计算的工作效率 3 , 5 1 , 5 2 】。 1 2 2 国内发展情况 我国有关科学工作者和工程技术人员结合大规模的坝工建设，对渗流理论、数值方 法和试验技术等都进行了广泛而深入的研究，并取得可喜的研究成果。南京水利科学研 究所的毛昶熙教授在渗流分析和控制领域开展了系统研究：河海大学的朱岳明教授提出 了排水子结构技术，准确模拟了排水孔幕的渗流行为：南京水科院的吴良骥同志对饱和 一非饱和的渗流计算做了重点研究，大连理工大学的刘迎曦教授则在混凝土重力坝的渗 透反演分析方面开展了探索研究 3 , 4 4 近几年来，随着计算机技术的发展和应用以及有限元法的迅速发展，有限元法在求 解渗流场问题方面取得了很大进展。渗流问题的理论相对完善，国内外对算法的研究也 很深入，但在实际工程应用中，一种相对简单可行、而又可以全面分析多种介质复杂边 界条件的渗流问题的有限元方法还有待于进一步研究，本文在这方面做了一些尝试。 1 3 土石坝中渗流计算的任务、内容及方法 1 3 1 土石坝中渗流计算的任务 渗流分析的主要任务是确定给定渗流场的水头、流速分布和渗流量等基本物理量， 并据此通过有关设计计算，对工程安全性和经济效益作出评价，选择合理的渗流控制措 施。 其次是通过分析找出各工程中渗流破坏易产生的部位和可能发生的隐患并进行针 对性的分析处理。在水利工程中易发生渗流问题的部位主要有：坝体，坝基，库岸。 坝体的渗流破坏，特别是均质土坝，由于施工控制质量不严，土石坝体碾压不密实， 坝体内留有隐患裂缝及通道或透水层，筑坝土料有溶蚀性，过多渗水或下游排水堵塞， 缺乏滤层保护，坝体浸润线抬高，填土湿度过大，填筑太快，以及运用管理时库水位放 空太快等原因，都可能导致事故的发生。 坝基渗流破坏一般与地基处理( 岩基或土基处理) 的质量有关，不均匀沉陷、软弱 夹层、滑动、超压、渗漏、管涌、液化等，都是大坝基础内部发生事故的主要因素。如 果在坝体内已经有了可见性裂缝，则可能表明坝基内部已发生过某种位移；地基内有软 弱夹层，下游又无减压排水措施，造成过高的扬压力，则会导致滑动；透水地基的持续 渗流和过量的渗漏容易导致管涌和溶蚀。 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 水库漏水的主要途径是库底或库边隐蔽的地下通道、断层带、裂隙组、溶沟、溶洞。 如果防渗措施不善，也会给水库蓄水带来困难或影响水库的安全。 渗流分析是土石坝设计的重要组成部分。工程实践表明，土石坝是否安全可靠和经 济合理，在很大程度上取决于能否正确进行渗流分析和选择合理的渗流控制措施。 1 3 2 石坝中渗流计算的内容 针对上述各种渗流问题，土石坝渗流计算包括如下内容： ( 1 ) 计算通过坝基坝体的渗水量。 ( 2 ) 计算坝基中的渗透力和水头分布，各关键点的水力坡降；土石坝体内的浸润 线；渗流出口坡降等，以便检验其稳定性。 ( 3 ) 预测、评价渗透变性( 或破坏) 1 3 3 土石坝渗流分析方法概述 研究土石坝的渗流问题的方法很多，归纳起来可分为三类：即基于各种公式计算渗 流特性的理论方法；基于绘制流网图确定渗流特性的图解方法；利用模拟原理确定渗流 特性的数值模拟和物理模拟方法，有些时候这些方法联合使用，便于互相校核。 下面简单介绍在渗流分析领域中应用较为广泛的三种方法：解析法、数值法和电拟 法。所谓解析法，系指利用有关数学手段直接定解基本微分方程的方法。通过解析解可 得到关于水头函数在所研究区域内分布的显式表达式。它既满足基本方程，又满足给定 的边界条件。一般地说，解析解是比较精确的，但其实用性差，这是因为到目前为止所 见到的解析解都是针对各向同性均质渗透介质和简单边界调节条件而建立的。 数值方法目前主要有有限差分法、有限单元法和边界元法。有限差分法是从微分方 程出发，将研究区域经过离散处理后，近似地用差分、差商来代替微分、微商。这样一 来，基本微分方程和边界条件地求解可归结为解一个线性方程组，所得结果为数值解。 有限差分法的优点是其原理易懂，算式简单，有较成熟地理论基础。其缺点是往往局限 于规则的差分网格，对曲线边界和渗透介质的各向异性模拟比较困难。有限单元法是对 近似计算的归纳和总结。它既吸收了有限差分法不离散处理的内核，又继承了变分计算 中选择试探函数，并对区域进行积分的合理作法，充分考虑了单元对节点参数的贡献。 有限单元法在模拟曲线边界和各向异性渗透介质方面比有限差分法具有较大的灵活性。 边界元法是通过把求解域边界剖分为若干个单元，化边界积分方程为线性代数方程组来 求解其数值解的。可见，建立所研究问题的边界积分方程是边界元法的基础。边界元法 与有限单元法相比，它便于处理无限或半无限渗透介质、渗流奇异( 如排水井点) 和自 由面等问题，而且由于它只对研究域的边界进行剖分，其数据信息量显著减少。一般说 大连理工大学硕士学位论文 其计算精度也高于有限单元法。边界元法的缺点是它所建立的系数矩阵满阵，而且是不 对称的，即使节点较少也将占据可观的内存。同时，它对三维非均质渗透介质问题的应 用尚存在相当大的困难。 电拟法是基于电场和渗流场符合同一形式的控制方程而进行求解的。电拟模型对渗 流来说是个数学模型，而不是物理模型。电拟法目前主要采用两种模型，即导电液模型 和电网络模型。由于导电液模型为连续介质模型，故它便于模拟急变渗流区问题，但用 它无法模拟非均质各向异性渗透介质，也不尽适应复杂的地质和边界条件。为了模拟更 加复杂的渗流场，逐步发展和研究了电网络模型，即电网络法。该方法既可基于差分原 理建立，也可基于变分原理建立，其基本原理是基于网络电路问题的解和渗流场的数值 解符合同一形式的差分方程和变分方程。由于基于变分原理而建立的电网络法吸收了有 限单元法的优点，故使该方法在模拟曲线边界和各向异性渗透性方面得到一定改进，电 网络法尽管在渗流分析中沿用已久，但由于它具有容量、稳定性基本不受限制和在解题 过程中不产生累积误差等特点，目前仍是求解大型复杂渗流场的有效工具。 本论文从实用出发，只应用了上述分析方法中的有限单元法。 1 4 本论文研究的指导思想和主要工作 在应用有限单元法分析土石坝渗流方面，前人已做过了许多有意义的事情，解决了 不少实际问题。但由于渗流问题本身的复杂性，采用数值模拟方法如何真实准确的模拟 渗流问题，是一个需要不断完善的过程。本文根据渗流基本理论和有限单元法的基本原 理，以大型通用有限元软件a n s y s 为平台，编写了相应的程序( s e e p a g e 程序) ，对土 石坝的稳定渗流场进行了仿真模拟计算，并利用其他渗流分析软件对a n s y s 计算的结果 进行了验证；然后在渗流分析的基础上，采用瑞典圆弧滑动法对土坝边坡的稳定进行了 分析计算。 本文具体工作内容为： 1 详细阐述渗流的理论基础及有限元分析的基本原理，并对求解渗流自由面的方法 进行评述。 2 简述大型通用有限元分析软件a n s y s 的基本功能，并阐明利用a n s y s 的热分析功 能模拟土石坝渗流场的基本原理和方法。 3 利用a n s y s 的a p d l 和u i d l 语言将整个渗流分析过程程序化、界面化，形成土坝 渗流分析模块。 4 利用a n s y s 对均质、非均质以及更复杂土石坝的二维稳定渗流场进行具体分析， 验证a n s y s 的热分析功能求解渗流问题的正确性和实用性。 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 5 在渗流分析结果的基础上，采用瑞典圆弧滑动法对土石坝边坡进行静力稳定分 析。 6 对利用a n s y s 计算三维渗流问题进行初步的探讨。 1 5 本章小结 本章针对目前水工结构中存在的主要渗流问题说明了渗流问题研究的重要性；介绍 了渗流分析的主要任务和国内外研究现状，总结介绍了渗流分析的基本方法，最后在总 结已有成果的基础上提出了本文的主要研究内容。 大连理工大学硕士学位论文 2 渗流有限元分析的理论基础及基本方法 2 1 渗流的理论基础 2 1 1 渗流的基本概念 地下水在孔隙、裂隙和溶洞中运动称为渗流。渗透水流所占有的空间区域称为渗流 场。 根据渗流所处的渗流场的不同，可将渗流分为饱和渗流和非饱和渗流。 从地表面向下，地下水一般可分为饱和带和非饱和带 j 4 1 。在非饱和带中介质的孔 隙中既有液相的水，也有水汽和其他气体，水的压力小于大气压力。非饱和带的下部是 毛细带，在毛细带中介质的孔隙逐步被水饱和，但其中水的压力仍然小于大气压力，可 视为非饱和带。水压力等于大气压力的界面称为自由面，是饱和带和非饱和带的分界面。 饱和带中水的压力大于大气压力。 在饱和带中的渗流称为饱和渗流，在非饱和带中的渗流称为非饱和渗流。 根据渗流的基本表征量，如水头、水力梯度、渗透流速的大小和方向是否随时间变 化，可将渗流分为稳定渗流和非稳定渗流。 当渗流的基本表征量中的任一或全部随时问而变化，则称此渗流为非稳定渗流，这 主要是由于天然或人为因素的影响，引起地下水的补给、径流和排泄条件的变化而造成 的，如降水、蒸发、地表水位变化、水井抽水等。 由于地下水位总是在不断变化着的，所以在多数情况下遇到的都是非稳定流。但当 地下水位变化不大时，可以将非稳定渗流当作稳定渗流考虑。稳定渗流模型具有计算量 小、简单的特点。 2 1 2 渗流的连续性方程 连续性方程是质量守恒定律在渗流中的具体应用，是渗流理论的基础，它表明，流 体在渗透介质中的流动过程中，其质量既不能增加也不能减少。 假定渗透水为不可压缩的均质液体，并且仅考虑垂直方向上的压缩，根据质量守恒 定律可以得到可压缩介质渗流的连续性方程【l 】： 冬+ 挈+ 誓：偌o + 叫罢 ( 2 1 ) 僦 o y 以口 式中，口为多孔介质压缩系数：卢为渗透水的密度；卢为水的压缩系数；v 。，v ，v 。渗 流沿坐标轴方向的分速度。 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 2 1 3 渗流的基本微分方程 根据达西定律，在非均质各向异性可压缩介质中有 8 h。8 h8 h ”x2 以i ； 7 ，瑚y 瓦。”：瑚：i 代入( 2 1 ) 式，得： 昙c t 争+ 斋c 。争+ 鲁限争= 偌o + 卅詈 设从= 偌缸+ 妒) ，则上式有： 昙( t 警) + 嘉( 等) + 昙( 屯警) = 段署 在q 内 ( 2 2 ) ( 2 3 ) ( 2 3 ) 式即渗流的基本微分方程，其中h = h ( x ，y ，z ，r ) 为水头函数；j | ，k y ，k ：为以五y ，z 轴为主方向的渗透系数；。为单位贮水量或贮存率；q 为渗流区域。 在不考虑土体压缩性或单位贮存率( 即，= o ) 时，上式变为拉普拉斯( l a p l a c e ) 方程形式的特殊情形： 昙( t 掣) + 导( 七，掣) + 昙( t 掣) ：0 ( 2 4 ) o x a yc r y 0202 目前，渗流分析有理论求解、物理模拟、数值模拟三种。 ( 1 ) 理论计算浸润线 该方法通过对工程作适当概化，采用相应近似公式求解，它只能适应于均质的简单 的工程，而对于有复杂边界条件的多介质渗流场，在计算理论未取得突破性进展前，该 法使用受限。 ( 2 ) 物理模拟试验模型模拟 用于研究渗流场的试验模型很多，如砂槽模型、粘滞流模型、电网络模型和导电液 模型( 即电模拟试验) 。从应用的广泛性看，当属导电液模型。由于该模型保持了原型 介质具有的连续性，并能精确的模拟原型的复杂边界和内部结构，因而是求解大型渗流 场的有效工具，但不足之处在于制作模型的费时费力。 ( 3 ) 数值模拟有限元数值分析 目前用于稳定渗流场的数值模拟方法中，应用较广的当属有限元法。有限元法主要 是求解渗流场内的水头函数，确定渗流场内的自由面和渗流量等渗流参数，有限元法即 把微分方程和边界条件按变分原理转变为一个泛函求极值的问题，首先把连续体或研究 域离散划分成有限个单元体，然后形成代数方程组，在计算机上求解。 大连理工大学硕士学位论文 2 2 有限单元法概述 有限元法是数值计算方法中应用最广的一种，它是在电子计算机的广泛应用和数值 分析方法发展的基础上发展起来的。有限( 单) 元法( f m i t ee l e m e n tm e t h o d ) 这一名称 是1 9 6 0 年c l o u g h 在计算结构分析论文中首先提出的，该方法很快就普及应用到整个固 体力学领域。随后有限元法逐渐被引用到流体力学领域，其中最早的是1 9 6 5 年 z i e n k i e w i e z 和c h e u n g 那篇求解拟调和微分方程的论文 3 0 , 4 8 1 。到7 0 年代，有限元法在稳 定渗流和非稳定渗流领域内得到了广泛的应用。有限元法应用广泛，发展之快，引人注 目，而且该法随着计算机的普及和完善仍在日益发展着。现在对各类坝型二维、三维稳 定渗流和非稳定渗流、饱和与非饱和渗流、非达西流及各向异性岩体裂隙渗流等都用有 限元方法进行分析，并已经应用于许多大型水利水电工程中。 有限单元法以剖分离散和分块插值为指导思想。其基本方法是将连续的求解区域离 散化为一组有限个、且按一定方式相互连接在一起的单元组合体，利用每一个单元内假 设的近似函数来分片地表达整个求解域上待求的未知场函数。由于单元能按不同的连接 方式进行组合，且单元本身又可以有多种形状，因此可以模型化几何形状复杂的求解域。 单元内的近似函数通常由未知场函数或其导数在单元各个结点的数值和其插值函数来 表达。这样一来，未知场函数或其导数在各个结点上的数值就成为新的未知量，从而使 一个连续的无限自由度问题变成离散的有限自由度问题。一经求出这些未知量，就可以 通过插值函数计算出各个单元内场函数的近似值，进而得到整个求解域上场函数的近似 值。显然随着单元数目的增加，或者随着单元自由度的增加及插值函数精度的提高，解 的近似程度将不断改进，只要单元满足收敛要求，近似解最后将收敛于精确解。 有限元方法的实施虽然也类似于有限差分法，但其实施方法不同。有限差分法是直 接从微分方程入手，以离散格式逐步近似逼近方程中的导数。有限元法的实施则相反， 按照变分原理求泛函积分找其函数值，即把微分方程及其边界条件转变为一个泛函求极 值的问题。有限元法是一种分块近似里兹( r i t z ) 法的应用，即首先把连续体或研究区 域离散划分为有限个、且按一定方式相互连接在一起的单元的组合体。再以连续的分片 插值函数建立一个个的单元方程后，依靠各结点把单元与单元连接起来，集合为整体， 形成代数方程组进行求解。 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 2 3 渗流分析的有限元法计算原理及实施步骤 2 3 1 渗流的基本方程和边界条件 对于符合达西定律，非均质各向异性土体，当坐标轴方向和渗透主轴一致时，其三 维渗流问题可归结为以下定解问题： 鲁( 也警) + 茜( 七，詈) + 0 也i o h ) 2 以詈在q 内 ( 2 s ) 初始条件： 日| f = o - - d o ( x ，y ，z ，o ) 在q 内 ( 2 6 水头边界： 驯l = 五( x ，y ，毛f ) 在r 1 上 ( 2 7 ) 流量边界： k 掣：五0 ，y ，z ，f )在r 2 上 ( 2 8 ) o n l ， 式中： q 一渗流区域，即为f 1 和r 2 所围成的的研究区域； 日= g ( x ，y ，= ，f ) 一域内各点水头值，为待求水头函数； l ( x ，y ，z ，o 一域内各点上的初始水头值； z ( x ，y ，z ，) 一r l 上的已知水头值； l ( x ，y ，z ，f ) 一r 2 的已知流量； 胛一r 的法线方向； l 一已知水头值的边界曲面； r 1 一已知流量值的边界曲面。 对于所研究的渗流场，根据变分原理，上述定解问题与下列泛函取极小值等价： 大连理工大学硕士学位论文 明) = 酬三卜嚅o h ，z + q ( 等) 2 + t 学) 2 日望a tj + f 2 h a l f ( 2 ” 2 3 2 渗流场的离散 在确定渗流场的边界后，将渗流场q 离散化，剖分为珊个互不相交的单元体p 。不 失一般性，设单元体的基函数m 是由单元体相应的m 个结点的位置坐标构成，则单元 体e 内任一点的水头表达式为： h = m 只 ( 2 1 0 ) 2 3 3 有限元计算公式 渗流场剖分成若干单元之后，渗流场就分解为各个单元之和。于是泛函式( 2 9 ) 相 应地分解为有关单元泛函之和，即 旭，= 薹赡i 筹n c 争2 吐c 参2 + 以日筹卜渺+ 圭j - 1 缈腑 ( 2 1 1 ) 将式( 2 1 0 ) 代入式( 2 9 ) 中，以，。( 日) 表示单元e 上的泛函，即： 以= 昧卜n q 学n 姒警，2 + 以日警肛+ 黟解 = i ；+ i ；+ i ； ( 2 1 2 ) 下面依次对( 2 1 2 ) 中各项求导数及其最小值。 第一项i ? ： = 鹏1 lk ( o 融h ) 2 坞学) 2 + 姒警) 2b 出 ( 2 1 3 ) 上式对单元各结点水头h i ，马，求导数，得： 盖= 者鹏卜警n qc 等n t 譬，2 一 = 三圳t 矗) 2 + 尼，者c 万o h ) 2 + 也高，2 j 蚴 ( 2 1 4 ) 将式( 2 1 0 ) 代入上式，得： 茎童塑! 茎王垒坚翌塑圭互塑整塞塑盟丝堡型 一一 令 k q - - 城t 警警+ k ，o 砂n fo 印n 嘶；警警 蔫卜，。 对单元8 的m 个结点水头求导数， 筹= 以刖导( 艺k = l 肌巩) ( 艺k = l = 从嘴巩争m 蚴 = 善争；s ：f ：n k n , d x d y d z 令s j = 如i 辩n i 捌皿| 鼬 8 s 1 1s 1 2 s 2 ls 2 2 s m ls 2 1 2 一 则： m )( 2 1 5 ) ( 2 1 6 ) ( 2 1 7 ) ( 2 1 8 ) 譬玑 嗡 忙 倦删 魄 渺 净筹 弘吲堕砂 叶 旦缸 盟融盟缸 唾孵坼缸 塑观 如；b l = 研一妈研一吗；研一吼 蝴 堕出 日 赡肛 牡 弘 出 妙 出 堕钟 有 旧r母 1 1 、rj 皿也： w i i u n i i 儿 鹏一嘲弼一；辑一吼 大连理工大学硕士学位论文 第三项e ：鬈是面积分，是表示r 2 边界的流量边界条件。将自由面边界l 看作流 量补给边界条件，即以= s 百o h ，s 是给水度，有： 弘rf 2 h d l + f 肛3 詈订2 蜉m 嘲兰。n k a h 。kd f + 黟荟m 崛订旺 e 对单元e 任一结点水头z 求导，有： 面o i 。2 黟j 仉bs n 罾m v t o h k 得 = 黟m 订+ 舻m 订，以订) f 姜 8 l ： o h l 8 i ： 8 h m 一时停卜坩 【优j 其中，【p 】。： 弓 ；弓：f f s n g n j a t ；舻) 这样，对任意单元e ，有 n 凡n i d f 驱氏n m d f ( 2 2 0 ) ( 2 2 1 ) 曷= 时埘+ i s 8 。+ p 】8 謦) 。俐8 z z ， 对所有单元得泛函求得微分之后，并使其等于零，就得到泛函对节点水头进行微分 的方程组 盖= 莩箐= 。卫，n ( 2 2 3 ) 式中，以为节点总数。 1 3 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 对上式汇总方程写成矩阵形式为 k m + 陋楞 + 【p 楞 = 式中扩) 为已知常数项。 2 3 4 单元渗透矩阵的形成 由式( 2 1 6 ) 得单元渗透矩阵为： k 】。= 蜀1世1 2墨 垦1如2 五2 k mk m l k w ( 2 2 4 ) ( 2 2 5 ) 其中，m 为单元结点个数，单元渗透矩阵的各元素为： 岛= 畎t 警警+ k ，o 砂n | o 砂n j m ；警警卜 c 2 筋， 以上各微分和积分都是在实际单元上对整体坐标进行的，可以利用复合函数的微分 法则和积分变量替换，变换为在标准单元上对局部坐标进行微分和积分a 由复合函数的微分公式： 柳： o x o n i 印 a ， 出 = 计1 8 n 鹫 o n i o r l 酬 西 其中m 为雅可比矩阵，即 纠= 苏 西 o x a 行 础 a f 砂 髫 砂 a 行 砂 a f 出 皙 出 a 刀 出 鸳 o n u a f o n m 口玎 8 n m d c h而 _ y 2z z ： y m2 把上面微分和积分变换公式代入式( 2 2 7 ) ，有 世，= r 。f ，f ，( t i o n , 百o n j + k ，o 砂n o 砂n j 嘶：i o n i o n j v , l 蜊朽 ( 2 2 9 ) 砭； hiiihiii二_j盟鸳毗一卸巩一西叭一鸳叭一却叭一鸳 一 壁型三奎堂堡圭堂垡塑 0 。 由于积分公式( 2 2 9 ) 比较复杂，故采用高斯积分公式把( 2 2 9 ) 化为下列形式： 2 喜耋砉他4c 砍警警+ q 警警+ t 警警黼孝r t , , g k ，眨，。， 式中：a i ，一。，a k 为三个方向相应的加权系数；孝7 ，可”，r 为积分点局部坐标；以为每 个方向取的积分点数。 单元右端项列阵各元素的计算式为： 群= l l 即i 拯+ 嗽n , d r 1i11 】( 2 3 】) = s l 卜n | j d 弘啊+ l l n ? j f d 弘q 式中，w 为入渗或蒸发水量； 为边界单位面积上流入空间渗流区域的流量；j 是蜕 化为二维情形的j a c o b i 矩阵。 2 3 5 总体渗透矩阵的形成 以上只是解决了单元渗透矩阵的问题，根据( 2 2 4 ) 式还要对共有节点的所有单元 组成的集合的每一个单元逐个计算，然后累加起来，对r 个节点都要如此，可以得到，2 个 方程： 医肛 = 妒 ( 2 3 2 ) 其中，医】为合成的整体渗透矩阵， 日 为结点水头列阵。 局= 聪：口j = 且 h l ： h ， h n + ： h w f = 甲 ( 2 t 3 3 ) 式中，是共有结点f ，的单元的个数；是结点总个数；甩是未知水头结点个数； n 一胛为已知水头结点个数。 对于总体渗透矩阵k 】，由于实际上每个结点的相关结点并不多，所以它是一个高 度稀疏矩阵，而且由于单元渗透矩阵的对称性，医】也是一个对称矩阵。扩) 是n 个元 蒙富强：基于a n s y s 的土石坝稳定渗流的数值模拟 素组成的列向量。结合边界条件求解式( 2 3 2 ) 就可以直接得到以结点水头值表示的近似 的稳定渗流场。 2 4 渗流自由面求解的基本方法及评述 在岩土边坡、土坝、地下洞室及地下水运动等渗流分析中，均存在有渗流自由面问 题。迄今，用有限元法、边界元法及离散元法等数值计算方法求解无压渗流时，最困难 最复杂的问题之一便是渗流自由面的确定，同时也是历年来国内外渗流工作者极感兴 趣的研究内容之一。由于渗流自由面的位置是预先未知的，属于混合边界问题，必须同 时满足d i r i c h l e t 边界条件( 第一类边界条件) 和n e u m a n n 边界条件( 第二类边界条件) ， 这使得求解渗流场问题很复杂。 用有限单元法求解有自由面的无压渗流问题，目前常用的方法有变动网格法和固定 网格法两种。 2 4 1 变动网格法 变动网格法是传统的求解渗流自由面的有限元方法。由于渗流自由面作为渗流域的 自然边界面是待定的，是一个非线性问题，因而需要迭代求解。迭代的具体步骤是： ( 1 ) 根据渗流概念和经验假定一条渗流自由面，以确定有限单元法的计算区域： ( 2 ) 将渗流自由面作为第二类边界，求解方程组得出渗流自由面节点的水头值日； ( 3 ) 比较假定渗流自由面节点的计算水头值日和其位置高程z ，视其是否满足条件 h = z ，若不满足，则用计算水头值日去改变节点的z 坐标，形成新的假定渗流自由面， 一般可以采用以下公式改变节点的z 坐标： z ：z ：+ 口：= z ：+ f ：一z ：1 t 2 3 4 ) 式中，z ? ，z ，分别为自由面上某一节点f 的第n 次和第月+ 1 次计算的自由面的z 坐标值； 日? 为第疗次计算出的自由面上节点f 的水头值：m 为修正系数，一般取为o 5 m 1 ： ( 4 ) 转入第( 2 ) 步，直到满足精度要求为止。 变动网格法最大的优点是渗流自由面和逸出点可以随着求解渗流场的迭代过程逐 步稳定而自行形成，迭代过程是收敛的。虽然移动网格法取得了许多成功的经验，但也 表现出方法本身的缺陷： ( 1 ) 当初始渗流自由面与最终稳定渗流自由面相差较大时，会使计算单元发生畸变， 乃至于相邻单元发生交替、重叠，以至于在计算过程中常需对渗流域进行重新剖分计算； ( 2 ) 当自由面附近渗流介质不均一，尤其是有水平介质层时变动网格会改变介质分 区： 大连理工大学硕士学位论文 ( 3 ) 当渗流域内有结构物时，变动网格常会改变结构的边界条件，计算精度受到影 响； ( 4 ) 在变动网格过程中，每一次迭代运算网格都随自由面的变动而变动，总体渗透 矩阵要重新形成，故需大量计时； ( 5 ) 在研究渗流与应力耦合作用中，由于应力分