有限元基础理论教程_lecture06

1、根据力的平衡，结点上的结点力与结点载荷或约束反力平衡。列出全部结点的平衡方程，并按照结点编号顺序排列，就得到了以全部结点位移分量为未知量的方程组，对应的系数矩阵K就是整体刚度矩阵， P为载荷列阵。 PK(2-29)2.7约束条件的处理由于整体刚度矩阵K是奇异的，在解方程组之前，要根据约束条件修改整体刚度矩阵。先考虑结点n有水平方向位移约束，与n结点水平方向对应的平衡方程为：nxnnnnnnnnPvKuKvKuK.2, 1212, 1212, 1211 , 120nu在K矩阵中，第2n-1行的对角线元素改为1，该行中全部非对角线元素改为0；在P中，第2n-1个元素改为0。为了保持K矩阵的对称性，

2、将第2n-1列的全部非对角元素也改为0。应该换成下面的方程：如果结点n在垂直方向有位移约束，则（2-29）中的第2n个方程修改为，0nv在K矩阵中，第2n行的对角线元素改为1，该行中全部非对角线元素改为0；在P中，第2n个元素改为0。为了保持K矩阵的对称性，将第2n列的全部非对角元素也改为0。yyxyxxPPPPPPvuvuvuvuvuvu644331665544222211000000*00*0000*00*0000*00*0000*00*0000*00*00*00*00*00*00*000000*000000*例2.9、结构的位移约束条件如图所示，结构平衡的方程组如下。修改整体刚度矩阵。根

3、据结点4的位移约束，修改该整体刚度矩阵的系数。yyxxPPPPvuvuvuvuvuvu633166554422221100000000*00*0000*00*0000*00*00*00*00000010000000000001000000*00*00*00*00*00*00*000000*000000*000000000010*0*0*000010000010*00*0*00*00*00*00*00*000000*00000000000133665544332211yxPPvuvuvuvuvuvu称对根据结点1和结点6的位移约束条件继续修改整体刚度矩阵可以得到以下的形式，如果结点n处存在一个已

4、知非零的水平方向位移，这时的约束条件为，*nnuu *12, 122, 1212, 1212, 1211 , 12.nnnnnnnnnnnuAKvKuAKvKuK在K矩阵中，第2n-1行的对角线元素乘上一个大数A，将方程修改为， 处理约束条件的意义在于，强迫单元内部的分片近似位移场满足整体结构的位移边界条件，得到整个结构的近似解。 在定义位移约束时，要消除结构的刚体位移。1）用相关单元的单元刚度矩阵计算结点力，再由力的平衡关系得到约束反力。与位移约束对应的约束反力如何计算？2）给矩阵相应的对角元素加上一个大数，将载荷列阵的对应元素置为零。相当于用一个刚度很大的弹簧代替位移约束。nnnnnnnn

