计算固体11_精品文档优质PPT.ppt

1、第十一章第十一章 代数方程组的解法代数方程组的解法 前面已经介绍了固体力学微分方程通过不同的方法可以化为代数方程组对于线性静力问题化为线性代数方程组，对于非线性静力问题则化为非线性代数方程组动力问题的特征值求解方法已有第七章中介绍过了，本章将介绍线性和非线性方程组的解法 11.1 线性代数方程组的解法线性代数方程组的解法 解线性代数组的方法有直接法和迭代法两种在有限元法发展的初期，由于当时计算机的内存不大，在求解大型问题时常常采用对存贮量要求较小的迭代法，但随着计算机硬件的发展，内存越来越大，而且，在直接解法中可采用内外存交换的方法，因此，现在的有限元程序中基本上均采用直接解法 但是，随着解题

2、的规模越来越大，直接解法的累积舍入误差已经引起人们的注意，因此，现在有的商用有限元分析软件中除了直接解法外，还增加了迭代解法，供用户选用因为，迭代解法计算时间虽然较长，但误差是可以控制的所以，下面既介绍直接法又介绍迭代法118,119 11.1.1 线性代数的一些基础知识线性代数的一些基础知识 在叙述代数方程的解法以前，先介绍下面要用的一些线性代数的基础知识（一）向量的范数 向量的范数是衡量向量大小的度量概念 它的定义为：#对任意向量 ，按一定的规律有一实数与之对应，记为 ，若 满足 （1）正定性或非负性,仅当时 ，；#（2）（2）齐次性,对任意的实数 ，有 ;#（3）（3）三角不等式,对任意

3、的 ，有 这时称 为向量 的范数，按照这个定义，向量的范数可以有很多种但常用的有以下三种（1）1-范数（2）2-范数（3）-范数（二）矩阵的范数 类似于向量范数的定义，可以定义 阶矩阵的范数 定义：#对于任意的 阶方阵 ，按一定的规则有一实数与之对应，记为 ，若 满足（1）,仅当 时，；#（2）对任意的实数 ，都有 ；#（3）对任意两个 阶方阵 ，都有 ,；#（4）相容性条件 则称 为矩阵 的范数 常用的矩阵范数有三种，是由三种常用的向量范数诱导出的矩阵范数，设 为 阶方阵（1）1-范数（2）2-范数 其中 为矩阵 的最大特征值（3）范数 从向量的范数式与矩阵范数式相比较可知，矩阵1-范数是矩

4、阵列向量的1-范数的最大值，故又称为矩阵的列范数，2-范数与 的特征值有关，故又称为谱范数，矩阵的 范数是矩阵行向量1-范数的最大值，故称行范数（三）方程的条件数与状态 对非奇异矩阵 ，称 为矩阵 的条件数，记为 由于矩阵的范数定义不同，因而其条件数也不同，但它们有等价的性质，即矩阵条件数的大小是衡量矩阵“好”“坏”的标志 考虑线性方程组 ，若系数矩阵有小扰动 ，即原方程的系数矩阵成为 ，这时方程的解也有扰动，即 ，下面讨论 如何受 的影响，由 和 相减得到 从而有 两边取范数 因此得到 如果 充分小，使得 ，则上式化为 上式表明，当系数矩阵有扰动时，解的扰动与 的条件数有关条件数越大，解的扰

5、动也越大 同样，常数项 的扰动也会引起方程解的扰动，由 及 可得到 故有 又由 ，可得到 这里，条件数表示相对误差的放大率 很自然，称 大的矩阵为“坏”矩阵更确切地说，若 为 阶矩阵，当 时，称 为病态矩阵 这时，在求解过程中，小的误差可能引起解的失真特别当 很小时，总是病态的这一点很容易理解，当 时，为奇异矩阵，这时，方程要么无解，要么有无穷多解，因此，当 的值很小时，必会带来解的不稳定性 11.1.2 直接法直接法 设 阶线性代数方程组为 其矩阵形式为 其中 称为线性代数方程组的系数矩阵，若系数矩阵的行列式不为零，则方程组有唯一解以下均假设方程组存在唯一解 直接解法有许多种，下面介绍常用的

