三维欧拉-伯努利梁与铁木辛柯-埃伦费斯特梁单元的显式方法：质量凝聚、临界时间步估计及流固耦合扩展

《European Journal of Mechanics - A/Solids》：An explicit approach for 3-D Euler–Bernoulli and Timoshenko-Ehrenfest beam elements: Mass lumping, critical time step estimates and extension to FSI

　　

•提出了适用于欧拉-伯努利梁单元和铁木辛柯-埃伦费斯特梁单元的统一方法。
•重点在于显式时间积分和一致性形函数。
•给出了根据转动质量矩阵凝聚方式的不同而变化的弯曲时间步估计。
•研究了剪力强系数对铁木辛柯-埃伦费斯特梁稳定性区域的影响。
•提出了扩展至充液管道系统中流体-结构相互作用的方法。

**引言**

梁理论已经发展完善并被广泛应用于许多工业领域（如航空航天或土木工程）的工程结构分析中，参见Reddy（2022）。两种主要且常用的梁理论是剪切约束的欧拉-伯努利（EB）理论和一阶剪切变形的铁木辛柯-埃伦费斯特（TE）理论。在EB理论中，基本假设是截面保持平面且垂直于梁轴，忽略剪切效应。Timoshenko（1921）对EB理论进行了改进，同时考虑了剪切变形和转动惯性，由此产生了铁木辛柯-埃伦费斯特（TE）梁理论（Elishakoff, 2019）。此外，值得注意的是，也提出了若干高阶剪切变形梁理论，例如Reddy提出的三阶剪切变形Reddy梁模型（Reddy, 1984；Heyliger和Reddy, 1988），综述参见Reddy（2022）。Ruocco和Reddy（2023）最近对EB、TE和Reddy梁理论进行了比较，展示了它们各自的行为和局限性。虽然TE模型代表了EB模型的有效改进，但它受恒定切向应变假设的限制，需要剪切修正（Cowper, 1966；Hutchinson, 2001）。此外，从数值角度来看，TE梁单元存在剪切锁定问题（Mukherjee和Prathap, 2001）。这影响了这些单元的精度和效率。已提出不同的方法以避免剪切锁定。一种常见方案是考虑缩减积分（Yokoyama, 1994），这在某些情况下可以缓解剪切锁定，但仍然影响单元的性能。事实上，剪切锁定是由于相互依赖的自由度（DOFs）之间插值的非一致性导致的。因此，在用于插值未知量的形函数构造中，必须考虑这种相互依赖性。为此，Reddy（1997）、Luo（1998）和Luo（2008）通过使用与静力平衡梁方程一致的形状函数，开发了不受剪切锁定影响的效率型TE梁单元。需要注意的是，Reddy（1997）和Luo（1998）中仅考虑了二维TE单元，该方法的扩展至三维TE梁单元是在Luo（2008）中完成的。在上述论文中，推导梁平衡方程、形状函数和刚度矩阵时仅考虑了静态条件。将动态运动方程与显式时间积分相结合的扩展具有重要价值，特别是对于柔性管道中的水锤现象和管道鞭甩等快速瞬态事件。

显式时间积分方法以条件稳定著称，计算时间步必须满足限制其大小的稳定性条件。此外，在结构动力学中使用显式时间积分时，通常伴随对角质量矩阵的使用以便于其求逆（Hughes, 2000；Wu, 2006；Wu和Gu, 2012；Belytschko等, 2014）。可以采用多种凝聚方案来获得这个对角质量矩阵。例如，经典的方法包括行求和法或对角缩放法（Key和Beisinger, 1971；Hinton等, 1976）。对角质量缩放已在理论上进行了分析，常见的缩放方案综述参见Voet等（2025）。在壳和梁单元的特殊情况下，需要考虑具有不同物理意义且单位不同的自由度（DOFs），即位移和转角。因此，必须使用适当的凝聚方案来处理每种类型的自由度（即平动和转动），以避免不一致性。显然，凝聚方案直接影响梁单元的精度（即动态响应、振动、频率方面）和稳定性（即临界时间步估计方面）。对于平动自由度，质量矩阵的常见凝聚方法基于质量守恒。相比之下，关于转动自由度尚无真正的共识。自Key和Beisinger（1971）、Krieg和Key（1973）、Hughes等（1978）以及Belytschko和Mindle（1980）的早期工作以来，一直有人在研究对应转动惯量的对角项缩放。例如，已证明增加转动惯量的缩放系数对临界时间步有显著好处（Key和Beisinger, 1971；Belytschko和Mindle, 1980）。Archer和Whalen（2005）提出了针对欧拉-伯努利和铁木辛柯-埃伦费斯特梁单元转动自由度的凝聚方案，使得对角转动质量项的总转动惯量与一致性矩阵相同。Zuo等（2014）提出的欧拉-伯努利梁单元凝聚方案基于最小化模态误差。Archer和Whalen（2005）及Zuo等（2014）所提出的对角质量矩阵的主要缺点是引入了转动自由度的负对角项，这在物理上是可疑的。

