马辉^1,2, 崔璨¹, 曾劲¹, 闻邦椿¹
1. 东北大学 机械工程与自动化学院, 辽宁 沈阳 110819;
2. 上海交通大学 机械结构强度与振动国家重点实验室, 上海 200240
收稿日期：2018-08-15
基金项目：国家自然科学基金资助项目(11772089)；中央高校基本科研业务费专项资金资助项目(N170308028, N180708009)；辽宁省高等学校创新人才支持计划项目(LR2017035)。
作者简介：马辉(1978-), 男, 河北安平人，东北大学教授，博士生导师;
闻邦椿(1930-)，男，浙江温岭人，东北大学教授，博士生导师，中国科学院院士。

摘要：利用解析方法对考虑弹性支承的凸肩叶片进行动力学建模, 在解析模型中考虑叶片的旋转效应、安装角、扭转角以及凸肩位置处截面变化的影响.利用ANSYS软件建立了相对应的有限元模型, 验证了所建立的解析模型的有效性.通过所建立的解析模型和有限元模型, 分析了叶根处支承刚度和凸肩位置对叶片动力学特性的影响.研究结果表明:随着叶根处支承刚度的增加, 叶片的固有频率升高, 但支承刚度增加到一定程度后, 支承刚度的增加不再引起固有频率的升高; 随着凸肩位置至叶尖距离的增加, 叶片的第1阶和第3阶固有频率有升高的趋势, 而第2阶固有频率则会随之降低.此外, 采用解析方法能够明显提高计算效率.
关键词：凸肩叶片弹性支承变截面动力学特性
Analysis of Dynamic Characteristics of the Rotating Shrouded Blade Considering Elastic Support
MA Hui^1,2, CUI Can¹, ZENG Jin¹, WEN Bang-chun¹
1. School of Mechanical Engineering & Automation, Northeastern University, Shenyang 110819, China;
2. State Key Laboratory of Mechanical System and Vibration, Shanghai Jiaotong University, Shanghai 200240, China
Corresponding author: MA Hui, E-mail: huima@me.neu.edu.cn
Abstract: An analytical model of the shrouded blade considering elastic support is established, in which the influence of the blade's rotational effect, setting angle, twist angle and variable cross-section are taken into consideration. ANSYS software is used to establish a corresponding finite element model to verify the effectiveness of the established analytical model. The influence of the support stiffness and shoulder position on the dynamic characteristics of the blade is analyzed through the established analytical model and finite element model. The results indicate that as the support stiffness increases, the natural frequency of the shrouded blade increases. However, after the support stiffness increases to a certain extent, the increase of the support stiffness will no longer cause the increase of the natural frequency. In addition, as the distance from the shoulder position to the blade tip increases, the first-order and third-order natural frequency of the shrouded blade tends to increase. In contrast, the second-order natural frequency will decrease. In addition, analytical methods can significantly improve the computational efficiency.
Key words: shrouded bladeelastic supportvariable cross-sectiondynamic characteristics
为了限制叶片强迫振动的振幅, 提高叶片适应冲击载荷的能力, 以及防止长叶片在流固耦合作用下产生的自激振动, 一种有效的方法就是在系统中采用特殊结构, 例如在航空发动机中采用叶冠、凸肩及缘板阻尼块等结构^[1].由于凸肩结构不受温度限制、结构简单有效而得到广泛应用, 因此, 深入研究凸肩结构对叶片振动的影响规律具有重要的理论意义和工程应用价值.洪杰等^[2]基于非线性接触有限元法, 研究分析了凸肩径向位置对接触状态和风扇叶片振动特性的影响规律.李红影等^[3]应用Donnell’s简化壳理论建立了壳板凸肩叶片的非线性振动微分方程, 并考虑几何非线性、阻尼、凸肩接触面正压力、摩擦力等因素研究了板壳凸肩叶片的静态和全局动态分岔.李健等^[4]用薄板-弹簧系统模拟组合凸肩叶片法.马静敏和任勇生^[5]利用扩展拉格朗日方程, 建立了翼型截面叶片的自由振动微分方程, 并通过Galerkin方法求解了不同翼型截面叶片的固有频率.
在实际工程应用中, 涡轮叶片一般带有一定的预扭角, 而该预扭角会对叶片的动力学特性产生一定的影响, 所以越来越多的学者开始在研究中考虑叶片预扭角的影响^[6-9].Yoo等^[7]提出一种针对带有集中质量点的旋转预扭叶片的建模方法, 并分析了集中质量点的质量和位置、叶片预扭角等参数对旋转叶片固有特性的影响.
通常来说, 固定支承是理想的约束条件, 所以对于具有弹性支承约束的叶片的振动分析更具有工程意义.Choi和Chou^[10]分析了在弹性支承条件下预扭叶片的振动响应, 并提供了一种计算组合凸肩叶片固有特性的方法.Choi等^[11]利用微分求积法研究了考虑旋转效应的Timoshenko梁的固有特性, Timoshenko梁的边界条件为自由、简支、弹性支承等不同边界条件的自由组合.叶茂等^[12]根据弹性支承点处的协调条件, 基于转换矩阵法, 给出了两端简支、两端固定、悬臂梁、两端自由边界条件下欧拉连续梁的特征方程.Naguleswaran^[13]研究了在线性变化的轴向力作用下的欧拉-伯努利梁横向振动, 给出了横向振动的振型函数求解方法, 并列举了多种经典边界条件下的频率方程.
在已发表的文献中, 很少有解析方法考虑变截面的扭形旋转凸肩叶片在弹性支承条件下的动力学特性.为了弥补这一不足, 本文通过解析方法建立了能够考虑以上各种效应的动力学模型, 并利用有限元软件验证了所建立模型的有效性.在此基础上分析了支承刚度和凸肩位置对叶片动力学特性的影响规律.
1 凸肩叶片动力学模型建立旋转凸肩叶片如图 1所示.图中, l₂是凸肩的厚度；h_s是凸肩的宽度；L, b, h分别为叶片的长度、宽度和高度, L=l₁+l₂+l₃.叶片在弯曲方向约束弹簧和约束扭簧的刚度分别为k_y, k_ry, 叶片在摆动方向约束弹簧和约束扭簧的刚度分别为k_z, k_rz, Oxyz为叶片叶根处的坐标系.叶根相对于叶盘的安装角为β₀, 叶尖处相对于叶根的扭转角为β_t, 叶片上任意截面相对于叶根处的扭转角为γ(x), 则γ(x)= β_tx/L, 所以可得叶片上任意截面相对于叶盘的扭转角为β(x)= β₀+γ(x).
图 1(Fig. 1)

