从一实例分析看弧长法与牛顿—拉普森法

魏鹏，李建光，李延强，陈超



从一实例分析看弧长法与牛顿—拉普森法

魏鹏，李建光，李延强，陈超

（西南林业大学土木工程学院，云南 昆明 650000）

在求解的非线性有限元方程中，牛顿—拉普森法和弧长法是两类重要的方法。牛顿—拉普森迭代法只能跟踪位移载荷曲线的上升段，但无法跟踪极值点以后的位移—载荷路径，而弧长法可以全程跟踪位移—载荷路径。对牛顿—拉普森法和弧长法的原理以及实施步骤进行了回顾，然后通过MATLAB对一则本构关系为非线性的算例使用牛顿—拉普森法和弧长法进行了计算并与理论解做了比较。数值结果表明，弧长法能很好跟踪全过程位移—载荷路径，并取得较好的效果，而牛顿—拉普森法在跟踪到极值点时产生发散而无法继续跟踪极值点以后的路径。

牛顿—拉普森法；弧长法；非线性有限元；MATLAB

1 牛顿—拉普森

对于非线性问题得到的非线性方程组，一般可表示为：

（）=0. （1）

写成平衡方程形式：

（）=. （2）

式（2）中：（）为切向刚度矩阵；为位移矢量；为施加的载荷向量。

平衡迭代的过程表示为：

{i+1}={i}[i]-1（－i）. （3）

式（3）中：[i]为方程的Jacobian矩阵（即刚度矩阵）；{}为结点力矢量；{i}为内力矢量；{－i}为不平衡力矢量。

牛顿—拉普森迭代步骤如下：①基于i时的结构构型计算i和i；②计算不平衡力矢量{－i}；③由i和不平衡力矢量{－i}计算位移增量；④更新位移向量△i；⑤重复②到④的过程直至计算收敛为止。

2 弧长法

弧长法的约束方程为：

{△}T{△}+2△2{}T{}=△2. （4）

式（4）中：为载荷比例系数；△为载荷增量；△为弧长半径。

根据的不同，可以分为不同类型的弧长法：=1，球面弧长法；=0，柱面弧长法；等于当前刚度参数值，椭球面弧长法。

弧长法的平衡迭代方程为：

其中：

以柱面弧长法为例推导，此时约束方程为：

{△}T{△}=△2. （9）

由位移增量关系可得：

为求得，把（6）式代入（10）式得：

把（12）代入（9）式可得关于一元二次方程：

式（12）中系数分别为：

对于收敛准则，一般采用不平衡力准则[3]。表达式为：

3 算例分析比较及讨论

本文的算例如下：

受拉杆施加的力=20 kN，杆长=50 cm，截面积=2 cm2，材料的应力应变关系为：

=0（1－/0）. （17）

式（17）中：0=0.002，0=21 000 kN/cm2。

分别用牛顿—拉普森法和弧长法计算位移—载荷曲线，并与理论解做比较。理论解通过把应力—应变曲线转化为载荷—位移曲线得到。过程如下：

=0（1－/0）；=/.

由以上两式，并代入相关参数可得：

=﹣8 4002+840. （18）

下面利用牛顿—拉普森迭代法和弧长法（分别如图1和图2所示），借助MATLAB对算例进行了实现，具体如图3、图4所示。

图1 牛顿—拉普森迭代法

图3 牛顿—拉普森法与理论解

图4 弧长法与理论解

由图1可知，用牛顿—拉普森迭代法计算得到的载荷—位移曲线与理论解的上升段是非常吻合的，但是由于在极值点附近出现发散，而不能继续跟踪曲线的下降段，但在跟踪载荷—位移曲线的上升段时还是有效的。由图2可知弧长法计算所得的载荷—位移曲线与理论解的全过程的逼近程度都是很好的。弧长法不仅可以跟踪曲线的上升段，还可以跟踪曲线中的下降段。因此，弧长法可以更好地跟踪载荷—位移曲线的特点。

［1］Carrera E.A study on arc-length-type methods and their operation failures illustrated by a simple model［J］.Computers & Structures，1994，50（2）：217-229.

［2］Zhou Z，Murray D W. An incremental solution technique for unstable equilibrium paths of shell structures［J］. Computers & Structures，1995，55（5）：749-759.

［3］刘国，卓家寿，夏颂佑.求解非线性有限元方程的弧长法及在工程稳定分析中的应用［J］.岩土力学，1993（4）：57-67.

［4］李元齐，沈祖炎.弧长控制类方法使用中若干问题的探讨与改进［J］.计算力学学报，1998，15（4）：414-422.

［5］Bellini P X，Chulya A.An improved automatic incremental algorithm for efficient solution of nonlinear finite element equations［J］.Computers & Structures，1987，26（1）：99-110.

2095－6835（2018）24－0005－02

TM402

A

10.15913/j.cnki.kjycx.2018.24.005

魏鹏（1991—），男，研究方向为非线性。

〔编辑：严丽琴〕