第4章 多自由度系统振动分析的数值计算方法(25页)

第4章 多自由度系统振动分析的数值计算方法

用振型叠加法确定多自由度系统的振动响应时，必须先求得系统的固有频率和主振型。当振动系统的自由度数较大时，这种由代数方程求解系统固有特性的计算工作量很大，必须利用计算机来完成。在工程中，经常采用一些简单的近似方法计算系统的固有频率及主振型，或将自由度数较大的复杂结构振动问题简化为较少阶数的振动问题求解，以得到实际振动问题的近似分析结果。

本章将介绍工程上常用的几种近似解法，适当地选用、掌握这类实用方法，无论对设计研究或一般工程应用都将是十分有益的。

§4.1 瑞利能量法

瑞利（Rayleigh ）能量法又称瑞利法，是估算多自由系统振动基频的一种近似方法。该方法的特点是：①需要假定一个比较合理的主振型；②基频的估算结果总是大于实际值。由于要假设主振型，因此，该方法的精度取决于所假设振型的精度。

§4.1.1 第一瑞利商

设一个n 自由度振动系统，其质量矩阵为[]M 、刚度矩阵为[]K 。多自由度系统的动能和势能一般表达式为

{}[]{}{}[]{}/2/2T

T

T x M x U x K x ?=?

?=??

(4.1.1)

当系统作某一阶主振动时，设其解为

{}{}(){}{}()sin cos x A t x A t ωαωωα=+??

?=+??

(4.1.2)

将上式代入式（4.1.1），则系统在作主振动时其动能最大值max T 和势能最大值max U 分别为

{}[]{}{}[]{}2max max /2/2T

T

T A M A U A K A ω?=?

?=??

(4.1.3)

根据机械能守恒定律，max max T U =，即可求得

{}[]{}{}[]{}

()2I T

T

A K A R A A M A ω== (4.1.4)

其中，()I R A 称为第一瑞利商。当假设的位移幅值列向量{}A 取为系统的各阶主振型

{}i A 时，第一瑞利商就给出各阶固有频率i ω的平方值，即

{}[]{}{}[]{}

2(1,2,,)T

i i i T

i i A K A i n A M A ω==

(4.1.5)

在应用上式时，我们并不知道系统的各阶主振型{}i A ，只能以假设的振型{}A 代入式（4.1.4），从而求出的相应固有频率i ω的估计值。从理论上讲，可用式（4.1.4）近似求解各阶固有频率，但由于对系统的高阶主振型很难作出合理的假设，所以，该式一般只用来估算系统的基频1ω。

§4.1.2 第二瑞利商

瑞利能量法也可以应用于由柔度矩阵Δ[]δ建立的位移运动方程。这时自由振动方程

{}[][]{}x M x δ=-

(4.1.6)

代入式（4.1.1），注意到[]δ、[]M M 是对称矩阵，以及[][][]K I δ=，则系统的势能为

{}[][][]{}2T T U x M M x δ=

(4.1.7)

由式（4.1.2）可得

2{}{}sin()x A t ωωα=-+

(4.1.8)

将上式代入式（4.1.7），系统势能的最大值为

4max {}[][][]{}2T T U A M M A ωδ=

(4.1.9)

由max max T U =可得

{}[]{}{}[][][]{}

2II ()T

T

A M A R A A M M A ωδ=

=

(4.1.10)

II ()R A 称为第二瑞利商。

可以证明，若所选假设振型{}A 很接近于第一阶主振型{}1A ，则由第一瑞利商和第二瑞利商计算出的2

ω值确实接近于21ω，而且比实际值稍大（所谓上限估计）。对于同一假设振型{}A ，第二瑞利商比第一瑞利商更接近真实值21ω，但其精确程度主要取决于假设振型

{}A 接近于第一阶主振型{}1A 的程度。

在图4.1.1所示三自由度系统中，试用瑞利能量法估算系统的第一阶固有频率。已知123m m

m m ===，123k k k k ===。

图 4.1.1

【解】 系统的质量矩阵为

[]001000001000001m M m m m ????

????==????????????

刚度矩阵为