图 1 旋转凸肩叶片示意图Fig.1 Schematic of the rotating shrouded blade

凸肩叶片上任意一点的位置在相对于叶盘的整体坐标系中的表达式见式(1), 其中, u, v, w分别为凸肩叶片上任意一点在径向、弯曲、摆动方向的变形；φ, ?分别为叶片弯曲方向和摆动方向的转角；θ是叶片绕旋转轴旋转的角位移；R_d为叶盘半径.

(1)

进而可得凸肩叶片总的动能:

(2)

式中m和ρ分别是凸肩叶片质量和密度；A(x)是凸肩叶片的截面面积:

(3)

考虑凸肩叶片在旋转过程中的离心刚化效应, 并考虑叶片由于剪切变形产生的剪切应变能，以及约束弹簧和约束扭簧的弹性势能, 凸肩叶片的整体势能如式(4)所示^[14], 其中, E, G分别为弹性模量和剪切模量；κ_y和κ_z分别是叶片对y和z的剪切系数；符号(′)表示一阶偏导.叶片所受离心力f_c(x)如式(5)所示，式中ρ是叶片材料密度，ω为气动载荷圆频率.

(4)

(5)

叶片截面绕y轴的截面惯性矩I_y(x), 叶片截面绕z轴的截面惯性矩I_z(x)的表达式如下:

(6a)

(6b)

作用在凸肩叶片上的外力所做的功W_non为

(7)

式中F_e是凸肩叶片所受的气动载荷, 是沿y方向的均布载荷:

(8)

式中:F₀是气动载荷的幅值; k_e是静子叶片数, 本文k_e=30;ω是气动载荷的圆频率, 与转速的关系为ω=2πΩ /60.
根据Hamilton能量方程：
可得到凸肩叶片的运动方程.
由于凸肩叶片运动方程是时变非线性耦合微分方程, 在直接求解过程中不易得到精确的解析解.在这种情况下, 常用方法是将连续方程进行离散化处理, 以此求解得到方程的近似解.本文采用Galerkin方法对运动微分方程进行离散化处理.
根据凸肩叶片截面的变化, 可将凸肩叶片分为三段等截面梁, 三段叶片在弯曲方向的振型函数可假设为式(9)形式.此外, 由于梁的连续性, 梁的第i段和第i+1段在连接点处的位移、转角、弯矩、剪力均相等, 所以可得到边界条件如式(10)所示^[15].将式(9)代入式(10)可得:P_iN_i=P_i+1N_i+1, 其中, N_i和N_i+1是振型函数中的待定系数, N_i=[A_i, B_i, C_i, D_i]^T, N_i+1=[A_i+1, B_i+1, C_i+1, D_i+1]^T, P_i和P_i+1的表达式见式(11).

