复合单元法在变截面梁自由振动分析中的应用*

吕中荣, 谢瑾荣,刘济科

(中山大学工学院，广东 广州 510275)

在工程结构中，有很多可以简化成梁结构的形式，或者本身就是梁。例如，土木工程中的跨度桥梁结构；机械工程中的转子轴承；航空工程中的天线等。因此，对梁的动力特性分析的研究在梁的设计和损伤检测中具有重要意义[1-2]。对梁的自由振动问题分析已有的方法很多，有理论方法，如Rayleigh-Ritz法[3], Component Modal Analysis 法[4]，Green 函数法等[5]；还有数值方法，如有限元法[6-8]。众所周知，利用有限元法要想得到较精准的高阶固有频率和振型，需划分较多的单元，因而导致计算量较大。

本文利用复合单元法建立了多阶梯变截面和渐变截面梁的有限元分析方程[9-11]，并进行了自由振动问题的分析。复合单元法的主要优点就在于其不但继承了有限单元法将结构离散化的思想，又引入经典的振动分析中的特征函数，二者相耦合。这样在保证一定计算精度的同时，能使划分的单元数量明显减少，从而可以减少计算量，提高计算效率。

1 复合单元法简介

复合单元法是一个相对较新的建模工具。这个方法将常规的有限单元法与高精度的解析理论相结合。用复合单元法解梁的问题时，位移函数是有限元法的位移函数和振动理论中的特征函数的耦合和扩大[9]。

1.1 复合单元法的位移函数

复合单元法的位移函数的一般表示为

uCEM(x,t)=uFEM(x,t)+uCT(x,t)

(1)

其中uFEM(x,t)为有限单元法(FEM 是finite element method 的缩写)的位移场

uFEM(x,t)=N(x)q(t)

(2)

这里，N(x)为有限单元法梁单元的形函数

(3)

q(t)为节点位移矢量

q(t)=[v1(t),θ1(t),v2(t),θ2(t)]T

(4)

上式中‘v’ 和 ‘θ’分别代表竖向位移和转角。

uCT(x,t)是解析法(CT 是classical theory 的缩写)位移函数

(5)

φi(i=1,2,…N)是解析法梁的特征函数，ci(t)是c-自由度。

由式(2)-(5), 式(1)可写为如下形式

uCEM(x,t)=S(x)Q(t)

(6)

其中S(x)为耦合后得到的复合单元法形函数

S(x)=[N1(x),N2(x),N3(x),N4(x),

φ1(x),φ2(x),...,φN(x)]

(7)

Q(t)为耦合后单元节点自由度

Q(t)=[v1(t),θ1(t),v2(t),θ2(t),

c1(t),c2(t),...,cN(t)]T

(8)

2 自由振动分析

2.1 整体刚度矩阵

假设矩形阶梯梁共有i个阶梯单元，假设每个单元长度为Li，截面高度hi分布为：

(9)

抗弯刚度矩阵可表示为：

(10)

复合单元法单元形函数组装成整体形函数的形式如下

(11)

其中，i代表阶梯数，j代表复合自由度级数。组装时，各个单元的有限单元法形函数N(x)须上下对应，解析法的特征函数φi不同单元的相同级数之间须上下对应。

在悬臂梁的自由振动问题中，特征函数φi为

(12)

其中，λi满足cos(λi)*cosh(λi)=1，且li-1

整体刚度矩阵的表达式为

(13)

其中：[kqq]对应FEM中的梁的元素刚度矩阵，[kqc]对应q自由度和c自由度的耦合，[kcq]是[kqc]的转置，[kcc]对应c自由度，是一个对角矩阵。

2.2 整体质量矩阵

各阶梯段的质量分布为

(14)

类似地，整体质量矩阵可表示为

(15)

其中，[mqq]对应FEM中关于梁的质量矩阵，[mqc]对应q自由度和c自由度的耦合，[mcq]是[mqc]的转置，[mcc]对应c的自由度。

得到整体刚度矩阵K和质量矩阵M后，梁无阻尼自由振动的特征方程表示如下

(K-ω2M)V=0

(16)

由此可确定系统的固有圆频率ω和相应的振型V。