5、nAuvKuKvKuK.2, 1212, 1212, 1211 , 12n结点在x方向位移为零，用刚度系数为A的弹簧代替，由于刚度很大，un是一个很小的值。方程修改为，0.)(.2, 1212, 1212, 1211 , 12nnnnnnnnvKuAKvKuK2.8整体刚度矩阵的特点与存储方法整体刚度矩阵具有以下几个显著的特点：对称性，稀疏性，非零系数带形分布。1）对称性由单元刚度矩阵的对称性和整体刚度矩阵的集成规则，可知整体刚度矩阵必为对称矩阵。利用对称性，只保存整体矩阵上三角部分的系数即可。2）稀疏性单元刚度矩阵的多数元素为零，非零元素的个数只占较小的部分。3）非零元素带形分布整体刚度矩阵

6、的非零元素分布在以对角线为中心的带形区域内，这种矩阵称为带形矩阵。在包括对角线元素的半个带形区域内，每行具有的元素个数叫做半带宽。最大半带宽用d表示，2) 1(最大植相邻结点编码的差值的d例2.9所示结构的最大半带宽为，82) 13(d二维等带宽存储设整体刚度矩阵K为一个n行、n列的矩阵，最大半带宽为d。利用带形矩阵的特点和对称性，只需要保存以d为固定带宽的上半带的元素，称为二维等带宽存储。整体刚度矩阵K每行中的上半带元素取出，保存在另一个矩阵K*的对应行中，得到一个n行、d列矩阵K*。把元素在K矩阵中的行、列编码记为r、s，在矩阵K*中的行、列编码记为r*、s*，对应关系如下： r*=rs*

7、=s-r+1K矩阵中的对角线元素保存在新矩阵中的第1列中，K矩阵中的r行元素仍然保存在新矩阵的r行中，K矩阵中的s列元素则按照新的列编码保存在新矩阵的不同列中。仍然会保存一些零元素，但是元素寻址很方便。10*0*0*000010000010*00*0*00*00*00*00*00*000000*000000000001称对整体刚度矩阵K10*0*0*000010000010*00*0*00*0*00*00*000*00000001等带宽矩阵K*二维等带宽存储所存储的元素数量取决于最大半带宽d的值，d的值则由单元结点的编码方式决定。相同的有限单元网格按照图2.13(a)的结点编码，最大的半带宽为

8、14；按照图2.13(b)的结点编码，最大的半带宽为18。图2.13(a)图2.13(b)一维变带宽存储 按行的一维变带宽存储，按照每一行的半带宽把半带宽内的元素保存到一维数组中。 使用辅助数组定义对角元素在一维数组中的位置。 需要存储的矩阵元素最少。 解方程时的地址计算比较复杂，会带来一些额外的计算量。2.9 线性方程组解法线性方程组的解法： 直接解法包括高斯消去法、等带宽高斯消去法、三角分解法，以及适用于大型方程组求解的分块算法和波前法等。 迭代解法雅可比迭代法、高斯-赛德尔迭代法、超松弛迭代法和共轭梯度法等。利用矩阵的对称、稀疏、带状分布等特点提高方程求解效率是关键。在方程组的阶数不是特

9、别高时，通常采用直接解法。当方程组的阶数过高时，为避免舍入误差和消元时有效数损失等对计算精度的影响，可以选择迭代方法。利用矩阵的对称、稀疏、带状分布等特点提高方程求解效率是关键。在方程组的阶数不是特别高时，通常采用直接解法。当方程组的阶数过高时，为避免舍入误差和消元时有效数损失等对计算精度的影响，可以选择迭代方法。高斯消去法的一般公式第m次消元时，以第m-1次消元后的第m行元素为主元行，对第i行元素(im)的消元公式为，)1()1()1()1()(mmjmmmmimmijmijKKKKK)1()1()1()1()(mmmmmmimmimiPKKPP),.,2, 1,(),.,3 , 2 , 1

10、(nmmjinm若原系数矩阵是对称矩阵，则在消元过程中的待消矩阵仍然保持对称。 )1()1(mijmjiKK)()1()1()1()1()(mijmmimmmmjmmjimjiKKKKKK如果雅可比迭代的原理bAx 0ADA)(11jijnijjiiiixabaxbxAxbxAxkk1A为非奇异矩阵，主元非零D为对角矩阵)(111kjijnjiiikikixabaxxADIA1bDb1通过迭代逐步逼近方程组的解， 波前法的特点刚度矩阵K和载荷列阵P不按照自然编号进入内存，而按照参加运算的顺序排列；以集成完毕的自由度作为主元对其它列的元素进行消元修正。不形成整体刚度矩阵。 共轭梯度法求解最优化问

11、题的搜索算法。nRxxf),(min迭代公式为，kkkkdxx111)(kkkkdxfd搜索方向，线性方程组，Axbxf)(ANSYS提供了多种求解器供选择，分为直接解法和迭代解法。直接解法包括： 波前法（Frontal Solver） 稀疏法（Sparse Direct Sovler）迭代解法包括：雅可比共轭梯度法（Jacobi Conjugate Gradient Solver， JCG），不完全共轭梯度法（Incomplete Cholesky Conjugate Gradient Solver， ICCG）预处理共轭梯度法（Preconditioned Conjugate Gradient Solver， PCG）代数多格法（Algebraic Multigrid Solver，AMG）区域分割法（Distributed Domain Solver， DDS）3 弹性力学轴对称问题的有限元法3.1 用虚功方程建立有限元方程3.2 三结点单元位移