传递矩阵法分析平行四边形板的自由振动问题

刘灿礼 袁丽芸 王俊鹏 施国强

摘 要：针对平行四边形薄板的自由振动问题，采用传递矩阵法，结合摄动理论分析平行四边形板，并通过与有限元法对比验证该方法的精确性.首先，利用坐标变换公式将平行四边形板的求解域转换到矩形板域，并将控制方程及其对应的边界条件变换到矩形板域内；随后，结合Fourier级数展开以及摄动理论，将平形四边行板的振动控制方程写为一階常微分矩阵的形式；再运用传递矩阵法对该一阶常微分矩阵方程进行求解，得到平形四边行板的振动固有频率.算例结果表明，在平行四边形板偏角较小的范围内，所提出的传递矩阵法可精确应用于平行四边形板的动力学特性分析.

关键词：平行四边形板；传递矩阵法；摄动法；一阶常微分矩阵方程

中图分类号：U462.31 文献标志码：A

0 引言

在弹性薄板的弯曲振动问题中，形状规则的矩形板和圆形板较易获得解析解.但对于平行四边形板的振动分析，因边界条件不易满足，其解析解极难获得，因而常采用数值解法或近似函数解法进行近似逼近.在数值解法中，有限元法边界适用性最佳，能得到较精确的结果，但进行参数化分析通常需要重复建模.在近似函数解法中，常用的有变分法.文献[1]用康托洛维奇变分法对平行四边形板进行了分析，但由于在近似函数中采用了均布载荷下的梁函数，因而仅适用于均布载荷的情况，其他载荷情形需要重新求解，且求高阶近似解的计算工作量也相当大.杨柳等[2]采用GD方法对平行四边形板的静力学问题进行了研究，该方法数学原理严谨，精度较高，可拓展到动力学问题的求解，但其本质仍然是用差分代替微分，存在一定的数值误差.

上个世纪20年代以来，由于具有力学概念清晰、精度高、易于编程计算等优点，传递矩阵法在结构振动中广泛应用[3] ，该方法的求解思路是将结构控制微分方程写成一阶常微分方程组的形式，则起点与终点的状态向量间即可用传递矩阵建立简单关系，利用边界条件对结构振动进行求解[4].由于传递矩阵法中要计算矩阵运动微分方程中系数矩阵的指数函数的幂，导致其精度不高，为此，向宇[5]运用微分方程与矩阵分析理论，结合精细积分技术，精确地计算了指数函数的幂，大大提高了传递矩阵法的计算精度.此后，传递矩阵法被广泛应用于复合多层梁[6]以及旋转壳[7]等规则结构的振动.该方法要求结构振动方程可降维写成一阶常微分矩阵方程的形式，因此在形状不规则的结构中的应用受到一定的限制.

基此，针对四边简支平行四边形板的自由振动问题，本文拟采用坐标变换将平行四边形板的求解域变换为矩形域，此后将借助摄动理论，忽略方程中的小项，采用傅立叶级数展开和无量纲化处理，得到平行四边形板自由振动的一阶常微分矩阵方程，并拟借助高精度的传递矩阵法来求解.此方法为该类结构的参数化研究提供一个有力的手段.

1 平行四边形板的动力学问题描述

考虑一厚度为h的各向同性均匀薄板，坐标系位于中面（x，y）上.板中面x，y，z三个方向的位移分量分别记为us，vs，ws；u，v，w分别为板内任意一点的3个位移分量.由薄板理论，可得其与中面位移分量的关系如下：

u=us-z■， v=vs-z■， w=ws （1）

与文献[8]不同的是，式（1）中薄板中面内的位移不为0，即薄板在弯曲振动动中，中面内的各点位移有平行于中面的位移，且变形后的中面不再是平面.

结合薄板的几何方程，以及广义胡克定律，将其内力向中面进行简化可得内力-位移关系：

（2）

其中，K=■，D=■ ；E，μ分别为板的弹性模量和泊松比，Nx，Ny分别为x，y方向的中面单位长度的薄膜内力，Nxy为中面单位长度的薄膜剪力；Mx，My分别为x，y方向的中面单位长度弯矩，Mxy为中面单位长度扭矩.

由板壳振动理论[9]，对薄板进行受力分析可得平衡方程：

（3）

式（3）中，Qx，Qy分别为沿x，y方向的中面单位长度横向剪力；pi（i=x，y，z）为作用在单元体中的单位面积体力.

引入等效剪力：

