固体力学学报

序号	支承情况	理论屈曲荷载（N）	本文	Yang和Kuo
（a）	两端铰接	246.45	247.21	0.31%	246.70	0.10%
（b）	一端固接，一端铰接	502.96	506.81	0.77%	504.72	0.35%
（c）	一端固接，一端自由	61.61	61.71	0.16%	61.68	0.11%
（d）	两端固接	985.80	994.95	0.93%	987.01	0.12%
（e）	两端铰接，中间铰支承	985.80	994.95	0.93%	987.01	0.12%

序号	支承情况	理论屈曲荷载（N）	本文	Yang和Kuo
（a）	两端铰接	246.45	247.21	0.31%	246.70	0.10%
（b）	一端固接，一端铰接	502.96	506.81	0.77%	504.72	0.35%
（c）	一端固接，一端自由	61.61	61.71	0.16%	61.68	0.11%
（d）	两端固接	985.80	994.95	0.93%	987.01	0.12%
（e）	两端铰接，中间铰支承	985.80	994.95	0.93%	987.01	0.12%

方法	增量步数	总迭代次数	计算时间（s）
本文	8545	13593	4.39
Yang和Kuo	8545	13612	4.97

方法	增量步数	总迭代次数	计算时间（s）
本文	8545	13593	4.39
Yang和Kuo	8545	13612	4.97

方法	增量步数	总迭代次数	计算时间（s）
本文	4282	8792	2.27
Yang和Kuo	4287	8803	2.73

方法	增量步数	总迭代次数	计算时间（s）
本文	4282	8792	2.27
Yang和Kuo	4287	8803	2.73

方法	增量步数	总迭代次数	计算时间（s）
本文	21329	43020	10.62
Yang和Kuo	21445	43327	11.73

方法	增量步数	总迭代次数	计算时间（s）
本文	21329	43020	10.62
Yang和Kuo	21445	43327	11.73

方法	增量步数	总迭代次数	计算时间（s）
本文	516	1035	1.22
Yang和Kuo	518	1049	1.36

方法	增量步数	总迭代次数	计算时间（s）
本文	516	1035	1.22
Yang和Kuo	518	1049	1.36

方法	增量步数	总迭代次数	计算时间（s）
本文	389	980	0.34
Yang和Kuo	389	1949	0.55

方法	增量步数	总迭代次数	计算时间（s）
本文	389	980	0.34
Yang和Kuo	389	1949	0.55

基于刚体准则的梁单元几何刚度矩阵推导及应用

PDF下载

陈安全 ¹^,²^,^** , 吴俊力 ¹ , 伏梨莹 ¹ , 范秀娟 ¹

固体力学学报 | 研究论文 2025,46(4): 546-559

收起

固体力学学报 | 研究论文 2025, 46(4): 546-559

基于刚体准则的梁单元几何刚度矩阵推导及应用

全屏

陈安全¹^,²^,^**, 吴俊力¹, 伏梨莹¹, 范秀娟¹

作者信息

¹内江师范学院建筑工程学院，内江，641100

²重庆大学土木工程学院，重庆，400045

通讯作者:

**E-mail：chen.anquan@foxmail.com.

Derivation and Application of the Geometric Stiffness Matrix for Beam Elements Based on the Rigid Body Rule

Anquan Chen¹^,²^,^**, Junli Wu¹, Liying Fu¹, Xiujuan Fan¹

Affiliations

¹School of Architecture and Engineering, Neijiang Normal University, Neijiang, 641100

²School of Civil Engineering, Chongqing University, Chongqing, 400045

出版时间: 2025-08-27 doi: 10.19636/j.cnki.cjsm42-1250/o3.2025.017

文章导航

摘要

收起

杆系结构的承载性能与其几何刚度密切相关，通常情况下梁单元几何刚度通过稳定函数或能量原理推导得出，推导过程较为复杂，并且对于位移高阶项的处理缺乏物理观. 研究首先基于三维实体梁线性化的增量虚功方程和空间力矩转动特性，分析得到梁单元节点保持平衡的关键特征-诱发矩矩阵，并证明了单元层面非对称的诱发矩矩阵组装到结构层面后将恢复对称性. 然后根据诱发矩矩阵和单元位移特征构造了含待定系数的三维梁单元几何刚度矩阵，并利用几何刚度矩阵的对称性和刚体准则推导出了简洁形式的三维梁单元几何刚度矩阵的显示表达，同时简化得到二维梁单元几何刚度矩阵. 通过对典型案例进行线性屈曲分析和非线性分析，结果表明所推导的几何刚度矩阵可以很好的用于梁类型结构的屈曲和后屈曲分析，并且具有较好的计算精度和较高的计算效率. 此推导梁单元几何刚度矩阵的方法物理意义清晰，推导过程简单，是推导单元几何刚度的一种新思路.

关键词

梁单元 / 几何刚度矩阵 / 对称性 / 刚体准则 / 线性屈曲分析 / 几何非线性分析

Abstract

收起