[]12

22

23

33

302101210

011k k k K k k k k k k k +--????

????=-+-=--????

????--????

柔度矩阵为

[][]11111122123K k δ-??

??==??????

粗略地假设振型为{}[111]T A =，从而得

{}[]{}3T A M A m = （1） {}[]{}T A K A k = （2）

{}[][][]{}

21001111001114(111)01012201010011230011T

A M M A m m m k k

δ????????????????==????????????????????????

（3） 式（1）、（2）代入式（4.1.4）得

{}[]{}{}[]{}

20.3333T

T

A K A k k m m A M A ω=== 式（1）、（3）代入式（4.1.10）得

{}[]{}{}[][][]{}

230.21414T

T

A K A k k

m m A M M A ωδ=== 系统的第一阶固有频率的精确值为2

10.198

k

m

ω=，显然第二瑞利商的结果较接近精确值，但误差还较大，这是因为假设振型与第一阶精确振型T

1{}[0.4450.8021]A =相差较远的缘故。

如果在图4.1.1的每一个质量上顺坐标方向分别作用一单位力，则以该静变形曲线作为假设振型，即取

{}[]3

56T

A =

则有

{}[]{}70T

A M A m = {}[]{}14T

A K A k =

{}[][][]{}2

353T

m A M M A k

δ=

由式（4.1.4）得

20.200

k m

ω= 由式（4.1.10）得

20.198

k m

ω= 可见，假设振型与第一阶主振型愈接近，则瑞利商结果愈接近于基频21ω。

如图4.1.2所示，已知梁的弯曲刚度为

EJ ，不计其质量，13m m m ==，

22m m =，求系统的第一阶固有频率。

图4.1.2

【解】 系统的质量矩阵为

[]0010002002000001m M m m m ????

????==????????????

柔度矩阵为

[]391171116117687119l EJ δ??

??=??????

粗略地选取假设振型为 {}[12

1]

T

A =，则 {}[]{}10T A M A m =

{}[][][]{}[]323

100911712912102011161127684800171191T l m l A M M A m EJ

EJ δ??

??????

????==??

????????????

????

代入式（4.1.10）得

{}[]{}{}[][][]{}2233

1016.55

29/48T

T

A M A m

EJ

m l EJ

ml A M M A ωδ===

ω=

系统第一阶固有频率的精确值为1 4.024ω=1%。在系统柔

度矩阵已知的情形下，若假设振型用{}[][]{}A M I δ=，则计算精度还可提高。

§4.2 邓克莱法

邓克莱(Dunkerley)法又称迹法。前述的瑞利能量法给出了系统最低阶固有频率的上限估计值，而邓克莱法则给出了系统最低阶固有频率的下限估计值。

如前所述，n 自由度系统的位移方程：

{}[][]{}0x M x δ+=

(4.2.1)

设其解为 {}{}()sin x A t ωα=+

代入式（4.2.1），并以2

ω除全式得主振型方程

[][][](){}{}2

0I M A δ-=

(4.2.2)

其特征方程为

21112

11112

1n 2

2122

22122

2n 21

2n1n2

nn 10001000

00

1n n n n nn m m m m m m m m m δδδωδδδωδδδ???????????

???????=?????

??

?????????

??－

当系统的质量矩阵为对角矩阵时，可展开为

()()21211122210n n nn n m m m ωδδδω--++++

=

由代数方程理论（多项式根与系数之间的关系）可知，上式中()211n -项的系数变号后

等于21

ω的n 个根22

2

12

n

111ωω、、之和，即

222

12n 111222111nn n m m m ωωωδδδ=+++＋＋＋

(4.2.3)

对等式（4.2.3）作如下处理： 等式左边，由于12n ωωω<<

<，

即1

2

n

1

1

1

ωωω>

>>

，故近似地只保留一项211

。

等式右边，令

[][][]D M δ=

(4.2.4)

[]D 称为动力矩阵（dynamic matrix ），则式（4.2.3）右边为动力矩阵的迹，记为[]tr D 。

因为()1,2,

,ii i n δ=是第i 个质量处作用单位力时系统在该处的柔度系数。设想系统只有

