基于位移控制新法的结构非线性有限元数值模拟

2021-01-21陈梦成李嘉钰

华东交通大学学报 2020年6期

陈梦成，李嘉钰

（华东交通大学土木建筑学院,江西南昌330013）

结构全过程分析是非线性分析理论的基础，它常常应用于大跨度空间结构的设计之中。在结构非线性全过程分析中，刚度矩阵的组装、计算方法及迭代方式的选取，都将对结果产生重大的影响。近年来，人们一直在对结构非线性全过程分析中的迭代控制算法进行研究，提出了压缩平衡迭代法、人工弹簧法、位移控制法、弧长法等多种方法。王志军等[1]在常规的位移控制法和弧长法中引入了简单适用的考虑路径搜寻和单因子加速的方法，能明显改善收敛情况；郑宏等[2]利用Sherman-Morrison 定理，对经典的荷载增量控制法程序做了少量的修补，构建了一个有良好数值特性且仅需很少运算量的位移控制法程序；杜柏松等[3]引入广义位移作为约束方程来确定荷载增量因子，从而控制结构荷载加载的全过程，形成了广义的位移控制法；许茂等[4]针对负刚度结构，采用了附加正刚度弹簧的方法进行了分析，解决了负刚度切线矩阵在有限元分析中造成的迭代求解不收敛的问题；李钢等[5]证明了拟力法在计算结构负刚度时的优势，它是一种高效的非线性分析方法；曲激婷等[6]基于向量式有限元的基本思想，提出了力-位移控制的复合模拟方法，此方法可以很好地求解结构的大变形问题；AL－Aukaily 等[7]采用了直接微分法和向前差分法，同时考虑结构的几何非线性和材料非线性，对求解结构响应所用位移控制法的相关参数进行了敏感度分析。

利用杜修力等提出的位移控制新方法[8]，编写相关计算程序，并结合实际算例，对结构的非线性全过程进行有限元分析，验证位移控制新方法的有效性、正确性和可行性。

1 U.L.列式下增量有限元方程

根据连续性介质力学和虚位移原理[9]，并以时刻结构的状态作为基准参考构型，应变、应力以及位移的增量关系为[10]

式中：KL为整体坐标下的弹性刚度矩阵；KG为整体坐标下的几何刚度矩阵，两者之和为单元切线刚度矩阵；fe为单元抵抗力向量。

2 位移控制新方法

2.1 传统位移控制法

传统的位移控制法（displacement control method，DCM）是求解材料、构件及结构包含极限强度后软化下降段的荷载-位移关系非线性全过程曲线的常用方法之一，它是通过引入给定力作用节点（位移控制点）的位移增量来反求荷载增量和其它未知节点的位移增量。其一般迭代方程为[11]

式中：[KT]为结构的切线刚度矩阵；Δλ 为荷载因子增量；{P}为参考荷载向量；[R]为不平衡力向量；j 为增量步数；i 为迭代次数。接下来根据需要控制的位移向量对结构的迭代平衡方程式（22）进行重新排序，即对结构的切线刚度矩阵、荷载向量以及不平衡力向量进行如下分块

至此，可以计算出每一步迭代中的荷载因子增量Δλ，最终可计算其余点位移分量。

2.2 位移控制新方法

传统位移控制法的求解过程需将有限元方程的刚度矩阵进行分块并重新排列，从而得到新的方程和系数矩阵，但新的系数矩阵将不再具有原刚度矩阵的对称性及带状性，会使计算量显著增加，且计算过程较为复杂，最终带来诸多计算上的不便。

一般地，有限元增量方程可写为

3 算例分析

3.1 几何非线性

选取William 肘式框架[12]进行计算分析，William 肘式框架是众多学者检验其所建立理论有效性及正确性的典型算例，其原因在于该结构具有显著的几何非线性特征，同时又具备可以用来与理论值对比的试验结果，如图1 所示。该结构两端固结，水平方向跨度为657.6 mm，竖直方向高为9.8 mm，杆件为细长矩形等截面直杆，截面宽19.1 mm，高6.17 mm，弹性模量为71 018.5 MPa，框架顶端作用一竖向集中力荷载。

在分析时，将该结构离散成20 个Euler-Bernoulli[13]梁单元进行模拟，取控制点的位移增量步步长u=0.1 mm，总位移增量步n=180 即分析结构位移控制点0 至18 mm 的荷载-位移关系曲线，允许位移收敛误差δ=10-4，在U.L.列式下，使用MATLAB 进行编程计算。

图2 给出了运用传统位移控制法、位移控制新方法得到的数值分析结果与William 试验结果之间的对比；图3 给出了运用传统位移控制法与位移控制新方法得到的结构最终变形曲线。可见，位移控制新方法能跨越荷载-位移关系中的极值点以及结构极限强度后的软化下降段模拟，且能准确地计算出结构的荷载-位移全过程曲线。两种方法所得到的结构最终变形曲线高度吻合，证明了位移控制新方法的可行性和正确性。

图2 William 肘式框架荷载-位移曲线Fig.2 The load-displacement curve for William frame

图3 William 肘式框架变形曲线Fig.3 The deformation curve for William frame

3.2 材料非线性

一混凝土棱柱体试件[14]尺寸为150 mm×150 mm×300 mm，混凝土等级为C30，考虑材料非线性，对该构件进行非线性有限元分析。

计算分析时，采用三维纤维梁单元，沿长度方向离散成3 个单元，截面划分成15×15 的纤维网格，取控制点的位移增量步步长u=0.1 mm，总位移增量步n=24，即分析构件应变为0～0.008 时对应的应力大小，允许位移收敛误差δ=10-4。混凝土的本构模型根据我国混凝土设计规范[15]，采用单轴受压应力-应变曲线。

图4 混凝土柱应力-应变曲线Fig.4 Stress-strain curve of concrete column

将运用传统位移控制法和位移控制新方法得到的数值分析结果与文献[14]结果进行对比，如图4 所示。可见，两种方法计算出的曲线与文献中曲线的吻合程度良好，且两种方法的结果高度吻合，证明了位移控制新方法的准确性和可行性；除此之外，本文在利用MATLAB 进行编程计算时，对两种方法的计算完成时间进行了测试，结果发现，传统位移控制法所需时间约为333 s，而位移控制新方法所需时间约为202 s，新方法近节约了40%的时间，优化了计算过程。