The load-bearing capacity of slender structural systems, particularly those comprising assemblages of beam elements, is critically governed by their geometric stiffness properties. Conventional methodologies for deriving the geometric stiffness matrix of beam elements are typically rooted in stability functions or variational energy principles and entail mathematically intricate formulations that frequently obscure physical interpretations, especially concerning the treatment of higher-order displacement terms. In marked contrast to these established approaches, the present study introduces a physically insightful framework by first investigating the induced moment matrix. This fundamental mechanical attribute rigorously preserves nodal moment equilibrium during finite three-dimensional rotations. Its formulation is established through the incremental virtual work principle, derived from consistent linearization of three-dimensional solid beam kinematics integrated with exact spatial rotational transformations. A pivotal theoretical finding demonstrates that while the induced moment matrix inherently exhibits asymmetry when examined at the individual element level, this asymmetry vanishes upon assembly into the global structural system, resulting in a symmetric structural-level geometric stiffness matrix. Building upon this foundation, a three-dimensional geometric stiffness matrix incorporating undetermined coefficients is systematically constructed by rigorously enforcing displacement compatibility conditions across the element. Subsequently, by exploiting the proven symmetry of the assembled global geometric stiffness matrix and strictly imposing the rigid body rule, which necessitates zero straining energy for arbitrary rigid displacements, a concise and fully explicit analytical expression for the three-dimensional beam element geometric stiffness matrix is derived. This expression is further simplified to yield its corresponding two-dimensional counterpart for planar analyses. Comprehensive numerical validations encompass eigenvalue buckling analyses of axially compressed prismatic members and geometrically nonlinear analyses of diverse structural configurations: cantilever beams, arches, spatial frames, and curved beams. These extensive simulations consistently confirm that the proposed geometric stiffness matrix achieves exceptional computational accuracy and numerical efficiency in predicting both critical buckling loads and complex post-buckling equilibrium paths. This work establishes a novel, mechanically transparent, and mathematically streamlined paradigm for deriving geometric stiffness matrices, offering a unifying perspective readily extensible beyond beam elements to potentially encompass plate and shell finite elements in future developments.

Key words

beam element / geometric stiffness matrix / symmetry / rigid body rule / linear buckling analysis / geometric nonlinear analysis

引用本文

陈安全, 吴俊力, 伏梨莹, 范秀娟. 基于刚体准则的梁单元几何刚度矩阵推导及应用. 固体力学学报, 2025 , 46 (4) : 546 -559 . DOI: 10.19636/j.cnki.cjsm42-1250/o3.2025.017

Anquan Chen, Junli Wu, Liying Fu, Xiujuan Fan. Derivation and Application of the Geometric Stiffness Matrix for Beam Elements Based on the Rigid Body Rule[J]. Chinese Journal of Solid Mechanics, 2025 , 46 (4) : 546 -559 . DOI: 10.19636/j.cnki.cjsm42-1250/o3.2025.017

正文

收起

0　引言

收起

柔性梁为众多工程与科学领域中的基本结构单元，例如钢框架结构、航天飞行器天线、机械臂、风力发电机等细长杆件类型的结构物均可视为由柔性梁组成，这类结构由于刚度较小，在正常服役阶段表现出明显的几何非线性特征^[1]. 近年来，有关梁结构的非线性分析数值方法得到了大力发展，由于计算便捷高效并且适用性广，使其逐渐成为了梁结构非线性分析的主要工具^[2].

对于梁结构几何非线性分析问题，主要的发展方向有梁柱理论和非线性有限元. 基于梁柱理论的稳定函数单元和稳定插值函数单元是框架几何非线性分析中一种常用的单元模式. 前者采用由梁柱理论建立的单元内力-变形关系导出单元结点力和刚度矩阵的计算方法，后者将梁柱理论和有限元方法相结合，根据平衡微分方程的解构造单元^[3,4]. 非线性有限元法是根据杆件单元的受力性状、变形特点等选择合适的位移函数、应变及应力表达式再利用虚功原理推导单元的切线刚度矩阵^[5]. 非线性有限元方法可方便地考虑结构的初始缺陷、截面翘曲等因素，能导出直梁、曲梁和薄壁梁等具有不同特性的梁单元. 因此，非线性有限元方法是当前杆系结构非线性分析的主要方法. 最简单的考虑几何非线性效应的有限单元模型是在使用三次Hermite插值函数的单元弹性刚度矩阵基础上，附加一个几何刚度矩阵来进行分析. 众多的研究者如Barsoum和Gallagher^[6]，Meek和Tan^[7]，Chan和Kitipornchai^[8]，Yang和Kuo^[9]等均对该模型进行了研究. 对于三次Hermite函数插值精度不够的问题，可以通过插入附加节点，构造更高次Hermite插值函数梁单元来提高数值计算精度^[10]. 夏拥军和陆念力^[11]、夏拥军和缪谦^[12]证明了该方法的有效性，其在梁单元中点处增加一个节点，通过静力凝聚法消去中间节点的自由度，得到了与传统两节点梁单元相同自由度数量的梁单元刚度矩阵. 除增加附加节点外，还可以通过约束单元位移场满足一定的条件来增加拟合项，进而提高梁单元精度. 例如，Chan和Zhou^[13,14]提出了一种PEP（Pointwise Equilibrating Polynomial）单元，该单元不仅要求形函数满足边界上的转角和位移条件，而且要求在单元的中点满足横向挠度微分方程及其导数，从而使形函数多项式提高到了6项. 聂国隽和钱若军^[15]、许红胜等^[16]采用类似的思路，基于Vlasov的薄壁梁变形截面假设，建立了一种新的空间薄壁梁单元刚度矩阵用于几何非线性分析. 郭鑫等^[17]基于空间Timoshenko梁的解析位移形函数，推导改进空间梁单元刚度矩阵. Zhao等^[18]考虑截面扭转和翘曲变形，推导出了14自由度的三维梁柱单元的二阶精确刚度矩阵. 邓继华等^[19]基于共旋法与子结构法推导了具有对称特征的几何非线性平面梁单元切线刚度矩阵用于非线性分析. 近期，Liu等^[20,21]考虑I型梁单元截面扭曲变形建立了一种18自由度的空间梁单元，并采用线性、三次多项式插值函数推导了单元弹性刚度矩阵和几何刚度矩阵，其推导过程特别注重结构转角节点的平衡.

