有限元分析：非线性分析求解算法——牛顿迭代法

来源：
熊库辛
https://zhuanlan.zhihu.com/p/417412391

结构有限元分析中，涉及到弹塑性材料、结构大变形大转动以及接触问题都是非线性问题，需要使用非线性工况进行求解。
Abaqus、OptiStruct及Ansys Mechanical等有限元求解器多使用牛顿迭代法，即Newton-Raphson迭代法，来求解非线性问题。
很多结构非线性培训资料都会用下面这张图来讲解牛顿迭代法，但说实话，笔者第一次看到这张图时，只知道是逐步逼近准确解，知其然但不知其所以然，不知道怎么来的。
image.png
牛顿迭代法有以下三种变形：常规牛顿迭代法，修正的牛顿迭代法以及牛顿下山法。有限元求解器一般会根据实际情况自动选择最合适的那种，帮助文档中一般不会讲得很细。
下面请听我一一道来。

1、常规牛顿迭代法

若不存在加速度和阻尼，有限元算法可以简化为求解平衡方程 非线性分析求解算法——牛顿迭代法 - 图2 ，其中刚度矩阵K和载荷F为已知条件，通常需要计算刚度矩阵K的逆矩阵从而求解位移矩阵X。
常规牛顿迭代法如下图所示，注意曲线L(u)并不是已知的，而是逐步求解出来的，下面详细说明求解过程。
image.png
假设分析步总载荷为100N，非线性算法中一般会分多次加载，比如分5次加载每次增加20N，5次加载过程称为5个增量步；每次加载时需要多次迭代才会达到平衡状态，称为迭代步，迭代时使用的算法就是牛顿迭代法。
假设在第n个迭代步时，需要加载大小为 f 的力，下面开始迭代，迭代顺序依次为：A>B>C>D>E。

  • 初始位移为上图中的A点， 对应横坐标为非线性分析求解算法——牛顿迭代法 - 图4 ，根据材料参数和结构形状可以计算出A点对应的刚度矩阵 非线性分析求解算法——牛顿迭代法 - 图5 ，根据方程 非线性分析求解算法——牛顿迭代法 - 图6 可得到位移 非线性分析求解算法——牛顿迭代法 - 图7 ，从而确定B点坐标( 非线性分析求解算法——牛顿迭代法 - 图8 ，非线性分析求解算法——牛顿迭代法 - 图9 )；
  • 根据B点位移 非线性分析求解算法——牛顿迭代法 - 图10 、材料参数以及结构形状可计算出此状态下结构的刚度矩阵 非线性分析求解算法——牛顿迭代法 - 图11 ，根据方程 非线性分析求解算法——牛顿迭代法 - 图12 求得C点坐标( 非线性分析求解算法——牛顿迭代法 - 图13 ，非线性分析求解算法——牛顿迭代法 - 图14 )以及残差 非线性分析求解算法——牛顿迭代法 - 图15 ;
  • 以C点为起点，对应刚度矩阵为 非线性分析求解算法——牛顿迭代法 - 图16 ，加载大小为 非线性分析求解算法——牛顿迭代法 - 图17 的力，可得到D点坐标( 非线性分析求解算法——牛顿迭代法 - 图18 ，非线性分析求解算法——牛顿迭代法 - 图19 )；
  • 根据D点位移 非线性分析求解算法——牛顿迭代法 - 图20 、材料参数以及结构形状可计算此状态下结构的刚度矩阵 非线性分析求解算法——牛顿迭代法 - 图21 ，根据方程 非线性分析求解算法——牛顿迭代法 - 图22 求得E点坐标( 非线性分析求解算法——牛顿迭代法 - 图23 ，非线性分析求解算法——牛顿迭代法 - 图24 )以及残差 非线性分析求解算法——牛顿迭代法 - 图25 ，若残差 非线性分析求解算法——牛顿迭代法 - 图26 小于设定误差值，则认为已经达到收敛，若残差 非线性分析求解算法——牛顿迭代法 - 图27 太大，则继续迭代。

觉得不错的话，点个赞吧，更多有限元相关硬核分享，请点击头像，查看我的主页
@熊库辛
从数学上来说，牛顿迭代法是数值计算方法，具体可翻阅数值分析教材，下面简单说明其数学原理。
设已知方程 非线性分析求解算法——牛顿迭代法 - 图28 有近似根 非线性分析求解算法——牛顿迭代法 - 图29 (假定 非线性分析求解算法——牛顿迭代法 - 图30 )，根据泰勒级数将函数 非线性分析求解算法——牛顿迭代法 - 图31 在 非线性分析求解算法——牛顿迭代法 - 图32 点展开，有
非线性分析求解算法——牛顿迭代法 - 图33
于是方程 非线性分析求解算法——牛顿迭代法 - 图34 可近似地表示为
非线性分析求解算法——牛顿迭代法 - 图35
根据该方程可求出迭代公式 非线性分析求解算法——牛顿迭代法 - 图36 ，其中k=0,1,…，取接近于方程根的一个初始值，然后开始迭代，一般经过几次迭代就能得到一个精度较高的解。
常规牛顿迭代法的优点是收敛速度快。

2、修正的牛顿迭代法

常规牛顿迭代法每次迭代都是求解新的刚度矩阵K以及K的逆矩阵，计算量较大，如果刚度矩阵总是使用初始的 非线性分析求解算法——牛顿迭代法 - 图37 ，并且保持不变，则可以不求解求逆矩阵，大大减少了计算量，这就是修正的牛顿迭代法，示意图如下。
image.png

3、牛顿下山法

常规牛顿迭代法要求初始值必须在准确值附近才会收敛，初始值不合适可能会导致结果不收敛。
为克服收敛问题，引入判断条件 非线性分析求解算法——牛顿迭代法 - 图39 ，满足此要求的算法称为下山法。
常规牛顿迭代法迭代公式为
非线性分析求解算法——牛顿迭代法 - 图40 (1)
若对计算结果进行加权平均作为新的改进值，
非线性分析求解算法——牛顿迭代法 - 图41 (2)
其中 非线性分析求解算法——牛顿迭代法 - 图42 ，称为下山因子。
将(1)式代入(2)式有
非线性分析求解算法——牛顿迭代法 - 图43 ， 非线性分析求解算法——牛顿迭代法 - 图44
选择下山因子时从 非线性分析求解算法——牛顿迭代法 - 图45 开始，逐次将 非线性分析求解算法——牛顿迭代法 - 图46 减半进行试算，直到能使下降条件 非线性分析求解算法——牛顿迭代法 - 图47 满足为止。
通过引入下山因子，保证了牛顿迭代法的收敛性。
因此牛顿下山法收敛速度快，一定收敛，而且对初始值要求不高，相比常规牛顿迭代法更稳健。

更多有限元相关硬核分享，请点击头像，查看我的主页
@熊库辛
码字不易，点个赞吧。
image.png
参考文献：
[1]李庆扬，王能超，易大义. 数值分析[M]. 第5版. 北京: 清华大学出版社, 2008:222-226.
[1]曾攀. 有限元分析及应用[M]. 北京: 清华大学出版社, 2004:325-327.

编辑于 2021-10-24 09:16