本文的目的是考虑EB梁单元，并将Luo（2008）中提出的不受剪切锁定影响的三维TE梁单元扩展至瞬态动力学。为此，分析了这些梁单元的显式时间积分。特别是，研究了转动自由度的凝聚对临界时间步的影响。此外，在铁木辛柯-埃伦费斯特梁单元的情况下，检验了剪力强系数对临界时间步的影响。

本文组织如下。在第2节中，使用虚功原理（VWP）推导了欧拉-伯努利和铁木辛柯-埃伦费斯特梁运动学的运动方程。第3节专门讨论有限元方法框架下的显式时间积分。特别是给出了考虑横向位移和弯曲转角之间相互依赖性的形函数。然后，表达了两种梁单元对应的质量矩阵和刚度矩阵。接着，重点讨论了质量矩阵凝聚、梁运动学以及剪力强系数对临界时间步（尤其是弯曲时间步）的影响。第4节分析了与本文理论结果一致的数值算例，既考虑了静态和动态算例，并将结果与参考解进行了比较。最后，第5节提出了扩展至流体-结构相互作用（FSI）的方法。在充液管道的背景下，将所描述的梁单元与使用面积平均欧拉方程建模的内部流体耦合，同时考虑管道的弹性和柔韧性。

**欧拉-伯努利和铁木辛柯-埃伦费斯特梁模型**

考虑一段长度为?的直管道，如图1所示。与未变形梁相关联的局部参考坐标系记为$E_x$、$E_y$、$E_z$，其中$E_x$对应管道轴。考虑未变形管道截面上的一点m，其在$E_x$、$E_y$、$E_z$坐标系中的坐标为x、y和z，以及同一截面上梁轴上的一点M，即坐标为$(x, 0, 0)$。从点M到截面上点m的相对位置向量为：

**一致性形函数**

为了数值逼近前述梁运动方程，考虑具有每节点六个自由度（三个平动和三个转动，对应位移场$d_x$、$d_y$、$d_z$、$\theta_x$、$\theta_y$和$\theta_z$）的三维二节点单元。在有限元方法框架下，使用形函数在计算域的每个单元上插值自由度（DOFs）。如在Daude等（2025b）和Luo（2008）中对铁木辛柯-埃伦费斯特梁所做的...

**数值结果**

为了验证所呈现梁模型的各种方面，考虑了静态和动态算例。一方面，静态结果允许将数值模型与解析解进行验证，并评估铁木辛柯-埃伦费斯特模型中考虑剪切变形的影响。另一方面，动态算例侧重于数值方案的稳定性以及质量凝聚对精度的影响...

**控制方程**

欧拉-伯努利和铁木辛柯-埃伦费斯特梁理论已被广泛应用于充液管道系统中的流体-结构相互作用（FSI），综述参见Tijsseling（1996）、Wiggert和Tijsseling（2001）、Li等（2015）以及Ferras等（2018）。对于充液管道，内部流体对横向位移起到附加质量的作用，从而得到以下直欧拉-伯努利充液梁的运动方程：

$$\rho_s A_s \frac{\partial^2 d_x}{\partial t^2} - \frac{\partial F_x}{\partial x} = -\pi d \tau_w + \rho_s A_s g \cdot E_x$$

$$\rho_s A_s + \rho_f A_f \frac{\partial^2 d_y}{\partial t^2} - \rho_s I_z \frac{\partial^2 \theta_x}{\partial x^2}...$$

**结论与展望**

本文在有限元方法和显式时间积分的背景下，为欧拉-伯努利和铁木辛柯-埃伦费斯特梁提出了一种通用的数值方法。使用了与梁平衡方程一致并因此考虑了自由度间相互依赖性的形函数。这些一致性形函数的主要效果是消除了铁木辛柯-埃伦费斯特梁单元的剪切锁定效应。相应的质量矩阵和...

**CRediT作者贡献声明**

F. Daude：撰写——初稿、验证、监督、软件、方法论、调查、形式分析、概念构思。
M. Sharifi：撰写——初稿、验证、监督、软件、方法论、调查、形式分析、概念构思。
P. Galon：撰写——审阅与编辑、监督、方法论、调查、形式分析、概念构思。
S. Potapov：撰写——审阅与编辑、方法论、调查、形式分析、概念构思。

**利益冲突声明**

作者声明，他们没有已知的竞争性经济利益或个人关系可能影响本文所述工作。

**致谢**

本工作是在EDF/CEA/Framatome三方研究所的"快速动力学"项目框架内完成的。计算设施由EDF提供。数值模拟使用Europlexus软件（http://www-epx.cea.fr/）进行。作者感谢审稿人的意见，这些意见有助于提高文章的质量。

F. Daude | M. Sharifi | P. Galon | S. Potapov