(9)

(10)

(11)

由此可得相邻两段梁振型函数系数之间的关系:

(12)

变截面梁作为一个连续体, 其三段梁的固有频率是相等的, 故有以下关系式:

(13)

为了使计算简便, 在每一段梁上将x作变量代换为(l_i-x), 根据弹性支承梁的边界条件可得式(14).其中, 由前两个条件可得:A₁=C₁, B₁=D₁.由后两个条件可得ΛN₃=ΛZ₂Z₁=Q=0, 其中Λ的表达式如式(15)所示.

(14)

(15)

由以上分析可得

(16)

其中, Q₁₁, …, Q₂₄是矩阵Q中各元素的值.由于A₁, B₁存在非零解, 所以其系数行列式等于0, 即

(17)

为了得到变截面梁的振型函数, 令A₁=1;根据三段振型函数系数之间的递推关系, 可得到变截面梁每一段的振型函数.利用同样的方法可求得凸肩叶片在摆动方向的振型函数.
利用Galerkin方法对凸肩叶片的运动方程进行离散化处理, 引入正则坐标U_i(t), V_i(t), W_i(t), Ψ_i(t)和Φ_i(t), 将凸肩叶片的径向、弯曲方向、摆动方向、弯曲转角, 以及摆动转角方向位移写为

(18)

式中, ?_1i(x), ?_2i(x), ?_3i(x), ?_4i(x)和?_5i(x)分别是凸肩叶片的径向、弯曲方向、摆动方向、弯曲转角及摆动转角方向的振型函数:

(19)

将式(18)和式(19)代入到凸肩叶片运动方程中, 并将其写成矩阵形式:

(20)

式中, M, G, D, K_e, K_c, K_s, K_acc, K_t, q和F分别是质量矩阵、科氏力矩阵、瑞利阻尼矩阵、结构刚度矩阵、离心刚化矩阵、旋转软化矩阵、加速度导致的刚度矩阵、弹簧刚度矩阵以及正则坐标系下的位移向量和外激振力向量.其中瑞丽阻尼矩阵的表达式如下:

(21)

其中,
式中：f_n1和f_n2表示凸肩叶片的前两阶固有频率, ξ₁和ξ₂是对应的模态阻尼比, 本文中分别取为0.02和0.04.
2 凸肩叶片动力学特性分析基于上文所建立的解析模型和有限元软件, 本节讨论凸肩叶片叶根处的支承刚度和凸肩位置对旋转凸肩叶片固有特性的影响.在利用有限元软件ANSYS所建立的模型中, 采用Beam188单元模拟凸肩叶片, Combin14单元模拟叶片弹性支承.凸肩叶片有限元模型如图 2所示.
图 2(Fig. 2)

图 2 凸肩叶片有限元模型Fig.2 Finite element model of the shrouded blade

2.1 支承刚度的影响叶片振动过程中以弯曲方向的振动为主, 因此摆动方向弹簧和扭簧刚度取固定值, 分析弯曲方向弹簧和扭簧刚度取不同值时凸肩叶片的固有特性.
取摆动方向弹簧刚度k_z=2×10⁸N/m, 扭簧刚度k_rz=2×10⁷ N·m/rad, 凸肩叶片的安装角β₀=30°, 扭转角β_t=10°.在弯曲方向上不同支承刚度条件下解析模型与有限元模型的动频曲线如图 3所示, 其中, 当弹簧刚度k_y=2×10⁶ N/m, 扭簧刚度k_ry=2×10⁵ N·m/rad时, 前两阶频率为叶片的弯曲频率, 第3阶频率为叶片的摆动频率; 在其他刚度条件下, 第1阶和第3阶频率是叶片的弯曲频率, 第2阶频率是叶片的摆动频率.各阶动频的最大误差如表 1所示.不同支承刚度下解析模型与有限元模型在叶尖位置处y方向的谐响应曲线如图 4所示.
图 3(Fig. 3)

图 3 不同支承刚度下的动频曲线Fig.3 Natural frequencies at different support stiffness

表 1(Table 1)

表 1 不同支承刚度下动频最大误差Table 1 Maximum error of natural frequencies at different support stiffness % 支承刚度 1阶 2阶 3阶 ky=2×106N/m, kry=2×105N·m/rad 0.51 0.11 0.66 ky=2×107N/m, kry=2×106N·m/rad 0.55 0.48 0.51 ky=2×108N/m, kry=2×107N·m/rad 0.09 0.58 0.51 ky=2×1010N/m, kry=2×109N·m/rad 0.12 0.58 0.51