6、几种 （一一）高斯消去法高斯消去法 高斯消去法是解线性代数方程组的直接解法中最基本的一种方法它的求解过程包括消元和回代两个过程 消元过程是依次按计算 回代过程是逐步解出 消元后，方程组化为同解的上三角形方程组 上述方程组的系数矩阵是对角线元素为1的上三角矩阵.高斯消去法简单易行，从计算过程可以看出，要求 （称为主元素）不为零，因此，这方法只适用于从1到 阶顺序主子式均不为零的矩阵 计算实践还表明，高斯消去法的稳定性差，当出现小的主元素时，会影响计算结果的精度，甚至得到错误的结果产生这种现象的原因在于舍入误差，为了在计算过程中抑制舍入误差的增长，应尽量避免小主元的出现基于这种想法就有人提出了列主

7、元消去法和全主元消去法 （二二）直接三角分解法直接三角分解法 直接三角分解法有Crout分解法和Doolittle分解法两种 Crout解法就是将方程组的系数矩阵分解成两个矩阵的乘积 其中 为下三角矩阵,为对角元素为1的上三角矩阵 根据矩阵乘法规则可以看出，对于任意的 有 方程组可写为 先由前一个方程组解出 ，再由后一个方程组求出 Doolittle分解也是将系数矩阵 分解成与前面相同的形式，但它的下三角矩阵 的对角元素为1，而上三角矩阵 的对角元素不一定为1根据矩阵乘法公式可得到三角矩阵 的元素计算公式，按前面的公式可求出y,x （三三）平方根法和改进的平方根法平方根法和改进的平方根法 前面

8、所述的解法适用于一般的线性代数组，即系数矩阵可以是不对称的用加权余量法、边界元法和无网格法离散后所得的线性代数方程一般不是对称矩阵，所以可用上述的解法而有限元法所得的线性代数方程组的系数矩阵是对称正定的，矩阵的这一特征使它的三角分解具有更简单的形式，这就是平方根法与改进的平方根法 平方根法也称Cholesky分解法，对于对称正定的矩阵可以分解为 其中 为下三角矩阵，为它的转置矩阵，由逐行相乘可得 即有于是方程组的解可写为 因此有 有了分解式后，可以求矩阵 的行列式值 在平方根法中，前面的公式需要有 次开方运算，为了避免开方运算，有人提出了改进的分解形式 其中 是对角元素为1的下三角矩阵，为对角

9、矩阵，即只有对角元素 ，其他元素均为零的矩阵矩阵 的计算公式为,对于 计算 上述的计算公式避免了开方运算方程组可写为 由这三个方程组可依次解 针对有限元法所形成的代数方程具有稀疏，对称正是和带宽小的特点，Irans120提出了波前法，这个方法是结合三角分解和计算机内外存交换的一种解法，特别适合于大型代数方程的求解，因此，已为大多数有限元软件所采用具体解法可参看该文 11.1.3 迭代法迭代法 迭代法119能保持矩阵的稀疏性，具有计算简单，程序容易编制的优点，并在许多情况下收敛较快 迭代法的基本思想是用构造一串收敛到解的序列，即建立一种从已知的近似解计算新的近似解的规则 迭代法的一般形式就是对于

10、线性代数方程组构造下列形式的同解方程组 其中 为阶 方阵，为 维实向量，任意取初始向量 代入迭代式 由此产生序列 ，当 充分大时，以 作为方程组的近似解上式的 称为迭代矩阵这种迭代方法称为单步定常线性迭代法下面介绍几种迭代法（一一）雅可比雅可比（Jacobi）迭代法迭代法 雅可比法中第 次的 值 代入下式（二二）高斯高斯-赛德尔赛德尔（Gauss-Seidel）迭代法迭代法 高斯-赛德尔迭代公式为 在雅可比法中要保存两个近似解向量 ，而在上式中，每算出新近似解的一个分量 ，在算下一个分量 时，用新的分量 代替老分量 进行计算，这样，在整个计算过程中，只需一个近似解向量通常认为新近似解比老近似解

11、更接近精确解，因此，可望这样迭代会收敛更快（三三）松弛法松弛法 为了加快迭代过程的收敛速度，在高斯赛德尔迭代法的基础上引入一个称为松弛因子的参数 ，采用下列迭代公式 一般 ,即为高斯-赛德尔法，时称为超松弛法（Successive Over-Relaction，简写为 SOR）则为低松弛法（Successive Under-Relaction）松弛因子的选取对收敛速度影响极大，但目前尚无可供实用的计算最佳松弛因子的方法实际计算时，通常根据系数矩阵的性质及实际计算经验，通过试算来确定松弛因子的值，在弹性力学中，松弛因子可取1.8左右收敛较快（四四）最速下降法和共轭梯度法最速下降法和共轭梯度法通常