在梁单元刚度矩阵推导过程中，除由单元弹性应变能推导出弹性刚度矩阵外，还可以由势能推导出几何刚度矩阵，其由单元初始节点力构成，代表的是非线性分析中单元初始节点力的不稳定效应. 对于平衡路径包含多个临界点的结构，采用增量迭代法进行非线性分析时必须同时将弹性刚度和几何刚度用于预测阶段，否则在路径曲率较大位置处发散^[22]. 传统方法根据插值函数推导梁单元几何刚度矩阵需考虑非线性应变的高阶项，推导过程繁冗，还可能因高阶项取舍不当而影响几何刚度矩阵的合格性. 本文将从三维实体梁单元理论出发，根据空间力矩节点平衡的特征，构建三维梁单元的几何刚度矩阵，然后利用对称性和刚体合格性检验，推导出刚体合格的三维梁单元和二维梁单元几何刚度矩阵的显示表达形式，最后结合数值案例验证该几何刚度矩阵的有效性.

1　三维梁单元非线性平衡方程

收起

三维实体梁在非线性变形过程中有三种位形状态，即C₀是初始（未变形）位形，C₁是上一步计算（已知）的位形，C₂是当前计算（未知）的位形. 根据更新的拉格朗日列式，参考于C₁位形，三维实体梁线性化的增量虚功方程可以写为^[23]：

(1)

(2)

(3)

其中¹V是梁在C₁位形处的体积，¹S是梁在C₁位形处的表面积，₁C_ijkl是本构张量，₁e_ij是梁从C₁位形到C₂位形的格林-拉格朗日应变增量张量的线性部分，¹τ_ij是梁在C₁处的柯西应力张量，₁η_ij是格林-拉格朗日应变增量张量的非线性部分，

和

分别表示作用在梁C₂和C₁位形处的表面力

和

所做的虚功. 本文中，物理量的左上标表示该物理量所处位形，左下标表示参考位形.

建立如图1所示的2节点12自由度三维实体梁单元，可以得出梁单元的增量刚度方程^[23]为：

(4)

其中[k_e]是弹性刚度矩阵，[k_g]是由势能推导出的几何刚度矩阵，[k_i]为诱发矩矩阵，{u}为单元位移向量，

}和

分别为梁单元在C₂和C₁位形处的节点力. 诱发矩矩阵[k_i]具体表达如下：

(5)

其中子矩阵[k_ia]定义为

(6)

[k_i]_b的表达形式与[k_i]_a相同，只是将表达式中的单元节点a的力矩转换成节点b的力矩即可.

诱发矩本质上是由三维梁单元的初始节点力矩经历三维旋转引起，并且由于旋转特性不同，圣维南扭矩

被称为半正切矩（Semitangential moments），弯矩

，

称为类正切矩（Quasitangential moments），这是在大变形分析中将梁单元所有变量都定义在C₂位形处节点能够保持平衡的关键，在后续几何刚度矩阵的推导中特别注意了这一特征. 应当注意的是，由于尚未考虑与其他相邻单元的兼容性，上述推导出的是离散的单个单元刚度方程，推导出的刚度矩阵在单元层面上是不对称的，如式（5）中的诱发矩矩阵[k_i]所表示的. 当单元与结构的相邻单元组装在一起时，刚度矩阵将在结构水平上恢复对称性，这在下一节中将得到证明. 可以看出[k_i]也是单元节点初始力矩的矩阵，也代表了节点初始力矩的不稳定特征，由势能推导出的刚度矩阵表示为几何刚度矩阵[k_g]是针对一般单元的叫法，例如桁架单元，二维梁单元. 对于三维梁单元来说，一般意义上的几何刚度矩阵应该是囊括了所有初始节点力（矩）表示的矩阵，故在后续的三维梁单元几何刚度矩阵的构造中特别关注了空间力矩转动产生的诱发矩特征.

2　三维梁单元刚体合格的几何刚度矩阵推导

收起

2.1　诱发矩矩阵在结构层面上的对称性

三维梁单元节点b端诱发矩矩阵[k_i]_b是一个非对称矩阵，可以分解为对称部分[S]和反对称部分[A]：

(7)

其中

(8)

(9)

考虑如图2中连接两个单元i和j的节点B，单元i在节点B（xⁱ=L）的反对称部分可以表示为：

(10)

同样地，单元j在节点B（x^j=0）的反对称部分可以表示为：

(11)

对于图2所示的结构节点B，应满足以下弯矩的平衡条件：

(12)

其中[γ]是转动矩阵^[24]，

(13)

其中l_p，m_p，n_p（p=1，2，3）是两组坐标系（x^j，y^j，z^j）和（xⁱ，yⁱ，zⁱ）之间的方向余弦，将式（13）代入到式（12）可得：

(14a-c)

在结构层面上，两个单元在节点B的诱发矩矩阵[k_i]的反对称部分[A]可以组装为：

(15)

将式（10），（11）和（13）代入到式（15），使用式（14）关系式，得到：

(16)

其中

(17)

式（17）中的系数表达如下：

(18)

根据转动矩阵[γ]的正交性和单位性^[24]，式（18）中所有的系数a₁₁，a₁₂，…，a₃₃都等于0，所以式（17）中A₁₂，A₁₃，A₂₃也都等于0，因此在节点B处相连的单元诱发矩矩阵的反对称部分在结构层面组装抵消，也就是