一个质量i m 存在，则系统成为单自由度系统，这时系统的刚度1i ii k =，固有频率i Ω为

21i i i ii i k m m δΩ==， 即21

ii i i

m δ=

Ω，于是有

[][][]()2

2

1

1

ii i i i m tr D tr M δδω≈===Ω∑∑

(4.2.5)

综上所述，式（4.2.3）可写为

2

222112

1

111

n

ω≈

+++

ΩΩΩ (4.2.6)

即系统的最低阶固有频率平方值的倒数，近似等于各质量i m 单独存在时固有频率平方值2i Ω的倒数之和。由于式（4.2.3）的左边舍去了一些正数值，从而所得的21ω值比真值小。式（4.2.6）称为邓克莱公式，计算出的结果为最低阶固有频率的下限估值。由于等式右边为动力矩阵

[]D 的迹[]tr D ，故邓克莱法又称为迹法，它只适用于[]M 为对角矩阵的系统。

邓克莱法在准确度上一般不如瑞利能量法，但由于它的计算较简单，且易考虑各质量或刚度的变化对最低阶固有频率的影响，故工程上仍经常应用它。

用邓克莱法计算例4.1中系统的基频。

【解】 由例4.1可知，系统的质量矩阵和柔度矩阵分别为

[]100010001M m ????=??????， []1

1111221

2

3k

δ??

??=??????

动力矩阵[]D 为

[][][]1111001111122010122123001123m D M m k

k δ??????

??????===??????????????????

其迹为

[][][]()6m

tr D tr M k

δ==

由式（

4.2.5）得系统的基频为

210.167

k m ω≈， 10.4ω= 上述结果与精确值相比误差较大，大约为8.08％。

已知一均质悬臂梁如图4.2.1所示，式中EJ 为抗弯刚度，M 为梁的总质量，l 为

梁长，其第一阶固有频率的平方2

231 3.515/EJ Ml Ω=。若在梁的自由端放置一激振器质量

为m ，设激振器质量与梁的质量之比/1/20,1/10,1/2,1m M =，试用邓克莱法估算系统的基频值，并说明激振器质量对均质梁固有频率的影响。

图 4.2.1

【解】 已知悬臂梁的固有频率的平方为

2231 3.515/EJ Ml Ω= （1）

由材料力学可知，其端点的柔度系数为3

/3l EJ δ=，激振器固有频率的平方为

23

2

3/EJ ml Ω= （2） 将式（1）、（2）代入式（4.2.6），得系统基频的平方为

()22

2121

22222

312121137.0661/ 3.515/3EJ m M Ml ω??ΩΩ??

≈+== ?

?ΩΩΩ+Ω+

??

?? （3） 由上式可知，系统的固有频率1ω与质量比/m M 值有关。将式（3）式改写为

=

（4） 对于不同的质量比，式（4）的值如表4-1所示：

表4.2.1 不同质量比的1/

ω

表4.2.1中的误差是与1Ω比较而言，可见，只有当激振器的质量为梁的质量的1/20以下时，激振器质量对梁的固有频率影响才可接受。

§4.3 李兹法

前述两种方法只限于估算振动系统的基频，但工程实际中往往需要求出前几阶的固有频率及相应主振型，应用瑞利能量法的困难在于较高阶固有频率的假设振型难于选择。李兹法在瑞利法的基础上较好地克服了上述困难，可计算系统的前几阶固有频率及主振型。

李兹（Ritz ）法不需要直接给出假设振型，而是将假设振型{}A 表示为有限（低维）个独立的假设模态{}

j φ的线性组合

{}{}[]{}1

s

j j j A c C φφ===∑

(4.3.1)

其中 []{}{}{}12s φφφφ=????， {}[]12

T

s C c c c =

{}j

φ为1n ?列阵，可预先选定，s n <，(1,2,

,)j

c j s =为待定常数。

将式（4.3.1）代入式（4.1.4），第一瑞利商为

{}[][][]{}{}[][][]{}

2I T

T T T

C K C R C M C φφωφφ==