12、把最速下降法和共轭梯度法归为极小化方法当方程组的系数矩阵 为对称正定时，方程组的解即为下列二次泛函的唯一极小点 求 的极小点问题的最简单方法是最速下降法，具体的做法是从其个初始点 出发，沿 在点 处的负梯度方向（称为搜索方向）求得 的极小点 ，即然而再从 出发，重复上述过程得到点 ，如此继续下去得到序列 ，使得 可以证明，从任一初始点 出发，用最速下降法所得到的序列 均收敛于问题的解其收敛速度取决于比值 ，其中 分别是矩阵 的最小和最大特征值但当 时，收敛十分缓慢，因此在实际中较少使用这个方法，而使用共轭梯度法 共轭梯度法（Conjugate Gradient,简写为CG）又称共轭斜量法，其基

13、本步骤是在点 处选取搜索方向 ,使其与前一次的搜索方向 关于共轭，即 然后从点 出发，沿方向 求得 的极小点 ，即由此得到方程组的近似序列 从前面的极小点公式不难求出上式的解为若取 由共轭的定义可得共轭梯度法的计算过程可归结为第1步，取初始向量 ，计算第 ,步，计算理论上可以证明，如果计算过程是精确的，则对任意的 值，用共轭梯度法求解时最多进行 步就能得到准确解 这表明共轭梯度法实质上是一种直接法但在实际使用时，由于有舍入误差，破坏了这种方法的有限步终止性，因此，实际上是将它作为迭代法使用 共轭梯度法的主要优点是存贮量小，计算简单，它能充分利用矩阵 的稀疏性，计算时只需存贮中 的非零元素并且此

14、方法具有超线性收敛性 大量的数值经验表明，当 的条件数很小，或 的特征值大部分集中在一点附近时，用此方法只需迭代很少几步就能得到高精度的解如果系数矩阵 是病态或特征值较均匀地分布在一个很长的区间内，则此法收敛得很慢 为提高收敛速度，近年来提出了许多有效的实用算法，如预优共轭梯度法，广义共轭梯度法，分块SSOR-CG法和不完全分解预优共轭梯度法等在一些大型有限元软件中有预优共轭梯度法供用户选用11.2 非线性代数方程的解法非线性代数方程的解法 不论是材料非线性问题还是几何非线性问题，经过离散之后，它们都归结为求解一个非线性代数方程组：#其中 是未知量，是 的非线性函数如果引用向量记号 上述方程组

15、可以用简单地用一个向量方程表示：#为了讨论方便，有时将这个向量方程改写成如下形式：#这 是未知向量 的一个向量函数，是一个 的矩阵，矩阵的元素 是向量 的函数，而 是一个已知的向量。#在全量问题的有限元的位移表述中，代表未知的节点位移向量，是内力的等效节点力向量，是载荷的等效节点力向量，而向量方程 表示关于节点的平衡方程其中每一个分量方程对应于一个自由度的平衡对于增量问题，我们将把 理解为位移增量向量，将 理解为增量载荷向量11.2.1 直接迭代法直接迭代法 前面提到非线性有限元方程最后可以化为非线性代数方程组假设我们有第 次叠代的结果 ，则可得到相应的 ，于是可以求得第 次的未知变量 ：#式

16、中右边第二项中括号内的力学意义为第 次迭代后不平衡力，当它等于零时，是精确解，当它不等于零时，必须叠代求解，以逼近精确值 为松弛因子，用以加速收敛当 时，上式化为：#这是最简单的叠代格式对于 的数值称为超松弛因子它的最佳大小随结构的刚度矩阵不同而不同，在实际使用时要依靠经验在第一次叠代时需取一个初值即 直接叠代法的收敛速度与初值的选取也有很大的关系这种直接叠代法在多自由度情况下可能出现不稳定现象，所以这种方法很少在非线性有限元法中采用 11.2.2 牛顿牛顿-拉甫逊拉甫逊（Newton-Raphson）法法 牛顿-拉甫逊法是一个最著名的求解非线性代数方程的方法，简称为牛顿法 在方程中，我们已知有一近似解，可以用截断的一阶泰勒（Taylor）级数得到一个改进的解：#由方程式可知：#称为雅可比（Jacobi）矩阵当外载荷 不随位移而改变时（称为保守系统），上式右边第二项为零 为载荷-位移（或载荷增量-位移增量）空间的切线刚度。#由前一式可以得到：#上式中右端括号中的项为第 次叠代的不平衡力上式的形式与前面的叠代式相似，上式中采用的是切线刚度该法的叠代过程如图11.1所示因为它每一次叠代步后都修改切线刚度，以适应问题的物理性