(19)

式（19）说明了在一个节点处所有单元的诱发矩矩阵的反对称部分将抵消，仅留下对称部分，因此结构的整体刚度矩阵（包括弹性，几何刚度矩阵和诱发矩矩阵）在结构层面上依然保持对称性，在后续的推导中将利用这一特性.

2.2　三维梁单元几何刚度矩阵的构造

在三维梁单元的增量刚度矩阵方程式（4）中，几何刚度矩阵[k_g]与诱发矩矩阵[k_i]均代表了几何非线性分析中初始节点力的不稳定效应，为统一表述，在本节中将两者统称为一般意义上的单元几何刚度矩阵，它应该包含诱发矩矩阵的特征，因为这是结构在变形后C₂位形处保持节点平衡的关键，所以在此处构造含待定参数的三维梁单元刚体合格的几何刚度矩阵为

(20)

式（20）[k_g]_r中的子矩阵[k_i]_a和[k_i]_b反映的是初始节点力矩经历三维旋转产生的诱发矩.

式（20）中还包含14个未知的子矩阵，首先考虑对称性. 在2.1节中已经证明[k_i]_a和[k_i]_b在单元层面上虽然不对称，但经过同一节点处单元组装在结构层面上将恢复对称性，故在此处不用考虑子矩阵[k_i]_a和[k_i]_b的对称性，但其他位置应该保持对称性. 故式（20）可以简化为

(21)

式（21）中已考虑了[k_g]_r的对称性，还包含有8个未知的子矩阵，在2.3节中将根据刚体准则来进一步求解.

2.3　三维梁单元几何刚度矩阵的求解

Yang和Chiou提出的刚体准则^[25]，可以用于测试非线性有限元几何刚度矩阵的合格性，其内容为：对于在C₁位形处有初始节点力作用而保持平衡的单元，在发生刚性位移后，该单元在新位形C₂处仍应保持平衡，并且作用在单元上的节点力大小不变，方向随刚性位移转动. 由此如果[k_g]_r是刚体合格的，则应满足下列等式

(22)

其中{u}_r是刚性位移，[T]是从C₁位形到C₂位形的转换矩阵，[I]是单位矩阵. 接下来将利用式（22）对[k_g]_r进行求解.

①首先考虑单元发生沿x，y，z轴的刚性平移，则

(23)

(24)

将式（24）代入式（22）中可得

(25)

将式（21）和式（23）代入式（25）中可得

(26)

故式（21）可以进一步简化为

(27)

式（27）中实际仅包含[k₁₁]，[k₁₂]，[k₁₄]，[k₂₄]4个未知子矩阵，下面将根据绕x，y，z轴的刚性转动来求解未知子矩阵. 为方便描述，假定这4个子矩阵分别为

(28)

需要注意的是，根据[k_g]_r的对称性，[k₁₁]必然是对称矩阵，则x₁₂=x₂₁；x₁₃=x₃₁；x₂₃=x₃₂. 整个式（28）共有33个未知量需要求解，下面将考虑发生微小刚性转动来求解这些未知量.

②单元绕x轴发生微小转动θ_r，此时刚性位移为

(29)

转动矩阵为

(30)

其中

(31)

当发生微小转动时，θ_r≈0，cosθ_r≈1；sinθ_r≈θ_r，式（31）化简为

(32)

单元初始节点力向量为

(33)

将式（27），（29），（30），（33），代入到式（22）中可得

(34)

③单元绕y轴发生微小转动θ_r，此时刚性位移为

(35)

其中L为梁单元长度，转动矩阵形式仍为式（30），只是子矩阵变为

(36)

当θ_r很小时，进一步化简为

(37)

将式（27），（30），（33），（35）代入到式（22）中可得

(38)

④单元绕z轴发生微小转动θ_r，此时刚性位移为

(39)

转动矩阵仍为式（30），此时子矩阵变为

(40)

当θ_r很小时，进一步化简为

(41)

将式（27），（30），（33），（39）代入到式（22）中可得

(42)

根据式（34），（38）和（42）中关系，可以得到求解[k_g]_r子矩阵独立的表达式23个，其余未定系数均取为0，即t₂₂=0，t₂₃=0，t₃₂=0，t₃₃=0，x₁₁=0，y₁₁=0，y₁₂=0，y₁₃=0，z₁₂=0，z₁₃=0，可得式（28）中的四个子矩阵为

(43)

将式（43）中4个子矩阵代入式（27）中可得三维梁单元刚体合格的几何刚度矩阵[k_g]_r的表达形式.

3　二维梁单元几何刚度矩阵

收起

忽略式（27）三维梁单元几何刚度矩阵中有关z轴方向线位移和x轴方向、y轴方向转角位移的行和列，可得6×6的刚度矩阵，即为二维梁单元的几何刚度矩阵.

(44)

其中

(45)

4　数值案例

收起

本节将对典型结构进行线性特征值屈曲分析和结构非线性分析来说明本文建立的梁单元几何刚度矩阵的有效性和适用性. 线性屈曲分析的有限元方程为（

，其中[K_e]是结构的弹性刚度矩阵，由梁单元弹性刚度矩阵组装得到，

是施加单位外部荷载所得到的结构几何刚度矩阵，{U}是结构的位移向量，P_cr为结构的屈曲临界荷载. 结构非线性分析采用增量迭代方法进行计算，预测阶段采用[K_e]+[K_g]_r，校正阶段单元恢复力计算采用

