Euler-Bernoulli梁作为工程结构中基本组成单元被广泛地应用在航空航天、基础建筑设施、工程机械等领域，例如飞机机翼中的桁架结构、卫星太阳能帆板、火箭发射台支架、大跨度桥梁、智能机器人的机械手臂。所以Euler-Bernoulli梁的振动分析是《机械振动基础》课程中的重要内容，精确的模态表达式是计算分析振动梁动态响应的保证。而在相关教材[1–4]给出Euler-Bernoulli梁的模态表达式中，除了简支梁，一般只能计算前12阶模态。有学者称之为12阶瓶颈问题(the 12th-mode bottleneck)^[5]。Tang^[6]通过提高浮点计算精度，从而得到前100阶振动梁模态函数的数值解。Goncalves等^[7]给出一种近似方法，可以计算振动梁无限阶模态，但该方法会对低阶模态形状带来误差。近来Goncalves等^[8]进一步给出Euler-Bernoulli梁在不同经典边界条件时的模态表达式，可用于高阶模态的数值计算，但是其模态表达式过于复杂，不便于课堂教学。

本文对Euler-Bernoulli梁模态的经典表达式进行了修正，可以精确计算前200阶模态。为便于读者参考借鉴，列表给出经典边界条件（如两端固支、固支-简支、固支-自由、自由-自由）下的修正模态表达式。

1.   计算Euler-Bernoulli梁高阶模态时存在的问题

由振动理论可知，等截面均质Euler-Bernoulli梁计算模态的常微分方程为

式中，EI和ρ分别为梁的截面抗弯刚度和密度，A为梁的横截面积，X(x)和ω分别为梁的模态表达式及其固有频率。

令$ {k^4} = \dfrac{{\rho A{\omega ^2}}}{{EI}} $，则式(1)可重新表示为

式(2)的解可表示为

式(3)中的系数A，B，C，D与梁的边界条件有关，除了简支梁，模态表达式中包含双曲正弦和双曲余弦函数，例如悬臂梁的第n阶模态形式表达式为

式中，L为梁的长度，k_n为第n阶波数，可通过如下频率方程计算得到^[1-4]

显然，当模态阶数n较大（n>4）时，式(5)中双曲正弦和双曲余弦函数值远大于三角函数值，所以

此时式(4a)也可简化为

从式(7)可以发现，随着模态阶数n的增加，当x趋于L时，$ {X_n}\left( L \right) \approx \cosh \left( {{k_n}L} \right) - \sinh \left( {{k_n}L} \right) $，这意味着需要通过两个非常大值相减得到一个较小的数值。由此引起的舍入误差会导致模态形状失真（特别是在x = L附近）。图1表示通过式(4)计算得到的悬臂梁第12~第14阶模态形状，从图中可以发现，当x/L > 0.9时，第12阶模态形状开始失真。并且随着模态阶数增加，模态形状的失真更加严重，并且失真向固定端扩展。两端固支、两端自由和固支-简支梁通过经典模态表达式^[1–4]计算时也存在类似问题，如图2所示。本文所有计算是在32位精度的个人计算机上完成。

图  1  通过式(4a)计算得到的悬臂梁第12~第14阶模态

下载: 全尺寸图片 幻灯片

图  2  不同边界条件下通过教材[1-4]计算Euler-Bernoulli梁的第12~第14阶模态

下载: 全尺寸图片 幻灯片

2.   模态表达式修正

注意到双曲余弦函数可表示为

还是以悬臂梁为例，将式(8)代入式(4b)，可得

式中，$ {\alpha _n} = \dfrac{{\exp \left( { - {k_n}L} \right) - \sin \left( {{k_n}L} \right) + \cos \left( {{k_n}L} \right)}}{{\sinh \left( {{k_n}L} \right) + \sin \left( {{k_n}L} \right)}} $。

将式(8)和式(9)代入式(4a)，整理可得

为了阐述本文方法可以避免数值计算中的舍入误差，将式(9)中的α_n代入式(10)，可得

从式(11)可以发现，随着模态阶数n的增加(n>4)，α_n趋于零，即$ 1 + {\alpha _n} \approx 1 $。并且$ \sinh \left( {{k_n}L} \right) \gg \sin \left( {{k_n}L} \right) $，式(11)中的双曲三角函数项$ \dfrac{{\sinh \left( {{k_n}x} \right)}}{{\sinh \left( {{k_n}L} \right) + \sin \left( {{k_n}L} \right)}} \approx \dfrac{{\sinh \left( {{k_n}x} \right)}}{{\sinh \left( {{k_n}L} \right)}} \leqslant 1 $，则式(11)可简化为

从式(12)可以发现，修正后的模态表达式在计算高阶模态时，从经典模态表达式中两个非常大值相减转换为两个非常大值相除，避免了经典模态表达式(4)在数值计算中的舍入误差，从而保证模态形状不失真。

对于两端固支、两端自由、固支-简支边界条件，通过式(8)代替文献[1-4]的经典模态表达式中的cosh(k_nL)，即可得到修正后的模态表达式，表1列出不同边界条件时的修正模态表达式。为了便于读者阅读，表中直接给出了k_nL的解。

表  1  不同边界条件时模态表达式

