岩土工程数值分析方法

1、第6章 岩土工程数值分析方法 有限元法 边界元法 有限差分法 离散单元法 1 有限元法 概述基本思路：将复杂的结构看成由有限个仅在结点处联结的整体，首先对每一个单元分析其特性，建立相关物理量之间的相互联系。然后，依据单元之间的联系再将各单元组装成整体，从而获得整体特性方程，应用方程相应的解法，即可完成整个问题的分析分析过程：结构离散化确定单元位移模式单元特性分析集成总体特性接方程求未知量1 有限元法方法运用的基本步骤：步骤1：剖分:将待解区域进行分割，离散成有限个元素的集合元素（单元）的形状原则上是任意的二维问题一般采用三角形单元或矩形单元，三维空间可采用四面体或多面体等每个单元的顶点称为节点

2、（或结点)，分析过程。步骤2：单元分析:进行分片插值，即将分割单元中任意点的未知函数用该分割单元中形状函数及离散网格点上的函数值展开，即建立一个线性插值函数 步骤3：求解近似变分方程有限个单元将连续体离散化，通过对有限个单元作分片插值求解各种力学、物理问题的一种数值方 法。有限元法把连续体离散成有限个单元：杆系结构的单元是每一个杆件；连续体的单元是各种形状（如三角形、四边形、六面体等）的单元体。每个单元的场函数 是只包含有限个待定节点参量的简单场函数，这些单元场函数的集合就能近似代表整个连续体的场函数。根据能量方程或加权残量方程可建立有限个待定参量的代数 方程组，求解此离散方程组就得到有限元法