，{Δu }为单元位移增量，所使用的弹性刚度矩阵[k_e]见文献[9]，路径追踪技术采用更新正交迭代方法^[26]. 为了对比说明本文所建立的梁单元几何刚度矩阵在数值分析中的计算效率，案例分析中将与Yang和Kuo^[9]建立的简化理论梁单元进行对比分析. 本文的案例分析基于MATLAB R2020a平台，计算机配备第12代Intel © Core i5-12400F处理器和32.0 GB内存，采用Windows11 64位操作系统.

4.1　轴心受压柱屈曲分析

分析如图3所示不同边界条件下轴心受压柱的屈曲载荷，柱的长度L=10 m，材料弹性模量E=1.2×10⁶ N/m²，柱截面惯性矩I=2.083×10^-3 m⁴，每根柱划分为20个二维梁单元，采用特征值屈曲分析计算结果如表1所示.

从表1可以得出，采用本文推导出的二维梁单元几何刚度矩阵进行不同支承条件下轴心受压柱的线性特征值屈曲分析获得的屈曲荷载与理论屈曲荷载P_cr=π²EI/（μL）²基本吻合，最大误差不超过1%，由此证明本文建立的二维梁单元几何刚度矩阵在线性特征值屈曲分析中的正确性和有效性. 同时与Yang和Kuo^[9]的简化理论梁单元特征值屈曲荷载对比，本文所得屈曲荷载略大一点，主要原因是本文所建立的几何刚度矩阵主要考虑的是单元宏观层面上应满足的物理特征，未考虑单元内部变形，造成结果略微偏大，但精度也足以满足实际工程使用.

4.2　悬臂梁非线性分析

如图4所示的悬臂梁受到自由端横向集中力P的作用，材料弹性模量E=1.2×10⁶ N/m²，梁横截面面积A=0.1 m²，截面惯性矩I=2.083×10^-3 m⁴，梁的长度L=10 m，悬臂梁划分为10个二维梁单元，结构非线性分析结果如图5所示.

从图5可以看出，采用本文建立的二维梁几何刚度矩阵对承受横向荷载悬臂梁进行几何非线性分析所得的非线性平衡路径与Yang和Kuo^[9]、Duy等人^[27]的计算结果非常吻合，也与解析解一致，证明本文所建立的几何刚度矩阵可以很好的用于二维梁结构的几何非线性分析. 同时通过表2中的计算数据可以发现，本文的梁单元相比于Yang和Kuo^[9]的简化理论梁单元在增量步数一致的情况下，总迭代次数略少，相应的计算时间也更少，具有更高的计算效率.

4.3　固接-铰接拱非线性分析

如图6所示深圆拱，左端固接，右端铰接，并在圆拱顶点处作用竖直集中荷载P. 圆拱半径R=100 cm，圆心角α=215⁰，拱截面采用正方形截面，边长为1 cm，弹性模量E=1200 kN/cm²，剪切模量G=600 kN/cm². 本文采用30个二维梁单元来求解深圆拱非线性行为，求解结果如图7所示.

从图7可以看出，本文所建立的几何刚度矩阵可以很好的模拟深圆拱的几何非线性行为，并且与Yang和Kuo^[9]、Wood和Zienkiewicz^[28]的分析结果以及解析解^[29]一致，说明了本文所建立的二维梁单元几何刚度矩阵对于复杂非线性问题同样具有很好的适用性. 同时通过表3中的计算数据可以发现，本文的梁单元相比于Yang和Kuo^[9]的简化理论梁单元，在设置相同的初始参数下，增量步数、总迭代次数略少，总体计算时间更少，计算效率更高.

4.4　空间直角框架线性屈曲分析和非线性分析

如图8所示的对称直角框架在a、c两个铰接端受到大小相同的弯矩作用，即M_za=-M_zc，框架在a端和c端处允许发生沿x轴平移和绕z轴转动的位移，而节点b仅能在yz平面内移动. 材料弹性模量E=71240 N/mm²，剪切模量G=27190 N/mm²，杆长L=240 mm，横截面面积A=18 mm²，截面惯性矩I_y=0.54 mm⁴，I_z=1350 mm⁴，扭转常数J=2.16 mm⁴. 在非线性分析时，在顶点b处引入大小为F_zb=5×10^-5M_za的干扰力. 此案例由于涉及较大的空间旋转与扭转变形，常用作于空间框架分析程序的测试案例.

首先进行线性屈曲分析，利用本文所建立的几何刚度矩阵，得到直角框架的屈曲临界弯矩值为624.76 N·mm，与理论屈曲弯矩^[30]

N·mm非常吻合. 非线性分析得到的直角框架荷载位移曲线如图9所示，与Yang和Kuo^[9]、Argyris等人的解^[30]也吻合较好，并且可以得到屈曲临界弯矩值为625.78 N·mm，与线性屈曲临界值和理论屈曲弯矩值均吻合很好. 此案例也说明了本文所建立的三维梁单元几何刚度矩阵的正确性和有效性. 同时通过表4中的计算数据可以发现，本文的三维梁单元相比于Yang和Kuo^[9]的简化理论三维梁单元，计算的增量步数、总迭代次数和总体时间更少，因此计算效率更高.

4.5　三杆空间梁结构非线性分析

如图10所示的三杆空间梁结构，在d端作用有大小相等的两个力F_xc、F_zc，a端固定在支座上，三根直梁彼此垂直连接，长度L=1 m，横截面面积A=0.01 m²，截面扭转常数J=1.406×10^-5 m⁴，截面惯性矩I_y=I_z=8.333×10^-6 m⁴，材料弹性模量E=1×10⁶ N/m²，剪切模量G=5×10⁵ N/m². 每根杆划分为10个三维梁单元，空间梁结构非线性分析结果如图11所示.