边界条件	模态表达式	αn	knL
两端固支	$ \begin{gathered} X_n^{\text{new}}\left( x \right) = \exp \left( { - {k_n}L} \right) - \cos \left( {{k_n}x} \right)- \\ \qquad {\alpha _n}\sinh \left( {{k_n}L} \right) + \left( {1 + {\alpha _n}} \right)\sin \left( {{k_n}x} \right) \\ \end{gathered} $	$ {\alpha _n} = \dfrac{{\exp \left( { - {k_n}L} \right) + \sin \left( {{k_n}L} \right) - \cos \left( {{k_n}L} \right)}}{{\sinh \left( {{k_n}L} \right) - \sin \left( {{k_n}L} \right)}} $	k1L = 4.73004 k2L = 7.85320 k3L = 10.99561 k4L = 14.13717 k5L = 17.27876 $k_{n} L = \left( {n + 0.5} \right)\pi \left( {n > 5} \right) $
自由-自由	$ \begin{gathered} X_n^{\text{new}}\left( x \right) = \exp \left( { - {k_n}L} \right) + \cos \left( {{k_n}x} \right)- \\ \qquad {\alpha _n}\sinh \left( {{k_n}L} \right) - \left( {1 + {\alpha _n}} \right)\sin \left( {{k_n}x} \right) \\ \end{gathered} $	$ {\alpha _n} = \dfrac{{\exp \left( { - {k_n}L} \right) + \sin \left( {{k_n}L} \right) - \cos \left( {{k_n}L} \right)}}{{\sinh \left( {{k_n}L} \right) - \sin \left( {{k_n}L} \right)}} $
固支-自由	$ \begin{gathered} X_n^{\text{new}}\left( x \right) = \exp \left( { - {k_n}L} \right) - \cos \left( {{k_n}x} \right)- \\ \qquad {\alpha _n}\sinh \left( {{k_n}L} \right) + \left( {1 + {\alpha _n}} \right)\sin \left( {{k_n}x} \right) \\ \end{gathered} $	$ {\alpha _n} = \dfrac{{\exp \left( { - {k_n}L} \right) - \sin \left( {{k_n}L} \right) + \cos \left( {{k_n}L} \right)}}{{\sinh \left( {{k_n}L} \right) + \sin \left( {{k_n}L} \right)}} $	k1L = 1.87510 k2L = 4.69409 k3L = 7.85476 k4L = 10.99554 k5L = 14.13717 $k_{n} L = \left( {n - 0.5} \right)\pi \left( {n > 5} \right) $
固支-简支	$ \begin{gathered} X_n^{\text{new}}\left( x \right) = \exp \left( { - {k_n}L} \right) - \cos \left( {{k_n}x} \right)- \\ \qquad {\alpha _n}\sinh \left( {{k_n}L} \right) + \left( {1 + {\alpha _n}} \right)\sin \left( {{k_n}x} \right) \\ \end{gathered} $	$ {\alpha _n} = \dfrac{{\exp \left( { - {k_n}L} \right) + \sin \left( {{k_n}L} \right) - \cos \left( {{k_n}L} \right)}}{{\sinh \left( {{k_n}L} \right) - \sin \left( {{k_n}L} \right)}} $	k1L = 3.92660 k2L = 7.06858 k3L = 10.21018 k4L = 13.35177 k5L = 16.49336 $k_{n} L = \left( {n + 0.25} \right)\pi \left( {n > 5} \right) $

下载: 导出CSV 

| 显示表格

为了进一步验证本文模态表达式在数值计算中的可行性，本文以两端固支梁为例，并与文献[8]结果进行对比。文献[8]给出两端固支梁的模态

式中$ {\sigma _n} = \dfrac{{1 + \exp \left( { - 2{k_n}L} \right) - 2\cos \left( {{k_n}L} \right)\exp \left( { - {k_n}L} \right)}}{{1 - \exp \left( { - 2{k_n}L} \right) - 2\sin \left( {{k_n}L} \right)\exp \left( { - {k_n}L} \right)}} $

定义模态形状的计算误差为

图3表示分别通过本文模态表达式（见表1）与文献[8]模态表达式所得两端固支梁高阶模态形状，为了便于对比，图3同时给出通过文献[1-4]中经典模态表达式的计算结果。图4表示两端固支梁前200阶本文模态表达式的计算误差。从图中可以发现，通过经典模态表达式计算得到的高阶模态形状存在失真，而本文结果与文献[8]结果完全吻合，前200阶模态的计算误差小于1.8 × 10^–16。这进一步验证了本文模态表达式的正确性，对比表1和式(13)也可发现，本文模态表达式更为简便。

图  3  两端固支梁的高阶模态形式

下载: 全尺寸图片 幻灯片

图  4  两端固支梁本文模态表达式与文献[8]对比的计算误差

下载: 全尺寸图片 幻灯片

图5表示经典模态表达式^[1-4]与本文模态表达式的绝对值之差$ \left| {X_n^{\text{new}}\left( x \right) - X_n^{}\left( x \right)} \right| $，从图中可以发现，当模态阶数n小于12时，本文结果与经典模态表达式结果之间的误差为零。而随着模态阶数n的增加，两者之间误差从x/L = 1开始向x/L = 0.06扩展。这主要是因为经典模态表达式计算所得的模态形状失真引起（如图2和图3所示）。

图  5  经典模态表达式^[1-4]与本文模态表达式的绝对值之差$ \left| {X_n^{\text{new}}\left( x \right) - X_n\left( x \right)} \right| $

下载: 全尺寸图片 幻灯片

3.   总结

注意到经典机械振动教材中Euler-Bernoulli梁的模态表达式并不适合高阶模态的数值计算，本文提出对模态表达式进行修正，从而可以精确计算至少前200阶模态。与前人结果相比，本文修正模态表达式相对简单，便于课堂教学。本文结果可作为现有机械振动教学内容的补充，并希望能对《机械振动基础》课程的教学实践有所帮助。