表 1 不同支承刚度下动频最大误差 Table 1 Maximum error of natural frequencies at different support stiffness

图 4(Fig. 4)

图 4 不同支承刚度下的谐响应曲线Fig.4 Harmonic responses at different support stiffness

由图 3和图 4可以看出, 随着叶根处支承刚度的增加, 凸肩叶片的各阶固有频率均增加, 从而导致了凸肩叶片的共振峰向右偏移且共振幅值降低.但当支承刚度增加到一定值时, 叶根处支承刚度的增加对叶片固有频率的影响不再明显, 这是因为当支承刚度增加到一定程度之后, 叶根处即相当于固支约束, 所以支承刚度继续增加几乎不影响凸肩叶片的固有频率, 即当支承刚度足够大之后, 谐响应曲线几乎重合.
此外, 为了比较解析法和有限元法的计算效率, 记录了两种方法在计算k_y=2×10⁶ N/m, k_ry=2×10⁵ N·m/rad时的计算时间.其中, 在计算固有频率时, 解析法耗时10.56 s, 有限元法耗时82.00 s; 在计算谐响应时, 解析法耗时12.94 s, 有限元法耗时182.12 s.由此可以看出, 解析方法的计算效率有明显优势, 尤其是谐响应的计算.
2.2 凸肩位置的影响摆动方向支承刚度、安装角、扭转角等参数与上节相同, 取弯曲方向弹簧刚度k_y=2×10⁶N/m, 扭簧刚度k_ry=2×10⁵N·m/rad, 改变凸肩位置, 比较其对叶片固有频率的影响.以l表示凸肩中心位置距叶尖的距离, 在不同凸肩位置下解析模型与有限元模型的动频曲线如图 5所示, 各阶动频的最大误差如表 2所示.不同凸肩位置下解析模型与有限元模型在叶尖位置处y方向的谐响应曲线如图 6所示.
图 5(Fig. 5)

图 5 不同凸肩位置下的动频曲线Fig.5 Natural frequencies at different shoulder positions

表 2(Table 2)

表 2 不同凸肩位置下动频最大误差Table 2 Maximum error of natural frequencies at different shoulder positions % 距离 1阶 2阶 3阶 l=L/4 0.57 0.17 0.56 l=L/2 0.51 0.11 0.66 l=3L/4 0.44 0.16 0.85

表 2 不同凸肩位置下动频最大误差 Table 2 Maximum error of natural frequencies at different shoulder positions

图 6(Fig. 6)

图 6 不同凸肩位置下的谐响应曲线Fig.6 Harmonic responses at different shoulder positions

由图 5和图 6可以看出, 随着l的增加, 第1阶固有频率和第3阶固有频率均随之增加, 其中, 第1阶固有频率和第3阶固有频率分别是叶片的第1阶弯曲频率和第1阶摆动频率.叶片的第2阶固有频率是叶片的第2阶弯曲频率, 并且叶片的第2阶固有频率会随l的增加而减小.由于在谐响应所取的转速范围内只激发出了凸肩叶片的第1阶固有频率, 所以随l的增加, 共振转速有一定程度的增加, 并且共振幅值有一定程度的降低.
3 结论考虑叶片在叶根处的弹性支承, 利用Hamilton能量原理和Galerkin方法建立了凸肩叶片的解析模型, 在所建立的模型中考虑叶片的旋转效应以及叶片的安装角和扭转角的影响.利用ANSYS软件建立了相对应的有限元模型, 验证了所建立解析模型的有效性, 并分析了叶根处支承刚度和凸肩位置对叶片动力学特性的影响规律, 得到如下结论:
1) 随着叶根处支承刚度的增加, 凸肩叶片的各阶固有频率均随之增加, 由此导致了叶片共振转速的增加和共振幅值的降低.但是当支承刚度接近于固支刚度时, 支承刚度的增加对叶片固有频率的影响变得不再明显.
2) 随着凸肩位置距叶尖处距离的增加, 叶片的第1阶弯曲频率和第1阶摆动频率均有增加的趋势.
3) 本文所建立的解析模型不仅在求解精度上与有限元模型接近, 并且在求解效率上相对于有限元模型有明显的优势.
参考文献