从图11可以看出，采用本文所建立的三维梁单元几何刚度矩阵计算所得的结构平衡路径与Yang和Kuo^[9]、Romero^[31]所得的结构非线性路径非常吻合，说明了本文推导的三维梁单元几何刚度矩阵可以很好的用于空间梁结构的非线性分析. 进一步通过表5中的数据发现，本文的三维梁单元相比于Yang和Kuo^[9]的简化理论三维梁单元，计算的增量步数、总迭代次数和总体花费时间更少，体现出计算效率更高.

4.6　45°空间悬臂曲梁非线性分析

如图12所示的45°悬臂曲梁在自由端作用有z向的集中荷载F，曲梁初始位于xy平面内，横截面是边长为1 cm的正方形，半径R为100 cm. 该案例涉及空间曲梁的大变形，具有较强的非线性特征，故常作为空间曲梁非线性分析的标准案例. 杆件截面扭转常数J=1.406×10^-5 cm⁴，材料弹性模量E=1×10⁷ N/cm²，剪切模量G=5×10⁵ N/cm². 曲梁划分为10个三维梁单元，非线性分析结果如图13所示.

从图13可以看出，本文计算所得的曲梁非线性路径与Yang和Kuo^[9]简化理论三维梁单元计算结果相同，与Simo和Vu-Quoc^[32]采用几何精确梁分析理论得到的结构非线性路径一致，证明了本文考虑节点平衡所建立的三维梁单元几何刚度矩阵同样适用于曲梁结构的非线性分析. 进一步通过表6中的计算数据可得，本文的三维梁单元相比于Yang和Kuo^[9]的简化理论三维梁单元，在增量步数相同的情况下，总迭代次数和计算花费时间更少，展现出更高的计算效率.

5　结论

收起

本文基于三维实体梁线性化的增量虚功方程，分析得到梁单元增量刚度方程中体现空间力矩转动节点保持平衡的诱发矩矩阵，并指出诱发矩矩阵虽然在单元层面上是非对称的形式，但在结构层面上考虑多个单元组装后将恢复对称性. 根据诱发矩矩阵，构造了含参数的三维梁单元的几何刚度矩阵，利用几何刚度矩阵的对称性和刚体合格性检测推导出了简洁的三维梁单元几何刚度矩阵的显示表达，同时简化得到二维梁单元几何刚度矩阵. 通过对典型案例进行线性屈曲分析和非线性分析，证明了本文所推导的几何刚度矩阵可以很好的用于梁类型结构的屈曲和后屈曲分析，并且具有较好的计算精度和较高的求解效率. 本文建立的推导梁单元几何刚度矩阵的方法注重非线性分析过程中单元的刚体行为，物理意义清晰，推导过程简单，避免了传统非线性单元模型中复杂的插值离散和推导，并且最终建立的几何刚度矩阵相较于传统几何刚度矩阵也更为简洁，是推导单元几何刚度的一种新思路，适用于增量小变形范围内的结构屈曲分析及非线性分析. 后续作者将根据类似思路尝试拓展建立起板、壳单元简洁形式的几何刚度显式表达形式，应用于更加复杂的结构非线性分析.

基金

收起

国家自然科学基金项目(52078082)
内江市基础研究与应用基础研究项目(NJJH202316)
内江师范学院科研资助项目(2023ZD03)
内江师范学院大学生创新创业训练计划资助项目(X2024077)

参考文献

收起

文献

收起

参考文献引证文献

排序方式：

[1]

刘小靖, 王记增, 周又和. 一种分析柔性梁超大挠度问题的小波算法[J]. 固体力学学报, 2013, 33(S1): 159-165.

(Liu

X J

, Wang

J Z

, Zhou

Y H

. A wavelet method for studying large deflection bending problems of flexible beams[J]. Chinese Journal of Solid Mechanics, 2013, 33(S1): 159-165. (in Chinese))

[2]

Ouyang

, Bai

, Liu

S W

, Chan

S L

. Refined three-dimensional pile element formulation for second-orderanalysis of pile-supported structures accounting for complex loading conditions[J]. Engineering Structures, 2024, 301: 117281.

[3]

Kassimali

. Large deformation analysis of elastic-plastic frames[J]. Journal of Structural Engineering, 1983, 109(8): 1869-1886.

[4]

Kassimali

. Large deformation analysis of elastic space frame[J]. Journal of Structural Engineering, 1991, 117(7): 2069-2087.

[5]

Yang

Y B

, Chen

A Q.

, He

. Research on nonlinear, postbuckling and elasto-plastic analyses of framed structures and curved beams[J]. Meccanica, 2021, 56: 1587-1612.

[6]

Barsoum

R S

, Gallagher

R H

. Finite element analysis of torsional-flexural stability problems[J]. International Journal for Numerical Methods in Engineering, 1970, 2: 335-352.

[7]

Meek

J L

, Tan

H S

. Geometrically nonlinear analysis of space frames by an incremental iterative technique[J]. Computer Methods in Applied Mechanics and Engineering, 1984, 47: 261-282.

[8]

Chan

S L

, Kitiporncai

. Geometric nonlinear analysis of asymmetric thin-walled beam-columns[J]. Engineering Structures, 1987, 9: 243-254.

[9]

Yang

Y B

, Kuo

S R

. Theory and Analysis of Nonlinear Framed Structures[M]. Singapore: Prentice Hall, 1994.

[10]