(4.3.2)

显然，由上式还求不出固有频率，但与j C 有关。由于瑞利法是固有频率的上限估计，故j C 的选择应当使上式给出的固有频率值最小，即上式对j C 的偏导数应等于零。 令 {}[]

[][]{}()T

T

K c C K C φφ=，{}[][][]{}()T

T M c C M C φφ=

于是由

I 0j R C ?=?可得 2()()

0(1,2,,

)j j

K c M c j s c c ω??-==?? (4.3.3) 而

{}[][][]{}()

{}[][][]{}{}[][][]{}{}[][][]{}[][]{}()22T

T

T T T

T

j j

j j

T

T

T j j C C K c C K C K C C K c c c c C K C K C c φφφφφφφφφφ????==+?????????==

??? ?

???

同理

[][]{}()

2T j j

M c M C c φφ???=??? 于是，式（4.3.3）可写为

[][]{}[][]{}20

(1,2,,)T

T

j j K C M C j s φφωφφ????-==??

??

这s 个方程可合并为一个矩阵方程

[][][]{}[][][]{}{}20T T

K C M C φφωφφ-=

(4.3.4)

[][][][][][]**

T T K K M M φφφφ???=??

????=????

(4.3.5)

上式中，*K ????和*

M ????为s s ?阶对称阵，分别称为广义刚度矩阵（generalized stiffness

matrix ）和广义质量矩阵（generalized mass matrix ）。这样，式（4.3.4）可改写为

()

{}{}2

0K M C ω

*

*

????-=????

(4.3.6)

这样，问题又归结为特征值问题，所不同的是，现在为s s ?阶矩阵的特征值问题，而

不是原系统n n ?阶矩阵的特征值问题。因而，李兹法是一种缩减自由度的近似解法。

由上式求得的s 个特征值22

212s ωωω、、

就是原系统前s 阶固有频率平方的近似值。将

解得的s 个特征矢量{}{}{}12s C C C 、、进行归一化，代入式（4.3.1）可求得原系统前s 阶主振型的近似值，即

{}[]{}

(1,2,,)j

j

A C j s φ==

(4.3.7)

由式（4.3.6）可知，各{}

j C 对矩阵[]M 是正交的，即有

{}[]{}0

()T

i j C M C i j =≠

所以

{}[]{}{}[][][]{}{}

{}0

()T

T

T

T

i j i j i j A M A C M C C M C i j φφ*

??===≠??

上式说明用李兹法求得的s 阶近似主振型对质量矩阵[]M 也是正交的，同时它们对刚

度矩阵[]K 也是正交的，因此，对它们可以用振型叠加法分析系统的各种响应。

同理，如果将式（4.3.1）代入第二瑞利商式（4.1.10），也可归结为减维特征值问题： ()

{}{}*

2*

0M

L C ω????-=????

(4.3.8)

这里

[][][][][]T

L M M φδφ*??=??

(4.3.9)

应当指出，由李兹法求得的s 个ω值中，前

2s 个或1

2

s +个ω值比较接近于真值，而后面的ω值误差比较大。因此，若想求前s 个固有频率及主振型的近似解，缩减的自由度数目

最好不小于2s 个，这样就能得到较精确的解。

图4.3.1(a)所示为一等直杆，杆长为l ，截面面积为A ，密度为γ, 试用聚缩质量的方法将其离散为有限自由度系统，并用李兹法求杆纵向振动时第一阶固有频率和主振型的近似解。

123

45

( a )

( b )

图4.3.1

【解】 将直杆等分为五段，每段的质量/5m Al g γ=等分为两半，各集中于每段的两端，然后将五段合并聚缩为5个质量1234m m m m m ====，52m m =，各聚缩质量之间由刚度为k 的5个弹簧相连接，如图4.3.1(b)所示。每段杆的拉压刚度确定为5/k EA l =。 这样，我们就得到五自由度的离散系统。

系统的质量矩阵和刚度矩阵分别为

[]100000100

00

010*******

000

0.5M m ?????

?

??=????????

， []210

0012100012100012100011K k -??

??--?

???=--?

?--????-??