Sx=Nxy，Vx=Qx+■ （4）

谐激励下，式（1）～式（4）中的各物理量可展开成幅值与谐激励因子ejwt乘积的形式，即：

f（x，y，t）=■（x，y，w）ejwt （5）

则可得谐激励下板的平衡方程如下：

（6）

同理可得， 内力-位移关系（式（2））和等效剪力（式（4））在谐激励下的形式，仅需将各物理量用其幅值替代即可.由于篇幅有限，不再一一列出.

为了便于使用边界条件，需将平形四边行板的求解域进行坐标变换成矩形板域.记α为平形四边行板的偏角，采用如图1所示坐标变换，斜坐标（x，y）与直角坐标系（ξ，η）之间的转换关系为：

x=ξ+ηsinα， y=ηcosα （7）

由链式求导法则，有：

■=■ ■， ■=■ ■+■ ■ （8）

将式（8）代入谐激励作用下的式（2）、式（4）、式（6）可得（ξ，η）坐标系下平行四边形板的内力-位移关系如下：

（9）

运动控制方程如下：

（10）

等效剪力为：

■ξ η=■ξη ，■ξ=■ξ -tanα■+secα■ （11）

由于式（9）～式（11）均为偏微分方程，无法直接转换为常微分方程，且在α较小时，tanα→0，可将上述方程中所含tanα项消除；因此，采用摄动理论，取tanα为摄动参数，将式（9）～式（11）中各物理量展开为W=W0+W1tanα+W2tan2α+…的形式，并只保留第一项W0作为主项，误差应该不大.

采用矩形板类似的求解方法，令ξ=Lc，（L=a），然后将各物理量沿η方向进行傅里叶级数展开并进行无量纲化可得：

（■s，■s ）=L■（■s n，■s n）sin■，（■s ）=L■（■s n）cos■，（■ξ， ■η， ■ξ， ■ξ ）=K■ （■ξ n，■η n，■ξ n，■ξ n）sin■，

（■ξ η，■η ξ，■η ，■ξ ）=K■（■ξ ηn， ■η ξn， ■η n， ■ξ n）cos■， （■ξ ， ■η ）=（LK）■（■ξ n， ■η n）sin■，

（■ξ η， ■η ξ）=（LK）■（■ξ η n，■η ξn）cos■， （■ξ ， ■z ）=■（■ξ n， ■z n）sin■， （■η ）=■■η ncos■，

（■ξ ）=■■x ncos■ （12）

把式（12）代入只保留主项后的式（9）～式（11）可得平形四边行板的内力-位移关系：

（13）

平行四边形板的平衡方程：

（14）

等效剪力：

■c n0=■c ηn0，■c n0=■c n0-■■cηn0 （15）

并引入无量纲转角状态变量的Fourier级数展开分量：

■c n0=-■ （16）

从式（13）～式（16）中消去其他中间变量，经推导整理后可得到平行四边形板的一阶常微分矩阵方程：

■=AZ0+F （17）

其中，Z0={■sn0，■sn0，■s n0，■c n0，■c n0，■c n0，■c n0，■c n0}T为板中面状态变量的Fourier级数展开分量，

A=

为系数矩阵，其非零元素如下：

A12=■，A15=1 ，A21=-■ ，A26=■ ，A34=-1 ，A43=-■ ，

A48=■ ，A51=-■ ，A56=■ ，A62=-■+■ ，

A65=-■，A73=-■+■ ，A78=■ ，A84=■ ，A87=1

F8×1=■[0，0，0，0，■c n0，■η n0，■z n0，0]T为所受外激励.

2 传递矩阵法的动力学问题求解

本文研究平行四边形板的自由振动特性不考虑外力，此时，式（17）中的F为零向量，由传递矩阵法，边界两端变量之间存在的传递关系：

Z1=TZ0 （19）

其中，Z0为起始状态变量，Z1为终端状态变量；T=eA為传递矩阵；引入文献[10]描述的精细积分法，当系数矩阵A已知时，可精确计算T的值.

通过坐标变化后，四边简支的平行四边形薄板的求解域可转换为一矩形域，其在矩形域两边的边界条件为：

c=0，1+■sinα： us=0， Ncn0≠0， vsn0≠0， Scn0=0， ws=0， Vcn0≠0， Mcn0=0， θcn0≠0 （20）

结合式（19）和式（20），可得平形四边行板自由振动时的特征方程：

（21）