Augarde

C E

. Generation of Shape Functions for Straight Beam Elements[J]. Computers & Structures, 1998, 68(6): 555-560.

[11]

夏拥军, 陆念力. 梁杆结构稳定性分析的高精度Euler-Bernoulli梁单元[J]. 沈阳建筑大学学报, 2006, 22(3): 362-366.

(Xia

Y J

, Lu

N L

. A new Euler-Bernoulli beam element with high accuracy for the stability analysis of beam structures[J]. Journal of Shenyang Jianzhu University (Natural Science), 2006, 22(3): 362-366. (in Chinese))

[12]

夏拥军, 缪谦. 一种新型空间梁单元及其在梁杆结构稳定分析中的应用[J]. 工程力学, 2009, 26(4): 86-91.

(Xia

Y J

, Miao

. A spatial beam element and its application to buckling analysis for framed structures[J]. Engineering Mechanics, 2009, 26(4): 86-91. (in Chinese))

[13]

Chan

S L

, Zhou

Z H

. Pointwise equilibrating polynomial element for nonlinear analysis of frames[J]. Journal of Structural Engineering, 1994, 120(6): 1703-1717.

[14]

Zhou

Z H

, Chan

S L

. Self-equilibrating element for second-order analysis of semirigid jointed frames[J]. Journal of Engineering Mechanics, 1995, 212(8): 896-902.

[15]

聂国隽, 钱若军. 薄壁梁元几何非线性分析模型的研究[J]. 空间结构, 2002, 8(3): 35-40.

(Nie

G J

, Qian

R J

. Research on the geometric nonlinear analytical model of thin-walled beam element[J]. Spatial Structures, 2002, 8(3): 35-40. (in Chinese))

[16]

许红胜, 周绪红, 舒兴平. 空间刚框架几何非线性分析的一种新单元[J]. 工程力学, 2003, 20(4): 39-44.

(Xu

H S

, Zhou

X H

, Shu

X P

. A new element for geometric nonlinear analysis of three-dimensional steel frames[J]. Engineering Mechanics, 2003, 20(4): 39-44. (in Chinese))

[17]

郭鑫, 李东升, 魏达, 张雨. 基于解析位移形函数的改进风力机叶片梁单元[J]. 太阳能学报, 2022, 43(4): 387-392.

(Guo

, Li

D S

, Wei

, Zhang

. Improved beam element for wind turbine blades based on analytical displacement shape functions[J]. Acta Energiae Solaris Sinica, 2022, 43(4): 387-392. (in Chinese))

[18]

Zhao

C H

, Pan

W H

, Luo

Y Z

. Exact matrix stiffness method for out-of-plane buckling analysis of funicular arches considering warping deformations[J]. International Journal of Structural Stability and Dynamics, 2023, 23(1): 2350010.

[19]

邓继华, 梁璐祎, 谭平, 周福霖. 共旋梁单元几何非线性分析的子结构方法[J]. 土木工程学报, 2024, 57(05): 65-75.

(Deng

J H

, Liang

L Y

, Tan

, Zhou

F L

. Substructure method for geometric nonlinear analysis of corotational beam elements[J]. China Civil Engineering Journal, 2024, 57(05): 65-75. (in Chinese))

[20]

Liu

Y Z

, Yang

Y B

, Liu

X H

, Guo

D Z

, Xu

. A novel beam element with nine DOFs per node for resolving the warping-distortion compatibility in analysis of frames and curved beams made of I-sections[J]. Thin-Walled Structures, 2023, 194: 111314.

[21]

Liu

Y Z

, Yang

Y B

, Liu

X H

, Guo

D Z

. Lateral-distortional buckling of frames composed of non-aligned I-members by a simple distortional beam element considering angling effect[J]. Thin-Walled Structures, 2024, 202: 112146.

[22]

陈安全, 林航葳. 不考虑几何刚度的框架结构非线性分析方法[J/OL]. 计算力学学报, https://link.cnki.net/urlid/21.1373.O3.20240718.1722.018.