复合单元法的振型函数可表示为

Φ=S(X)V

(17)

3 数值算例

算例1：运用复合单元法与文献[12]中的应用实例分析结果比较, 其中的悬臂梁结构和尺寸如图1所示。梁的物理参数为：密度ρ=2 664 kg/m3, 弹性模量E=60.6 GPa。

图1 多变截面悬臂阶梯梁(单位：mm)

采用本文方法计算的频率与文献[12]的比较见表1所列。由表可以看出，采用本文方法208个自由度的复合单元计算结果与文献[12]中500个自由度的理论值精度相当。表明本文方法具有自由度少，精度高的优点。

算例2：进一步采用本文方法对一渐变截面梁的自由振动进行分析。计算的悬臂梁结构如图2所示。梁的高度沿着长度方向的变化假定为：

(18)

图2 渐变截面悬臂梁模型

梁的物理参数：密度ρ=7 800 kg/m3, 弹性模量E=210 GPa, 梁的宽度b=20 mm, 高度h0=20 mm，总长L=1 200 mm。采用4个有限单元加30个c-DOF 的复合单元进行计算。并与ANSYS解进行比较。前六阶频率的比较如表2所列，由此可见，本文方法和ANSYS解符合得非常好，表明本文方法的正确性和有效性。

表1 固有频率与文献[12]的比较1)

1)ω1B,ω2B分别代表第一阶和第二阶面外弯曲模态频率；ω1c代表第一阶面内弯曲模态频率；

CEM(4×200c)表示采用4个有限单元，200个c-DOF，其他类推。

表2 固有频率与ANSYS的解比较

4 结 论

利用复合单元法对变截面梁的自由振动进行了分析，其主要优点就在于不但继承了有限单元法将结构离散化得思想，又引入级数形式的解析法，二者相互耦合。这样在保证一定精度的同时，能使划分的单元数量明显减少，从而大大减少计算量。本文的两个算例表明，通过复合单元法进行自由振动分析时，将大大减少计算量，同时保证一定的计算精度。

参考文献：

[1] SINHA J K, FRISWELL M I,EDWARDS S. Simplified models for the location of cracks in beam structures using measured vibration data [J]. Journal of Sound and Vibration, 2002, 251(1): 13-38.

[2] LU Z R, LIU J K, HUANG M,et al. Identification of local damages in coupled beam systems from measured dynamic responses [J]. Journal of Sound and Vibration, 2009, 326(1/2): 177-189.

[3] BHAT R. B. Nature of stationarity of the natural frequencies at the natural modes in the Rayleigh-Ritz method [J].Journal of Sound and Vibration, 1997, 203(2):251-263.

[4] 李惠彬. 振动理论与工程应用 [M].北京：北京理工大学出版社, 2006.

[5] KUKLA S,ZAMOJSKA I. Frequency analysis of axially loaded stepped beams by Green’s function method [J].JournalofSoundandVibration, 2006, 300:1034-1041.

[6] 李景涌.有限元法[M].北京：北京邮电大学出版社.2002.

[7] 吕中荣, 罗绍湘, 刘济科. 预应力对预应力梁响应的影响[J]. 中山大学学报：自然科学版, 2006, 45(2): 119-120.

[8] 徐伟华, 吕中荣, 刘济科. 弱耦合系统模态局部化和频率曲线转弯现象[J]. 中山大学学报：自然科学版,2007, 46(4): 9-12, 20.

[9] 曾攀. 计算力学中的高精度数值分析新方法-复合单元法[J]. 中国科学：E辑, 2002, 30(1)：39-46.

[10] LU Z R，LAW S S. Discussions on “Composite Element Method for Vibration Analysis of Structure” [J].Journal of Sound and Vibration, 2007, 305(1/2):357-361.

[11] LU Z R，LAW S S. Dynamic Condition Assessment of a Cracked Beam with the Composite Element Model [J].Mechanical Systems and Signal Processing, 2009, 23(2): 415-431.

[12] JAWORSKI J W，DOWELL E H. Free vibration of a cantilevered beam with multiple steps: Comparison of several theoretical methods with experiment [J]. Journal of Sound and Vibration, 2008, 312(4/5): 713-725.