式（21）中，矩阵TT中的元素Tij表示传递矩阵T的第i行第j列元素.

要使特征方程（式（21））有非零解，则特征方程系数矩阵的行列式│TT│=det（TT）=0，这样就可以得到一个关于板的固有圆频率ω的方程，其求解过程可以利用MATALAB编程实现.此后，可得平行四边形薄板的振动固有频率如下：

f=■ （22）

3 算例分析

为了验证本文方法的正确性，下面将取一个偏角较小的平形四边行薄板为例，进行动力学分析.

算例1 考虑尺寸为2.00 m×1.00 m×0.02 m，偏角α=π/36 的四边简支平行四边形薄板，其物理参数如表1所示.

通过本文方法和有限元仿真求得的（1，1）阶振动固有频率的解如图2和图3所示.

其它阶次固有频率的本文解和有限元解如表2所示.从表2可以看出，随着阶数的提高，两种方法所得固有频率的差别逐渐增大，但本文提出的传递矩阵法求得的解和有限元解的误差均低于2%，表明传递矩阵法具有较高的求解精度.误差存在的主要原因可能在于忽略了偏角带来的其他小项；因此，需要进一步讨论偏角的影响.

算例2 考虑平形四边行板的 a， b，h 不变，但偏角发生变化时的情形（偏角分别取为α=π/36，α=π/18，α=π/9，α=π/6）；因此可得不同角度下本文解和有限元解，如表3所示.

算例2表明，随着角度的增大，传递矩阵法的解和有限元解的误差总体上是逐渐增大的，并且阶数越高，误差越大.当角度小于π/9 （20°）时，两者的误差可控制在2%左右，为了节省篇幅，其余角度进行的算例验证不再一一列表给出.通过不同角度下本文解和有限元解的对比表明：角度一旦超过π/9 （20°），误差较大.

4 结语

采用传递矩阵法和摄动理论，本文提出了一种分析平行四边形板的动力学特性的新方法.算例表明，在偏角较小（0～20°） 范围内，该方法和有限元法求得的固有频率基本吻合.说明该方法可应用于平行四边形薄板动力特性问题的求解，为其参数化研究提供了一个便利的手段.

参考文献

[1] 王磊.平行四边形弯曲板的康托洛维奇法[J]. 固体力学学报，1983，4（3）：420-426.

[2] 杨柳，彭建设.解平行四边形板弯曲问题的GD法[J]. 成都大学学报（自然科学版），2014，33（3）：230-233.

[3] 朱镜清.结构振动计算中传递矩阵法的应用问题[J]. 地震工程与工程振动，1984，4（1）：36-44.

[4] 万浩川，李以农，郑玲.改进的结构振动传递矩阵法[J]. 振动与冲击，2013，32（9）：173-177，82.

[5] 向 宇.分析结构自由振动的传递矩阵精确形式[J]. 振动与冲击，1999，18（2）：69-74.

[6] 吴寒剑，向宇，袁丽芸.PCLD梁的振动和阻尼特性分析[J]. 广西科技大学学报，2014，25（4）：40-45，77.

[7] 张又军，陆静.环向加肋旋转壳自由振动特性分析[J]. 广西科技大学学报，2014，25（1）：54-58.

[8] 徐芝纶.弹性力学[M]. 北京：高等教育出版社，2006.

[9] 曹志远.板壳振动理论[M]. 北京：中国铁道出版社，1989.

[10] 向宇，黄玉盈，曾革委.精细时程积分法的误差分析与精度设计[J].计算力学学报，2002，8（3）：276-280，319.

Abstract：The transfer matrix method is used to analyze the free vibration of parallel quadrilateral plate， and the accuracy of the method is verified by comparing with the finite element method. Firstly， using coordinate transformation formula， the solution domain with the parallel quadrilateral shape is converted into the rectangular one， and the governing equations and the corresponding boundary conditions in the transformed rectangular domain are obtained； then， combined with the Fourier series expansion and perturbation theory， the governing equation for the plate is written in the first order ordinary differential matrix form； next， the transfer matrix method is used to solve the matrix equations to get the natural frequency of the parallel quadrilateral plate. Example results show that when the parallelogram plate's angle is in a small range， the transfer matrix method proposed in this paper can be applied to the analysis of dynamic characteristics of a parallel quadrilateral plate.

Key words：parallel quadrilateral plate； transfer matrix method； perturbation method； the first order ordinary differential matrix equation

（學科编辑：张玉凤）