(Chen

A Q

, Lin

H W

. Nonlinear analysis method for frame structure with no consideration of geometric stiffness[J/OL]. Chinese Journal of Computational Mechanics, https://link.cnki.net/urlid/21.1373.O3.20240718.1722.018. (in Chinese))

[23]

Yang

Y B

, Chen

A Q

, He

, Wu

Y T

. Rational nonlinear analysis of framed structures and curved beams considering joint equilibrium in deformed state[J]. International Journal of Non-Linear Mechanics, 2020, 125: 103538.

[24]

McGuire

, Gallagher

R H

, Ziemian

R D

. Matrix Structural Analysis[M]. New York: John Wiley & Sons, Inc., 2000.

[25]

Yang

Y B

, Chiou

H T

. Rigid body motion test for nonlinear analysis with beam Elements[J]. Journal of Engineering Mechanics, 1987, 113(9): 1404-1419.

[26]

陈安全, 范秀娟. 几何非线性分析中的更新正交迭代策略[J/OL]. 计算力学学报. https://link.cnki.net/urlid/21.1373.O3.20250108.1356.034.

(Chen

A Q

, Fan

X J

. Updated orthogonal iteration strategy in geometric nonlinear analysis[J/OL]. Chinese Journal of Computational Mechanics, https://link.cnki.net/urlid/21.1373.O3.20250108.1356.034. (in Chinese))

[27]

Duy

, Pruettha

, Tinh

Q B

, Jaroon

. Locking-free isogeometric Timoshenko-Ehrenfest beam formulations for geometrically nonlinear analysis of planar beam structures[J]. Mechanics of Advanced Materials and Structures, 2022, 31(2): 356-377.

[28]

Wood

R D

, Zienkiewicz

O C

. Geometrically nonlinear finite element analysis of beams, frames, arches and axisymmetric shells[J]. Computers & Structures, 1977, 7: 725-735.

[29]

DaDeppo

D A

, Schmidt

. Instability of clamped-hinged circular arches subjected to a point load[J]. Journal of Applied Mechanics, 1975, 42(4): 894-896.

[30]

Argyris

J H

, Hilbert

, Malejannakis

G A

, Scharpf

D W

. On the geometrical stiffness of a beam in space—a consistent v.w. approach[J]. Computer Methods in Applied Mechanics and Engineering, 1979, 20: 105-131.

[31]

Romero

. A comparison of finite elements for nonlinear beams: the absolute nodal coordinate and geometrically exact formulations[J]. Multibody System Dynamics, 2008, 20(1): 51-68.

[32]

Simo

J C

, Vu-Quoc

. A three-dimensional finite-strain rod model. Part II: Computational aspects[J]. Computer Methods in Applied Mechanics and Engineering, 1986, 58(1): 79-116.

2025年第46卷第4期

PDF下载

引用本文

BibTeX

文章信息

doi: 10.19636/j.cnki.cjsm42-1250/o3.2025.017

接收时间：2025-06-27
首发时间：2026-03-20
出版时间：2025-08-27

补充材料

相关文章

文章信息

作者

出版历史

收稿日期：2025-06-27

基金

国家自然科学基金项目(52078082)

内江市基础研究与应用基础研究项目(NJJH202316)

内江师范学院科研资助项目(2023ZD03)

内江师范学院大学生创新创业训练计划资助项目(X2024077)

作者信息

¹内江师范学院建筑工程学院，内江，641100

²重庆大学土木工程学院，重庆，400045

通讯作者:

**E-mail：chen.anquan@foxmail.com.

参考文献

分享链接

https://castjournals.cast.org.cn/joweb/gtlxxb/CN/10.19636/j.cnki.cjsm42-1250/o3.2025.017

分享至

全文二维码

扫描看全文

引用本文

BibTeX

本文的引用情况

2种不同金属材料的力学参数

科 Family	属数 Number of genus	种数 Number of species	占总种数比例 Percentage of total species (%)	属 Genus	种数 Number of species	占总种数比例 Percentage of total species (%)
鹅膏菌科Amanitaceae	2	11	5.26	鹅膏菌属 Amanita	10	4.78
小菇科 Mycenaceae	2	12	5.74	丝盖伞属 Inocybe	5	2.39
多孔菌科 Polyporaceae	8	14	6.70	蜡蘑属 Laccaria	5	2.39
红菇科 Russulaceae	3	23	11.00	小皮伞属 Marasmius	6	2.87
				小菇属 Mycena	11	5.26
				光柄菇属 Pluteus	5	2.39
				红菇属 Russula	17	8.13
				栓菌属 Trametes	5	2.39

关闭全屏

BibTeX
EndNote
RefWorks
TxT

序号	支承情况	理论屈曲荷载（N）	本文	Yang和Kuo
（a）	两端铰接	246.45	247.21	0.31%	246.70	0.10%
（b）	一端固接，一端铰接	502.96	506.81	0.77%	504.72	0.35%
（c）	一端固接，一端自由	61.61	61.71	0.16%	61.68	0.11%
（d）	两端固接	985.80	994.95	0.93%	987.01	0.12%
（e）	两端铰接，中间铰支承	985.80	994.95	0.93%	987.01	0.12%

序号

支承情况

理论屈曲荷载（N）

本文

Yang和Kuo

曲荷载（N）特征值屈

误差

特征值屈曲荷载（N）

误差

（a）

两端铰接

246.45

247.21

0.31%

246.70

0.10%

（b）

一端固接，一端铰接

502.96

506.81

0.77%

504.72

0.35%

（c）

一端固接，一端自由

61.61

61.71

0.16%

61.68

0.11%

（d）

两端固接

985.80

994.95

0.93%

987.01

0.12%

（e）

两端铰接，中间铰支承

985.80

994.95

0.93%

987.01

0.12%

方法	增量步数	总迭代次数	计算时间（s）
本文	8545	13593	4.39
Yang和Kuo	8545	13612	4.97

方法

增量步数

总迭代次数

计算时间（s）

本文

8545

13593

4.39

Yang和Kuo

8545

13612

4.97

方法	增量步数	总迭代次数	计算时间（s）
本文	4282	8792	2.27
Yang和Kuo	4287	8803	2.73

方法

增量步数

总迭代次数

计算时间（s）

本文

4282

8792

2.27

Yang和Kuo

4287

8803

2.73

方法	增量步数	总迭代次数	计算时间（s）
本文	21329	43020	10.62
Yang和Kuo	21445	43327	11.73

方法

增量步数

总迭代次数

计算时间（s）

本文

21329

43020

10.62

Yang和Kuo

21445

43327

11.73

方法	增量步数	总迭代次数	计算时间（s）
本文	516	1035	1.22
Yang和Kuo	518	1049	1.36

方法

增量步数

总迭代次数

计算时间（s）

本文

516

1035

1.22

Yang和Kuo

518

1049

1.36

方法	增量步数	总迭代次数	计算时间（s）
本文	389	980	0.34
Yang和Kuo	389	1949	0.55

方法

增量步数

总迭代次数

计算时间（s）

本文

389

980

0.34

Yang和Kuo

389

1949

0.55