[1]	吴绵绵.涡轮叶片摩擦阻尼减振设计方法研究[D].南京: 南京航空航天大学, 2010. ( Wu Mian-mian.A vibration-reduced design method for turbine blade with friction damper[D].Nanjing: Nanjing University of Aeronautics and Astronautics, 2010.http://cdmd.cnki.com.cn/Article/CDMD-10287-1011252176.htm)
[2]	洪杰, 文敏, 马艳红, 等. 凸肩径向位置对风扇叶片振动特性的影响[J].航空动力学报, 2015, 30(12): 2817–2823. ( Hong Jie, Wen Min, Ma Yan-hong, et al. Influence of the shroud's radial position on vibration characteristics of the fan blade[J].Journal of Aerospace Power, 2015, 30(12): 2817–2823.)
[3]	李红影, 谢里阳, 郭星辉. 壳板凸肩叶片的分岔特性研究[J].振动工程学报, 2012, 25(3): 260–267. ( Li Hong-ying, Xie Li-yang, Guo Xing-hui. A study on the bifurcation characteristics for shrouded blade of cylindrical panel[J].Journal of Vibration Engineering, 2012, 25(3): 260–267.DOI:10.3969/j.issn.1004-4523.2012.03.006)
[4]	李健, 郭星辉, 魏伟明. 组合凸肩叶片的固有特性分析[J].工程力学, 2007, 24(8): 151–155. ( Li Jian, Guo Xing-hui, Wei Wei-ming. Analysis of inherent property on assembling shroud blades[J].Engineering Mechanics, 2007, 24(8): 151–155.DOI:10.3969/j.issn.1000-4750.2007.08.027)
[5]	马静敏, 任勇生. 风力机叶片的自由振动特性分析[J].振动与冲击, 2015, 34(17): 105–110. ( Ma Jing-min, Ren Yong-sheng. Free vibration characteristics of wind turbine blades[J].Journal of Vibration and Shock, 2015, 34(17): 105–110.)
[6]	Carnegie W. Vibrations of pre-twisted cantilever blading allowing for rotary inertia and shear deflection[J].Journal of Mechanical Engineering Science, 1964, 6(2): 105–109.DOI:10.1243/JMES_JOUR_1964_006_019_02
[7]	Yoo H H, Kwak J Y, Chung J. Vibration analysis of rotating pre-twisted blades with a concentrated tip mass[J].Journal of Sound and Vibration, 2001, 240(5): 891–908.DOI:10.1006/jsvi.2000.3258
[8]	汤凤, 孟光. 带冠涡轮叶片的接触分析[J].噪声与振动控制, 2005, 25(4): 5–7. ( Tang Feng, Meng Guang. Contact analysis of shrouded turbine blade[J].Noise and Vibration Control, 2005, 25(4): 5–7.DOI:10.3969/j.issn.1006-1355.2005.04.002)
[9]	Xie F, Ma H, Cui C, et al. Vibration response comparison of twisted shrouded blades using different impact models[J].Journal of Sound and Vibration, 2017, 397: 171–191.DOI:10.1016/j.jsv.2017.02.056
[10]	Choi S T, Chou Y T. Vibration analysis of elastically supported turbomachinery blades by the modified differential quadrature method[J].Journal of Sound and Vibration, 2001, 240(5): 937–953.DOI:10.1006/jsvi.2000.3267
[11]	Choi S T, Wu J D, Chou Y T. Dynamic analysis of a spinning Timoshenko beam by the differential quadrature method[J].AIAA Journal, 2012, 38(5): 851–856.
[12]	叶茂, 谭平, 任珉, 等. 中间带弹性支承各种边界条件连续梁模态分析[J].工程力学, 2010, 27(9): 80–85. ( Ye Mao, Tan Ping, Ren Min, et al. Modal analysis of multi-span beams with intermediate flexible constraints and different boundary conditions[J].Engineering Mechanics, 2010, 27(9): 80–85.)
[13]	Naguleswaran S. Transverse vibration of an uniform Euler-Bernoulli beam under linearly varying axial force[J].Journal of Sound and Vibration, 2004, 275(1): 47–57.
[14]	Sabuncu M, Evran K. Dynamic stability of a rotating pre-twisted asymmetric cross-section Timoshenko beam subjected to an axial periodic force[J].International Journal of Mechanical Sciences, 2006, 48: 579–590.DOI:10.1016/j.ijmecsci.2006.01.010
[15]	崔灿, 李映辉. 变截面铁木辛柯梁振动特性快速计算方法[J].动力学与控制学报, 2012, 10(3): 258–262. ( Cui Can, Li Ying-hui. A solution for vibration characteristic of Timoshenko beam with variable cross-section[J].Journal of Dynamics and Control, 2012, 10(3): 258–262.)