1、1第七章 材料非线性问题的有限元法7.1 引 言前面各章所讲述的问题，都属于线性变形体系。所谓线性变形体系是指位移与载荷呈线性关系的体系，而且当载荷全部撤除后，体系将完全恢复原始状态。这种体系也称为线性弹性体系，它需满足下列条件：（1）材料的应力与应变关系满足虎克定律：（2）位移是微小的；（3）所有约束均为理想约束。在分析线性弹性体系时，可以按照体系变形前的几何位置和形状建立平衡方程，并且可以应用叠加原理。根据这种理论建立起来的方程是线性的，对于小应变和小位移的情形这种分析是适用的。实际结构的位移与载荷可以不呈线性关系，这样的体系称为非线性变形体系。如果体系的非线性是由于材料应力与应变关系的非
2、线性引起的，则称为材料非线性，如材料的弹塑性性质、松驰、徐变等。如果结构的变位使体系的受力发生了显著的变化，以至不能采用线性体系的分析方法时就称为几何非线性，如结构的大变形、大挠度的问题等。还有一类非线性问题是边界条件非线性，或状态非线性，如各种接触问题等。但本书只讨论前两类非线性问题的有限元解法，即材料非线性和几何非线性问题的有限元解法，对接触问题的有限元解法，读者可参考其它书籍。材料非线性问题的处理相对比较简单，通常不必修改整个问题的表达式，而只需将应力应变关系线性化，求解一系列的线性问题，并通过某种校正方法，最终将材料特性调整到满足给定的木构关系，从而获得了问题的解。对于几何非线性问题，
3、那就需要对公式进行根本的修改，这个问题将在后面详细讨论，不过应该指出，用于求解材料非线性问题的基本迭代方法也同样适用于几何非线性问题的求解。事实上，有些工程结构问题同时具有这两类非线性性质，它们可以统一地加以处理。本章将首先介绍用有限元方法处理非线问题的一般方法，然后讨论这些方法在非线性弹性、弹塑性和蠕变问题中的应用。在介绍弹塑性问题的处理方法前，为便于讨论，需扼要叙述一下 Mises 屈服准则和 Prandtl-Reuss 塑性流动理论，并据此写出弹塑性矩阵表达式。最后对平面刚架的极限分析做了简要介绍。7.2 非线性问题的一般处理方法非线性问题用有限元法离散化应得到如下形式的一组代数方法：
4、FKfP或写成（7.1）0f其中 。虽然线性方程组 直接求解并无困难，但对于方Ff程组（7.1） ，单元刚变矩阵是单元节点位移向量的函数，直接求解就行不通。然而，下面介绍的非线性方程组的各种解法，仍以反复地求解线性方程组去获得满足一定精度要求的非线性方程组的解答。7.2.1 直接迭代法2对于方程（7.1）（7.2）0fKfP最简单的求解方法是直接迭代法。开始求解时先假定一组初始值 代入上式的中，可求得改进了的一次近似值KfK101式中00重复上述过程，将迭代格式写成 fKnn1迭代一直进行到误差的某种范数小于预先规定的容许值 er，即满足ne1则停止。可以看出，该法的每一次迭代都需形成一次系数
5、 矩阵，并求解一次线性代数方1nK程组。这里还隐含着一个假定，系数矩阵 可以表示成 的显函数，因此该法只适用于与变形历史无关的非线性问题，例如非线性弹性问题及可利用形变理论分析的弹塑性问题。而对于依赖于变形历史的非线性问题，直接迭代法是不适用的，例如加载路径不断变化或涉及卸载及反复加载等必须利用增量理论分析的弹塑性问题。图 7.1 显示了单变量问题中这种迭代过程收敛和发散的可能性。通常，如 曲线P是凹的，则迭代发散。图 7.1 直接迭代法7.2.2 Newton-Raphson 方法（简称 NR 方法）若已获得方程（7.1）的第 n 次近似解 ，为了求得改进的近似解，可利用仅n保留线性项的 T
6、aylor 级数展开式（7.3）01 nnnd3则有（7.4）nnn1前式中 为切线矩阵，即d（7.5）TKdP于是从式（7.3）可以得到式（7.4）中的 n（7.6）fPKnTTn 11重复上述迭代过程，直至达到所要求的精度。Newton-Raphson 方法的迭代过程示于图 7.2 中，通常迭代过程是收敛的。但当所选取的初始值偏离真实解较大时，正如图 7.2（b）所表示的那样，发散也是可能的。图 7.2 Newton-Raphson 迭代法由式（7.3）看出，该法在每项迭代中必须重新形成一个系数阵并求解一次线性代数方程组。应该指出，如果原始的离散化方程组是通过变分原理导出的，则切线刚度矩阵
7、总是对称的，而利用直接迭代法，系数矩阵的这种对称性不一定能保持。TK7.2.3 修正的 Newton-Raphson 方法（简称修正的 NR 法）为克服 Newton-Raphson 方法中每次迭代都需形成一次系数矩阵和求解一次方程组的缺点，通常切线矩阵总是采用它的初始值，即令（7.7）0TnK因此式（7.6）现在变成（7.8）fPnTn10这样，每次迭代求解的是系数矩阵相同的方程组。对于这种系数矩阵不变的线性方程组，如果我们在迭代的开始就将 求逆或分解，则以后每次迭代只须进行一次回代过程。因TK而每次迭代的计算工作量大大地减小了，但收敛速度却变得较慢。不过从总的效果看，还是合算的。图 7.3
8、 表示了这种方法的迭代过程。一种改进方案是在迭代若干次后就将切线矩阵修正一次，修正到当前的值，这一改进方案有时是非常有效的。4图 7.3 修正的 Newton-Raphson 迭代法7.2.4 增量法为了便于理解，假定方程（7.1）式表达的是结构应力分析问题，其中 代表结构的位移， 代表结构的载荷。所谓增量法首先将载荷分为若干步：f，相应的位移也分为若干步： 。每二步之间的增,210f ,210长量称为增量。增量解法的一般做法是假设第 m 步载荷 和相应的位移 为已知，mfm而后让载荷增加为 ，再求解 。如果每步)(1mfff )(1载荷增量 足够小，则解的收敛性是可以保证的。同时，可以得到加
9、裁过程的中间结m果。为了说明这一方法，可将式（7.1）改写成（7.9）0fP其中 是描述载荷变化的参数。将上式对入求导，得到00 fdKfdT由此得出（7.10）01)(fT式中 仍为切线矩阵。TK有多种方法可用来求典型常微分方程（7.10）的解。其中最简单的是 Euler 法，它可表达成（7.11）mTmmTm fKfK 1011)(其中下标指示增量加载的次数，即 1或 mmfff1在实际的数值计算中，由于并不用式（7.9） ，而是用其增量形式（7.10） ，所以解常会发生漂移现象。为了克服这一缺点，可将 N-R 法或修正的 N-R 法用于每一增量步。如采用N-R 法，则对于 的 m+1 次
10、增量步，和 N-R 法的第 n+1 次迭代表达式可写成50)()( 10110111 nmTmnmnnm KfPfP （7.13）式中 是 时第 n 次改进的切线矩阵。由上式解出 ，则nTK11 n1（7.14）nmnnm111开始迭代时，取 ，然后连续地进行迭代，直至方程（7.9）在 时能够m01 m在规定误差范围内被满足。由式（7.13）可以看出，当采用 N-R 法进行迭代时，每次迭代都需重新形成 ，nTK1并求解一次代数方程组，致使计算工作量很大。因此通常采用修正的 N-R 法，这时式（7.13）中的 应代之为nmTK1（7.15）)(011mTmTnTKK如果每一增量步采用 N-R 法
11、进行一次迭代，由式（ 7.13）则有（7.16）)(01011 fPTm且假定在前一增量步结束时支配方程（7.9）是精确满足的，即 0fm则有（7.17）mTfK011这实际上就就是 Euler 法，即式（7.11） 。不过采用式（7.16）时，能将前一步支配方程的误差 在本次增量步中加以校正，这是带自校正的增量法。采用此法解方程0fPm时，解的漂移不严重。本节所讨论的上述算法是目前用于求解离散的非线性方程组的常用算法。由于用有限元分析非线性问题计算工作量很大，且有时收敛很慢甚至会导致解的发散，因而引起许多计算工作者的关注。一些加速收敛的措施和方法，好的修正计算方案已相继提出。读者如有需要可查
12、阅有关文献。在以上介绍的各种解法中，很难说哪一种算法最好。因为在某种情况下最经济有效的方法，在另一种情况下则不然，甚至解收敛很慢或不收敛。不过，要为一个通用程序编入一种解法时，增量法是合宜的。因为只要选择足够小的增量步，解总是收敛的，如可能，还可在每一增量步中采用 N-R 法或修正的 N-R 法，使计算结果满足一定的要求。7.3 非线性弹性力学问题如果我们只考虑小变形，则平衡方程在整个求解域内可写成（7.18）vTFdB上式中的积分运算实际上应逐个单元进行，并按单元集成法把它们对节点的平衡的贡献进行叠加。右端节点力向量也由单元集成法形成。6几何关系可写成（7.19）B但此时物理关系是非线性的，
13、一般可写成（7.20）0),(f由式（7.18）至（7.20）可以导出与非线性方程组（7.1）相同的表达式 )(FP因而上节中所讨论的各种解法，原则上都可在此应用。根据非线性弹性力学问题中非线性应力一应变关系的不同表达方式，可以产生下列几种求解方案。7.3.1 直接迭代法如果材料的应力一应变关系能表达成 )(D于是由式（7.19）可将应力写成 B)(将上式代入式（7.18）后，得（7.21）FK)(其中 vTdvBD)(用迭代法求解方程（7.21）时，总是首先取 ，然后求出 ，由000)(K求出位移的第一次近似值。重复这一过程，迭代格式为FK10（7.22）FKn17.3.2 切线刚度法（即
14、N-R 法）如果材料的应力一应变关系能表示成增量关系（7.23）dDdT)(式中 为切线弹性矩阵，可将式（7.18）改写为)(TD（7.24）vTFB0)(将 对 求导，注意到关系式（7.23） ，则有（7.25）TKd7式中（7.26）vTTdvBDK)(是切线刚度矩阵。在 N-R 法迭代中，可首先取 ，由式（7.19）求得 ，再由式（7.23）确00定 。将 代入式（7.26）求出 ，并由式（7.24）求出 ，由)(0TD)(00TK0。解出 ，于是得到位移的第一次近似值为 。1K01重复上述过程，可将迭代公式表示为（7.27）nnnT11式中（7.28）vnTnFdvB我们知道，式（7.
15、28）右端的第一项是与 相等价的节点载荷，因而 可理解为n结构上的不平衡节点载荷向量，每迭代一次相当于对不平衡节点载荷作一次矫正。7.3.3 初应力法初应力法实质上是修正的 N-R 法的具体应用。如果材料的应力一应变关系可以表示成（7.29）)(f即应力分量能由给定的应变分量确定，我们就可以用具有初应力的线弹性应力一应变关系（7.30）0D去代替上式。其中D 采用非线性应力一应变关系式在 时的切线弹性矩阵。这就是说，可用调整初应力值 的方法，使在给定应变 下，用式（7.30）算得的应力 与0用式（7.21）算出的结果相同。于是有（7.31）elDf )(0其中 。在一维情况下初应力的定义见图
16、7.4Del图 7.4 初应力法8如令（7.32）vTTdvBDK0为结构的起始切线刚度矩阵，则将式（7.30）代入（7.18） ，得（7.33）vTT PFF00式中（7.34）vTdvBP0为与初应力等价的节点载荷。根据式（7.33） ，可以采用如下迭代步骤：首先由下式计算第一次近似位移 FKT101由 求出 ，用式（7.31）求出初应力111110)(Df由式（7.34）求出 dvBPTv101用 对位移进行一次矫正1P101PK于是位移的第二次近似值是 112继续进行这种矫正，直至收敛。由此可见，用初应力法迭代时，每一步都是根据真实应力 与弹性应力 之nnD差决定初应力 的，并从而决定
17、不平衡的节点载荷 ，依次进行一次调整，使位移n0P进一步逼近真实值。初应力法迭代过程在单向应力状态下的应力变化示于图 7.5 中。图 7.57.3.4 初应变法9对于某些问题，例如蠕变问题，其应变由应力值决定(7.35)(f这时，上式可用具有初应变 的线弹性应力一应变关系0（7.36）)(0D来代替，即调整初应变值 ，使在给定应力 下，用式（7.36）求得的应变 与用式0（7.35）求得的结果一样。于是有（7.37）elff )()(110式中 。在单向应力的情况下初应变的定义见图 7.6。1Del图 7.6将式（7.36）代入（7.18） ，则有（7.38）vTT RFdvDBFK00其中（
18、7.39）vTR0为与初应变等价的节点载荷。 见式（7.32） 。0T根据式（7.37） ，可以采用与初应力法类似的迭代过程。先计算位移的第一次近似值 FKT101由 计算应力 ，由式（7.37）确定初应变 ，再由式（7.39）计算R 1，然后用11101R求出 ，对位移进行一次矫正1112重复上述过程，直到收敛。可见，初应变法的每一步迭代是根据真实应变 与弹性应变 之差，去对位n1D10移进行调整的。在单向应力状态下初应变迭代过程示于图 7.7 中。图 7.77.4 弹塑性应力一应变关系在弹塑性小变形情况下，弹性力学中的平衡方程和几何关系仍然成立，但描写应力一应变关系的物理方程却是非线性的。
19、7.4.1 材料的塑性性质对于大多数金属材料，单向拉伸试验的应力应变曲线如图 7.8 所示。当应力小于屈服极限 时，应力和应变呈线弹性关系（一般假定比例极限与屈服极限重合） ；而应力T达到屈服极限后，应力一应变显示出非线性关系。图 7.8（a）所示为应变硬化材料，图7.8(b)所示为理想弹塑性材料。由实验知，对于硬化材料，应力超过屈服极限 时，要增T加应变，应力也需增加。而对于理想弹塑性材料，应力达到 时，材料则发生流动变形。T然而无论哪种情况，当应力达到 后，卸载都是弹性的，卸载路经沿直线，其斜率大T致与初始加载路径相同。完全卸去载荷后，试件会留下残余变形。因此，弹塑性材料当应力达到或超出
20、时，应力一应变之间并无一一对应关系，即应变不仅决定于当时T的应力，而且还与整个加载的历史有关。这是与上节讨论的非线性弹性材料的区别所在。还应当指出，材料进入塑性后，经卸载再加载时，应力一应变关系沿 CBD 变化。对于硬化材料显然屈服应力提高了。材料进入塑性后，经卸载并反向加载，材料会再次屈服，对于各向同性材料，再次屈服的应力的绝对值大体上与开始卸载时的应力相等。(a) (b)图 7.87.4.2 Mises 屈服准则Mises 屈服准则假定，材料在复杂应力状态下的形状改变能达到单向拉伸屈服时的形状改变能时，材料就开始屈服。于是，Mises 屈服条件可写成11T213221 )()()(式中 为
21、主应力， 是单向拉伸屈服极限。321,T如果定义等效应力 213221 )()()( 则 Mises 屈服条件为（7.40）T在一般应力状态下，等效应力可表示成 212222 )(6)()()( zxyxxzzyyx 引进应力偏量（7.41）cpzycx zxyx式中 是平均应力，则等效应力可用应力偏量表示为3/)(zyxcp（7.42）21222 )(zxyxzyx 若记Tzxyzxyzyx 则可将等效应力简洁地表示成（7.43）21)3(T对于大多数金属材料，实验已经证明，Mises 准则能较好符合实际情况。因此，下面将结合 Mises 准则来讨论。首先，讨论复杂应力状态下的强化规律。假定
22、材料进入屈服后，载荷按小增量方式逐步加载，在一个载荷增量作用下应力和应变都会增加一微小的增量 和 。此时总d应变增量可分成弹性的和塑性的两部分（7.44）Ped弹性应变增量 在卸载后可完全恢复，而塑性应变增量 在卸载后则不能恢复，如ed d图 7.9 所示。12图 7.9与等效应力对应，定义等效应变（7.45）212222 )(3)()()()12 zxyxxzzyyx 对于单向拉伸 ，于是等效应变恰等于拉伸应变0, zxyxzyx 。对应于塑性应变增量的等效应变称为塑性等效应变增量，记作 。由于塑性变形的pd泊桑比 ，于是有5.0 212222 )(3)()()(32 zxpyzxpxpzz
23、pyypxp ddddd (7.46)注意到塑性变形中的体积应变等于零，上式还可表达成（7.47）21222 )(13 zxpyzxypzpyxpp dddd 若记TzxpyzpxypzpypxpP d 212（7.48）则等效塑性应变增量可改写成（7.49）21)(3PTpdd对于一般应力状态，实验资料证明了如下的应变强化规律：材料进入屈服以后进行卸载或部分卸载然后再加载，其新的屈服应力值仅与卸载前的等效塑性应变总量有关。这就是说，只有当等效应力适合（7.50）)(pdH时，重新屈服才会发生。这里的函数 H 反映了新的屈服应力对等效塑性应变总量的依赖关系。13由于单向拉伸时， 就是拉伸应力
24、， 就是拉伸塑性应变增量 ，所以式pdpd（7.50）中的函数关系可以通过单向拉伸的 和 之间的关系来确定，如图 7.10 所示。式（7.50）反映了屈服和强化之间的关系，称为等向强化材料的 Mises 准则。 图 7.107.4.3 Prandtl-Reuss 塑性流动理论如果将等向强化 Mises 准则的式（7.50）写成（7.51）0)(),(PPdHF则 F 可以看成 n 维应力空间的一个曲面，称为屈服面。屈服面的位置决定于当时材料中的等效塑性应变总量。对于金属一类的材料，理论和实验的广泛研究表明，塑性应变增量和屈服面之间存在如下关系（7.52）FdP其中 是一个待定的比例常数。上式可
25、以解释为塑性应变增量“向量”垂直于 n 维应力空间的屈服面，如图 7.11 所示。因此，式（7.52）称为法向流动法则。图 7.11现在来决定式（7.52）中的常数 。为此先计算 。由式（7.42） ，并注意到，就有0 zyxzxzxyyzxxy zzyyx 3,32,2 将上式代入（7.52） ，并注意到式（7.48） ，得1423Pd上式等号两边的模应相等，于是有 249)(TPT由式（7.49）和（7.43） ，从上式可得出 2Pd由于在加载时 应取正值，于是有P这样，向流动法则最终化为（7.53）Pd7.4.4 应力一应变关系如前所述，当应力产生一无限小增量时，总应变增量可分解成弹性的
26、和塑性的两部分，即 Pedd弹性应变增量与应力增量是线性关系，可写成（7.54）)(PeD式中D e 是弹性矩阵。用 前乘上式的两边，得T)(PeTT dDd如果将强化材料的 Mises 准则（7.50）写成微分的形式(7.55)pTH并利用式（7.53） ，则上式可化为PeTeTp dDdDdH 由此可得到等效塑性应变增量 和总应变增量 的关系式p15dDHdeTeP将式（7.53）代入（7.54） ，并应用上式，得到dDHdeTeee 记（7.56）eTeePDHD（7.57）Peep于是得到增量形式的弹塑性应力一应变关系（7.58）dep通常称为弹塑性矩阵。epD顺便指出，对于理想弹塑性
27、材料，上述应力一应变关系（7.58）仍可应用，不过此时取 。0H7.4.5 弹塑性矩阵表达式上节已导出了弹塑性矩阵 的一般表达式，为便于应用。下面对三维空间问题，轴epD对称问题和平面问题分别写出它的显式。1三维空间问题的弹塑性矩阵利用空间问题的弹性矩阵表达式（1.7） ，并注意到 ，容易看出0 zyx TzxyzxyzyxeeD 2223= TzxyxzyxG其中 G 为剪切模量。于是注意到16eTTeDD)(由式（7.56）容易得 2zxyzzxyzzxyzx 22xz2zyx222 )3(9xPGHD(7.59)将式（1.7）和（7.59）代入式（7.57） ，得到空间问题的弹塑性矩阵为
28、 222222 111211 zxzxyzxyzxzxyzx yy xyxzxyxz zzzyyxepED （7.60）式中（7.61）)3(29GH2轴对称问题的弹塑性矩阵我们记 Tzrzr ddd在式（7.60）中划去最后二行二列，可以得到轴对称问题的弹塑性矩阵17（7.62） 222 11211 zrzrzrzr zzzrrrepED式中 仍由式（7.61）决定，但等效应力为2122)(3zrzr3平面问题的弹塑性矩阵对于平面应力问题，应力增量和应变增量分别为 Txyyxdd等效应力为 212)3(xyxyx由此 Txyyx23从而应用平面应力问题的弹性矩阵（2.6）式，有 xyyxGD
29、)1(3将上两式代入（7.56） ，得 2222 )1()(1)(1)( xyxyzyxyzzyzxPQED （7.63）式中（7.64）GHxyyxyx 9)(2)1(222由式（7.57） ，经整理，将平面应力问题的弹塑性矩写成18 222 )1(9)1(211 EHRPPQEDxyyxyxyep（7.65）其中（7.66）2219yxxxyREHP而 Q 的表达式（7.64）可改写为Q)(2对于平面应变问题，材料进入塑性后， ，等效应力可写成12124)(3xyyx这时的弹塑性矩阵可直接由式（7.65）得出，只要将近其中 E 换成 换成 。,1217.4.6 切线模量 的计算H上节导出了
30、几种情况下的弹塑性矩阵 ，它们的表达式中都含有切线模量 ，它epDH是表示材料硬化性能的参数。由式（7.55） ，微分形式的 Mises 屈服准则给出 ，Pd由此可以看出 是等效应力相对于塑性等效应变的变化率。它必须通过单向拉伸试验所H给出的应力一应变由曲线来确定。下面推导 的计算公式。H设单向拉伸试验时，材料进入塑性后的应力一应变关系以某种函数的形式给出，如图 7.12所示。（7.67）)(f图 7.12由于单向拉伸试验给出的是全应变，为了建立与塑性应变之间的关系，必须将它分解为弹19性部分与塑性部分之和，且弹性部分服从虎光定律 ，Ee故有 Pe将其代入按试验曲线给定的函数关系式 )(EfP
31、上式两边取微分，则有 dffdfdPP)(式中 代表应力一应变曲线在 点的斜率。上式经整理可得dEf EfdP1对于单向拉伸，显然 即 ，故 正是对应于复杂应力状态下的 ，代PPHdP入上式则可得到（7.68）EfH1这就是计算 的公式。 的计算步骤是，首先求出等效应变 ，然后在给定材料的H 应力一应变曲线上求得该点所对应的曲线斜率 ，最后代入式（7.68）就可求得 。f H这里有几种情况应该注意，一种情况是，在弹性范围内不能计算 的值，因此时它没有塑性变形，故 应取 值，在编写程序时，可根据计算机的情况取一个适当大的数即可表示 。另一种情况是，对理想弹塑性材料，当材料进入屈服时，由 ，此时可
32、0f由式（7.68）知， 应取 O 值。最后，如果材料硬化程度不剧烈，即材料的应力一应变H曲线较平缓时，则可考虑用以下的近似公式计算 ：H（7.69）f求出后，就可以利用 7.4.5 节的有关公式，计算不同情况下的弹塑性矩阵。7.5 弹塑性问题的求解方法弹塑性问题的有限元求解方法完全是在线弹性问题有限元求解方法的基础上发展起来的，所以本章前面各章所介绍的线弹性问题有限元公式仍然适用。只是在塑性区内应力和应变不再成线性关系，必须采用应力微分和应变微分的关系式（7.58） ，由于式（7.58）右端的系数矩阵 与当时的应力水平有关，所以这个关系是非线性的。epD用有限单元法求解弹塑性问题，也是用适当
33、的方法将问题线性化，然后用一系列的线性解去逼近一个非线性问题的解。而每一个线性解的求解方法和过程和线弹性问题的处理20非常类似。为达到线性化的目的，可以采用增量法（或称增量加载法） 。一般的做法是按比例施加载荷，将结构的弹性极限载荷作为第一个增量，其余的载荷再分成若干等分。如果实际载荷不是按比例施加的，可根据实际情况确定载荷增量。当材料进入塑性后，只要载荷增量适当地小，由式（7.58） ，应力增量和应变增量的关系可近似地表示成（7.70）epD且可认为 仅与加载前的应力水平有关，而与应力和应变的增量无关。这样式（7.70）epD就可视为线性的。7.5.1 增量切线刚度法首先，可以进行一次线弹性
34、的分析，得到弹性极限载荷 下结构的位移，应变和应0F力，分别记为 ， 和 ，在此基础上将载荷分为 n 个增量，并相继作用于结构。00作用载荷增量 时，对于应力处于弹性状态的单元，单元刚度阵用1F(7.71)dvBDkeTe计算。而对于进入塑性的单元，单元刚度阵用(7.72)VepTe计算。 中的应力取增量加载前的应力 。epD0结构总刚度矩阵的形成与前几章所介绍的线弹性问题 方法完全相同，也是利用单元定位向量，直接由单元刚度矩阵进行集成，由此求解 101)(FK可得到 ，从而求出 和 ，这样第一次增量加载后的位移，应变和应力是1101在加载第 k 次载荷增量 后，可通过求解kF（7.73）kk
35、FK1)(得到 ，进而得k21kkkk1重复上述过程，直至第 n 个载荷增量。最后的位移、应变和应力就是弹塑性分析的结果。增量切线刚度法，由于每个载荷增量步都需要重新计算一次刚度矩阵，所以这个方法也称为变刚度法。增量切线刚度法在实施过程中，随着不断地增量加载，解会越来越偏离真实解，即出现所谓解的漂移现象。产生这种现象的原因除了应力一应变关系用增量形式代替微分形式外，还来自过渡单元的不准确的刚度矩阵的计算。必要时可以用下面介绍的方法加以校正。结构在逐步加载过程中，塑性是不断扩展的。有一些单元在某一增量载荷步前处于弹性区内，而在该增量载荷步后可能进入塑性区，这些单元称为过渡单元。对于过渡单元，其刚
36、度矩阵如果简单地按式（7.71）或（7.72）计算都会带来较大误差。通常在计算其刚度矩阵时，应该用下述带权平均弹塑性矩阵 。epD如用 表示过渡单元达到屈服时所需的等效应变增量，而用 表示该单元下次增e s量加载所引起的等效应变增量，如图 7.13 所示，记 sem图 7.13显然，对过渡单元有 。定义加权平均弹塑性矩阵10mepeepDmD)1(则对过渡单元应按下式计算刚度矩阵（7.74）vepTedvBk如果采用数值积分计算 ，则在过渡的积分点上应作上述考虑。ek通常对 的估计，开始往往不够精确。一般第一次估计是根据上次增量载荷的结果s推出，然后用解算的结果来修改 ，经过二三次的迭代就可得
37、到比较精确的结果。s应该指出的是，上述用迭代处理和校正过渡单元刚度矩阵的计算，虽然提高了求解的22精度，但在迭代中却付出了成倍的计算时间。下面介绍另一种较好的求解方案。上述增量切线刚度法，对每一增量加载步，被求解的线性代数方程组（7.73）可以作如下的变动：（7.75）)()( 11vkTkkk dvBF式中 kk FF210是本次增量加载后结构中承受的总载荷向量。而右端最后一项积分是前一增量加载步末与结构中的应力等价的结点载荷。很显然，这是我们在 7.2 节中所提到的带自校正的增量法的方程。由于这一方程能对以前的载荷不平衡作一校正，使求解结果具有较好的精度。因每次增量载荷步，方程（7.75）
38、都具有“自校正”的作用，有人建议可以不去考虑上述过渡单元的迭代计算。实践证明，带自校正的增量切线刚度法（或称一阶自校正法）是一种较好的求解方案。7.5.2 初应力法对于弹塑性问题，增量形式的应力一应变关系可以定义为 0dDde而 P0其中 由式（7.56）计算。对某单元开始进入塑性以后的每次加载，用增量代替微分PD（7.76）PeD00位移增量 应满足的平衡方程是（7.77）)(0RFK由于式（7.76）第一式中的系数阵 为弹性矩阵，式（7.77）中的刚度矩阵 就e 0K是弹性计算中的刚度矩阵。而 vPTdvDBR)(是与初应力 等值的节点载荷，它是不平衡的矫正载荷。0可以看出，式（7.77）
39、中的矫正载荷决定于应变增量 ，而 在求解前是未知的，因此对每个增量载荷步，一个迭代过程是必须的，以便同时求出位移增量和应变增量。第 k 次增量载荷的迭代公式是（7.78）),210(10 jRFKjkjk23第一次迭代是在 下作纯弹性计算。以后的迭代是根据前次000迭代求出的 和加载前的应力水平计算 ，然后按式计算 。逐次进行迭代，jkjkR1jk直至收敛。值得注意的是，对于过渡单元，初应力的计算不应计及总应变增量 中在进入屈服之前的部分（参看图 7.13） ，即矫正载荷应用下式计算（7.79）vPTdvmDBR)1(7.5.3 初应变法对于弹塑性问题，增量形式的应力一应变关系可定义为（7.8
40、0）)(0dde而 P0由式（7.53）和（7.55） ，有（7.81） dHdd TP10用增量代替微分，将式（7.80）和（7.81）线性化（7.82）TPeD1)(00此时位移增量应满足的平衡方程为(7.83)(0RFK式中 仍为弹性计算中的刚度矩阵，而0KdvDBRpvT)(v TeTdvDBH1（7.84）是与初应变等价的节点载荷，或称为矫正载荷。由于矫正载荷与应力增量 有关，而 本身又是待确定的量，因此求解方程组（7.83）必须采用迭代方法。第 k 次增量载荷步的迭代公式为（7.85）),210(10 jkkjkRFK24每一次加载时，首先取 作一次弹性计算，以后的迭代用上次迭代的
41、结果和0加载前的应力水平计算 ，逐次进行迭代，直至收敛。jkR7.5.4 方法的比较1增量切线刚度法是在每次增量加载中用调整刚度矩阵的办法求非线性问题的近似解的。初应力法实质上是在每一增量步上确立结构中弹塑性的应力值和弹性解应力值之差，并不断按弹性方式重新分配这个差值，以恢复平衡。而初应变法则是不断调整初应变的过程。2增量切线刚度法由于在每次增量加载中都形成一次刚度矩阵，并求解一次，故计算时间一般比初应力法和初应变法长。因初应力法和初应变法在迭代中应用的是弹性刚度矩阵，如果求解一开始就将刚度矩阵进行三角分解，则在每次迭代中只需进行右端项的约化和回代，计算时间较节省。3可以证明，对一般强化材料，
42、初应力法的迭代过程总是收敛的，而对初应变法，收敛的充分条件是 。对理想塑性材料，初应力法和初应变法的迭代过程都是发散1/HG的。4如果要为一个通用程序配置一个求解方案的话，建议采用带校正的增量切线刚度法。通常该方法能得到较好的求解精度。如果将增量切线刚度法与初应力法相结合，则效果更好。7.5.5 算例带孔平板的拉伸 问题承受单向拉伸的开孔板条的几何尺寸和载荷情况如图 7.14 所示。利用对称性，计算区域取板条的四分之一部分。采用 Mises 屈服准则，应变硬化材料采用的是常斜率的单轴硬化曲线。图 7.14(b)和(c) 分别表示理想塑性和应变硬化情况下不同载荷水平塑性区的扩展。图 7.14图
43、7.14（d）还示出了一次加上全部载荷，用初应力法得到的解。令人感兴趣的是，用这种单步长解得到的塑性区却与增量方法的结果非常相似。7.6 杆系结构的塑性分析杆系结构的塑性分析相对连续体的塑性分析来说比较简单。例如对超静定桁架来说，若材料是属于理想弹塑性，则可以认为杆件截面应力达到屈服后，该杆件在继续加载时将退出工作，这样就应该修改结构的刚度矩阵，以作新一轮的加载分析。本节主要讨论平面刚架的塑性分析问题。平面刚架塑性分析的直接目标是寻求结构丧失承载能力的极限状态，与确定极限载荷，故又称极限分析。采用以下假设：（1）当出现塑性铰时，假设塑性区退化为一个位于塑性铰处的截面，而其余部分仍为弹性区。25（2）载荷按比例增加，即作用于平面刚架上的全部载荷都可以用一个载荷参数 P 表示，且均为节点载荷，因而塑性铰只能出现在节点处。如有非节点载荷，则可在载荷作用截面处设置节点。（3）每根杆件的极限弯矩相同，但各杆的极限弯矩可不相同。（4）假定结构的位移是微小的，忽略剪力和轴力对极限弯矩的影响。 （5）不考虑弹塑性变形的发展过程。极限载荷的求解过程是通过弹性分析不断探索新的塑性铰形成时对应的载荷增量，并不断修改结构的刚度矩阵和不断施加新的载荷增量，从而最终求得刚架的极限载荷。这