3、的数值解。有限元法被用于求解线性和非线性问题，建立了各种有限元模型，如协调、不协调、混合、杂交、拟协调 元等。有限元法有效、通用性强应用广泛，许多大型或专用程序系统供工程设计使用。 有限单元法的理论基础虚位移原理：受给定外力的变形体处于平衡状态的充要条件是，对一切虚位移，外力所作总虚功恒等于内力总虚功ATidxdyWieWWaSTSATbeldsFldxdyFW最小势能原理定义1：外力从位移状态退回到无位移的初始状态时所作的功称为外力势能定义2：应变能和外力势能的和称为总势能最小势能原理：位移状态l为真实位移状态的充要条件是，对应l的势能一阶变分为零。也即对应位移l的势能取驻值(对线弹性问题为

4、最小值)(*aSTSATbpldsFldxdyFE)(21aSTSATbATpldsFldxdyFdxdyE0pE 有限元法的基本方程单元位移函数),( | ),( | ) 1 ,(1111vuyxi),( | ),( | )3 ,(3333vuyxk),( | ),( | )2 ,(2222vuyxjyaxaayxvyaxaayxu654321),(),(假设任意点的位移：单元位移函数：332211332211),(),(),(),(),(),(),(),(vyxNvyxNvyxNyxvuyxNuyxNuyxNyxu Nvu或：321321000000NNNNNNN Tvuvuvu33221

5、1插值函数(形函数)321321000000NNNNNNN)()(),(312312133221yxyxyxyxyxyxycxbayxNiiiijiiyxNiyxNijiiii,3 , 2 , 1, 0),(3 , 2 , 1, 1),(311),(iiyxN形函数特点：单元应变矩阵 BNLvuxyyxxyyx00 xNyNxNyNxNyNyNyNyNxNxNxNNLB332211321321000000单元应变矩阵(几何矩阵)： 321332211321321000000BBBxNyNxNyNxNyNyNyNyNxNxNxNB 3 , 2 , 1,002100ibccbxNyNyNxNBii

6、iiiiiii单元应力矩阵 SBDDxyyx 321321SSSBBBDBDS单元应力矩阵： iiBDS 单元刚度矩阵岩土体或结构体发生虚位移，单元结点的虚位移为 ，相应的虚应变为 ，则根据虚功原理有： * tdAtdsPNtdAFNTAATTATTnn* * eeTeTeATTeTkFtdABDBFn* tBDBtdABDBkTATen单元刚度矩阵：总体刚度矩阵由于虚位移 的任意性，等式两边与其相乘的矩阵相等，则：设结构体剖分成n个单元，根据虚功有：总体刚度矩阵 * eekF PUKkFnieeTnieT1*1* 荷载列阵结点位移列阵总体刚度矩阵PUK等参元分析平面任意四边形单元 Tmmvu

7、vuvu2211结点位移矩阵：01234mmNNNNNNN0000002121插值函数：)1)(1 (41),()1)(1 (41),(43NN)1)(1 (41),()1)(1 (41),(21NN mixNyNyNxNBiiiii, 3 , 2 , 100几何矩阵： mBBBB21 iiiiNNJyNxN1单元刚度矩阵： mmpllpiTiTTAeWWJtBDBddJtBDBtdABDBkn111111)( 模型范围与边界效应模型范围岩土工程涉及无限域或半无限域，但处理问题时只能对有限域进行离散化模型范围可取结构体轮廓尺寸的34倍边界效应以小变形理论为基础的有限元法中，力与变形的影响范围是

8、无限域，因此，设定有限域，并假定模型边界的位移为零或为受力边界就会带来误差，靠近边界越近误差越大，靠近边界越远误差越小 初始地应力场与释放荷载初始地应力场自然状态的岩体处于一定的初始地应力状态，在结构荷载作用下，岩体内的应力为荷载产生的应力与初始地应力之和释放荷载由于初始地应力的存在，开挖将导致部分岩体卸荷，通常采用沿开挖面作用着与地应力等价的“开挖释放荷载” 施工建造过程的模拟开挖释放荷载1i1ii1ii1ihvihPivP空单元施工过程模拟Hh第1步开挖第2步开挖第3步开挖Hv空单元空单元空单元 节理及不连续面的模拟平面问题节理单元-GoodmanGoodman单元是无厚度4结点单元结点传

9、递切向力与法向力节理应力-应变关系nsnsnsuuKK00 xynsns节理单元本构模型nFnutannksFsutannsckmaxsF00节理单元刚度矩阵(假定位移沿单元长度线性变化) nnnnssssnnnsssnnssnseKKKKKKKKKKKKKKKKKKKKtlk200020200202202002020200202026称）（对考虑嵌入的节理单元模拟(考虑转动) nnnnssssnnnsssnnssnseKKKKKKKKKKKKKKKKKKKKtlk200020200202202002020200202024称）（对 nsnsnsuuKKKM0000000nKLK341变厚度节

10、理单元平面六结点变厚度节理单元相当于四边形等参元yx02134566161iiiiiivNuNvu位移函数：41)1)(1 (41iNiiii6 , 5)1)(1 (212iNii形函数:具有一定厚度的单元可按四边形等参元处理当厚度很小时按等厚度或无厚度节理单元处理 1111ddJBDBdVBDBkTVTensnsnsuuKK00 多节理岩体的模拟等效连续体利用节理单元来模拟密集分布或随机分布的节理与裂隙是不适宜的，会给离散化与计算带来诸多困难与麻烦从总体上考虑节理裂隙对岩体的影响，将岩体视为等效的正交异性、各向异性或各向同性体，即等效连续体将岩体中的节理裂隙当成存在于岩体材料内的一种损伤，因

11、如损伤力学的原理建立损伤模型来考虑节理裂隙的影响层状岩体均质各向同性岩体受彝族节理(层理)切割形成层状节理岩体，由于节理弱面的影响使岩体具有横观各向同性的特征，可按一般横观各向同性的连续体来建立有限元模型岩层走向与纵轴平行且该纵轴为一应力主轴时xyyxxyyxmGnvvvvvnvM0001)1 (0)1 (122121222)1 (2/)21)(1/(1111222111vEGEEnnvvvEMyxz岩层走向与纵轴正交，计算平面平行于层面xyyxxyyxnvvvnvvnvvnvM2210001012212222222122)1 (2/)21)(1/(1111222111vEGEEnnvvvEM

12、yxz 岩体工程中的弹塑性问题非线性分析的基本方法分段线性增量法：将总荷载分成若干增量nPPU1P0P0UnU误差 iiiPUK1荷载增量与位移增量的关系： niiUUU10总位移： niiPPP10总荷载：误差修正方法nPPU1P0P0U误差 11iiiiPPUK一阶自校正法：iiiUUU11UPU0P0U误差11iAiiPPUK牛顿迭代法：2UAPAU1P 有限元法的实现模型建立(范围及参数)前处理(模型剖分)形状函数几何矩阵本构关系单元刚度矩阵总体刚度矩阵单元结点位移结点位移列阵边界条件结点位移列阵模型应力与应变场2 边界元法 概述边界元法是同有限元法并行发展的另一类数值方法,该方法在岩

13、石力学中的应用自20世纪70年代以后有了较大的发展边界元法通常只须在边界上进行离散化,因而具有数据处理工作量小、占内存小、速度快等优点，但在处理多介质问题、复杂的非线性问题时效率低它以定义在边界上的边界积分方程为控制方程，通过对边界分元插值离散，化为代数方程组求解。它与基于偏微分方程的区域解法相 比，由于降低了问题的维数，而显著降低了自由度数，边界的离散也比区域的离散方便得多，可用较简单的单元准确地模拟边界形状，最终得到阶数较低的线性代数 方程组。又由于它利用微分算子的解析的基本解作为边界积分方程的核函数 ，而具有解析与数值相结合的特点，通常具有较高的精度。特别是对于边界变量变化梯度较大的问题

14、 ，如应力集中问题或边界变量出现奇异性的裂纹问题，边界元法被公认为比有限元法更加精确高效。由于边界元法所利用的微分算子基本解能自动满足无限远处的条件，因而边界元 法特别便于处理无限域以及半无限域问题。边界元法的主要缺点是它的应用范围以存在相应微分算子的基本解为前提，对于非均匀介质等问题难以应用，故其适用范围远不如有限元法广泛，而且通常由它建立的求解代数方程组的系数阵是非对称满阵，对解题规模产生较大限制。对一般的非线性问题，由于在方程中会出现域内积 分项，从而部分抵消了边界元法只要离散边界的优点。 边界元法有两种直接法：直接建立关于边界未知量的积分方程，通过离散化求得边界未知量，并进而求域内任一

15、点的场函数值间接法：设定一个在域内满足支配方程但包含若干未知系数的解，在边界上强迫其满足边界条件，求得该系数，进而求得边界上及域内各点的场函数值域边界元 直接边界元法基本方程*,ubtubt,相同结构第一状态下体积力、边界力与位移场：相同结构第二状态下体积力、边界力与位移场：由功的互等定理：dbudstuudbudst*外域外域*klP*kluklPklukb分布体力第一种情况第一种情况单位力单位力il方向第一种情况外力：在无限域 上i点，沿l方向施加单位集中力 内力：在轮廓线 上，k方向的应力 位移：在 内及 上，任一点在k方向的位移第二种情况体力：在无限域 上沿k方向有分布体力表面荷载：在

16、轮廓线 上，沿k方向荷载位移：在 内及 上，任一点在k方向的位移il*klP*klukbkPku由功的互等定理：dbudsPududsuPkklkklkilkkl*函数，具有如下性质：被称为DeltaDiracil 内点在外点在外点在ixfixfidxfiiil,21, 0上点在区域，上点在边界，式中：iiCuCduiilikil121)(2),( 1),(2),( 1(*方向方向方向方向yxlyxkdsPudsuPuCdbudsPududsuPkklkklilikklkklkilkkl当不考虑体力时：将边界离散成n个线段单元并假设 与 沿边界均匀分布： dsuUdsPPIPUuPuICPds

17、uudsPuCTijTijnjjijnjjijiinjjkklnjjkklilijj)()(*111*1*；二阶单位矩阵；式中：kPku直接边界元法边界支配方程： PUuH 边界位移列阵u 边界应力列阵P nnnnnPPPPccH1111100边界应力影响系数矩阵：边界位移影响系数矩阵： nnnnUUUUU1111 间接边界元法基本方程(不连续应力法)外域),(yxpj内域dsyx上分布的荷载沿闭合边界曲线),(yxp dsPdsPPFFFyxyx作用在微段ds上的荷载为：曲线上所有荷载在j点产生的位移： dsPuuj*曲线上所有荷载在j点产生的应力： dsPSj*将边界离散成n个线段单元并假

18、设 在单元内均匀分布： niijiniLijijPUPdsudsPuui11* niijiniLijijPSPdsSdsPSi11*曲线上所有荷载在j点产生的位移：曲线上所有荷载在j点产生的应力：间接边界元法边界支配方程： 应力影响系数矩阵量边界元表面应力总体向量边界元表面位移总体向式中： 1 SuKSUu边界元法求解平面问题的步骤模型建立(范围及参数)将边界划分成单元将原岩应力反作用在单元上利用基本方程求解边界单元上的作用力与位移利用开尔文基本解与功的互等定理求解内部点的应力与位移3 有限差分法 概述有限差分方法是将所有研究区域内的基本控制微分物理方程与边界条件近似差分方程表示，而将求解微分

19、方程的问题变成在研究区域内特殊点上求解代数方程的问题 差分公式有限差分离散化基础：以增量之比代替连续导数xuxudxdux0limyxxxi1i1iuxu有限差分网格连续函数 的泰勒展开yxu,3033202200! 31! 21xxuxxuxxuuu30332022003! 31! 21xxuxxuxxuuu30332022001! 31! 21xxuxxuxxuuu以0点为x，y坐标原点：203102231022xuuuxuxuuxu略去高阶微量：同理：204202242022yuuuyuyuuyuxuuxu2310差分公式(一阶二阶)：20310222xuuuxuyuuyu2420204

20、20222yuuuyu 75860241uuuuyxyxu)22(2111931032uuuuxxu差分公式(三阶四阶)：11930140444641uuuuuxxu876543210220224222241uuuuuuuuuyxyxu)22(21121042032uuuuyyu121042040444461uuuuuyyu应力函数的差分解xyxyxyyx22222为应力函数，则：设 8675220200312022004220220412121xyxyxxyyxyyx利用差分公式：4 离散单元法 概述离散单元法是20世纪70年代初兴起的一种数值计算方法，适合节理岩体的应力分析离散元法也将模型

21、划分成刚性单元，单元之间可以相互叠合，也可以相互分离单元之间相互作用的力可以根据力和位移的关系求出，而个别单元的运动则完全根据该单元所受的不平衡力与不平衡力矩的大小按牛顿运动定律确定单元之间不需满足变形协调方程 离散单元法的基本方程物理方程-力和位移的关系nFsFnFsFsssnnnUKFUKF切向刚度系数法向刚度系数snKK运动方程-牛顿第二运动定律iiFFMFiFFijijFxeMmFu/ IM / 合力：合力矩：加速度：角加速度：tUtUtUtUtUtU )()()()(0101加速度、速度、位移关系： 离散单元法的计算机实现力和位移的计算循环力-位移关系运动方程力边界条件位移边界条件力位移已知力求位移MIIFumumFumumyyyxxx 考虑阻尼的运动方