韩 伟,齐 龙,肖中云
(中国空气动力研究与发展中心,四川 绵阳 621000)
旋翼是直升机不可或缺的重要部件,其动力学特性对于直升机的动力稳定性和安全性具有不可忽视的影响,对于旋翼,一些看起来不重要的设计变化都可能导致动力学特性的显著变化[1]。为了保证满意的动力学特性,就必须在研制的过程中予以控制,而这些特性主要指其模态特性:固有频率和振型等。
桨叶在实际振动中,三个自由度通常是耦合的,其传统的固有振动特性计算方法主要包括:假设模特法、传递矩阵法和有限元法。最早使用的假设模态法是假设一组满足相应边界条件的函数序列,去逼近要求的模态函数,应用能量原理求得桨叶固有频率和振型的方法,此法主要包括拉格朗日法、瑞利法和伽辽金法,这种方法很不方便,目前已经很少使用。文献[2-3]基于传递矩阵法分析了旋翼的固有振动特性,虽然该法简单方便,计算效率较高,但这种方法容易遗漏固有频率,特别是当两个固有频率很接近时。文献[4-5]基于有限元方法分析了旋翼的固有振动特性,但这种方法主要适用于单体问题。文献[6]基于多体动力学程序建立了旋翼的结构动力学模型,文献[7]基于多体动力学方法计算了铰接式旋翼的固有频率,但是这种基于都是基于变分原理的完全积分方法,存在不可避免的剪切自锁现象,有限体积方法具有积分守恒的特点,能有效避免剪切自锁现象[8-9]。目前少有文献基于有限体积法对旋翼系统进行详细的固有振动特性分析。
基于有限体积法结合多体动力学方法,编写相关程序建立了旋翼的多体模型,其具有非线性负扭、后扫掠角并计及桨毂的扭转刚度,首先计算了旋翼在工作转速下的固有频率和振型,对比分析表明有限体积法能够对旋翼进行准确的模拟并展现此方法不存在剪切自锁的现象;同时还研究了转速、后扫掠结构和桨榖刚度对旋翼固有频率的影响。
有限部分梁的两端为点a和点b,梁上任意一个截面的内力ϑ为:
式中:F、M—梁截面的内力和内力矩;ϑa、ϑb—两端所受到的内力和内力矩记;—端点a到端点b的梁域,在梁域上有如下的平衡公式:
其中,力臂矩阵U定义如下式所示:
式中:ra、rb—梁的两个端点的位置坐标向量,上式中利用支点代替梁截面的剪切中心来简化计算平衡方程中的力矩项,r0是支点的位置坐标向量,支点一般取梁的两个端点。外部力定义如下:
其中,τ是广义上的外力和外力矩,其包含外部荷载τout和惯性力τin两部分,其中惯性力可表示如下:
式中:M—梁的惯性矩阵;r̈—线加速度;ω̈—角加速度。
桨叶模型采用几何精确梁模型[10],此模型为三维模型,对于梁的状态用四个变量(轴向变形u、摆振变形v、挥舞变形w和描述截面转动变形量的Gibbs-Rodriguez参数)描述。
采用三节点有限体积梁单元对梁模型进行离散,如图1 所示。每个三节点梁单元有两个控制中心节点Ⅰ和Ⅱ,中心节点采用2个高斯积分点使得多项式的积分达到了三阶精度,每个三节点有限体积梁单元中心节点具有四个变量(轴向变形u、摆振变形v、挥舞变形w和截面转动Gibbs-Rodriguez 参数g),有限体积梁模型的两个控制中心节点位于参考系的坐标和sII=,控制体中心节点之间的边界面可以通过两个中心节点的变量插值得到。
图1 三节点有限体积梁单元受内力和外力示意图Fig.1 The Three-Node Finite Volume Beam Element Under Internal and External Loads
存在如下所示的二次多项式组合函数(Nj(s)):
于是,梁单元上的任意一点的位置r(s)和方向参数g(s)都可以表示为节点1、2、3相应参数的线性组合。参考线上任意一点的位置为:
式中:Nj—第j个节点的形函数;rj—第j个节点的位置向量,上式表示的是以j为指针的求和运算,其中j=1,2,3。同理,任意一点的方向参数可表示为:
式中:gj—第j个节点的方向参数,当s=-1,0,1分别对应的是3个节点在局部坐标系中的位置。
离散后,三节点有限体积梁单元受力,如图1所示。该单元的离散方程包含三段梁的平衡方程,其公式如下所示:
式中:U—力臂矩阵;ϑI、ϑII—两个中心节点Ⅰ和Ⅱ的内力和内力矩;F—广义上的外力项。
列出所有离散单元的平衡方程,并进行方程组装,可得到一般矩阵形式的方程,其压缩形式如下:式中:M(q,t)—整体的质量矩阵;Ks(q,t)—结构刚度矩阵;自由度向量;—空气动力或者其它的外部力以及惯性力。
对于上述的振动方程,去掉外力项,即τout=0,其动力学变为如下的一般形式:
上式还可以进一步写成如下的简化形式:
式中:K(q,t)=Ks(q,t)+Ω2KCF(q,t),由梁的理论可假设位移响应满足如下关系:
将式(12)写成如下矩阵形式:
式中:ωn—系统的固有频率,ωn对应的解qn是系统的振型向量。考虑相应的边界条件以后,考虑相关位移约束方程后,求解上述频率方程,即可得到系统的固有频率和振型。例如对于固定约束,该处线位移和角位移均为0,而自由端的弯矩和剪力都为0。
基于有限体积C0梁单元建立了现代先进直升机UH-60A旋翼的多体模型,该旋翼具有许多先进直升机旋翼的几何特征,如非线性负扭转、后掠和不同的翼型组合等,桨叶的非线性扭转角分布,如图2所示。文献[7,11]给出了较为详细的结构和材料参数。该旋翼由四片桨叶组成,桨叶的展弦比15.3,后掠位置位于桨叶的93%R 处,后掠角为20°。桨叶由两种翼型构成,两段为SC1095翼型,中间部分为SC1095-R8翼型。桨叶的具体原始数据,如表1所示。
图2 UH-60A桨叶的扭转角分布Fig.2 The Torsion Angle Distribution of UH-60A Blade
表1 UH-60A旋翼的相关参数Tab.1 The Correlation Parameter of UH-60A Blade
桨叶模型采用三节点有限体积梁单元建模,一共才采用了23个梁单元,直梁部分具有19个梁单元,后扫掠结构具有4个梁单元,桨毂采用孤立的刚性节点模拟,如图3所示。支撑点采用固定约束,旋翼根部采用变距铰、挥舞铰和摆振铰约束,并在变距、挥舞和摆振方向上施加变形弹簧。为了考虑旋翼的变距系统的刚度对其固有频率的影响,将其等效为一常数刚度值,施加在旋翼根部的变距铰上,刚度值K0=84708N·m/rad[12]。桨叶的惯性通过集中质量单元来考虑,集中质量单元附加在梁单元的节点上。旋翼桨叶的后扫掠结构,采用真实拓扑建立模型,考虑其几何拓扑结构对其固有特性的影响。该模型为这里的第Ⅰ种模型。
图3 UH-60A旋翼的有限体积梁模型Fig.3 Finite Volume Beam Model of UH-60A
共建立三种不同的旋翼模型,以详细分析不同的边界条件、模型参数和拓扑结构对旋翼固有特性的影响,第Ⅱ种旋翼模型,除桨叶无后扫掠结构外,与模型Ⅰ完全相同。其中第Ⅲ种旋翼模型的桨叶根部采用摆振铰、挥舞铰,无扭转弹簧,即桨叶根部节点的扭转变形为零,其扭转自由度被完全锁死,如表2所示。
表2 三种模型的对比Tab2. The Comparison of Three Modals
在上述有限体积多体模型的基础上,基于非对称特征值QZ分解法计算了UH-60A旋翼在总距角为14.5°、工作转速状态下的扭转、挥舞和摆振的前几阶固有频率。这里的频率计算结果与参考值的对比[7],如表3所示。表中的频率值按工作转速无量纲化,并给出了相应模态的固有振型。
表3 计算结果与参考值的对比Tab.3 The Comparison of Computation Values and Reference Values
由上表可知,计算的频率值与参考值的误差均在工程计算误差要求范围以内,最大误差不超过7%,由此验证了这里的模型是可靠的。需要说明的是文献[7]以直梁模型对桨叶模型进行了等效化处理,也导致了这里计算值与参考值存在一定的差别。
计算所得的各阶模态对应的前五阶挥舞、前两阶摆振和前两阶扭转振型,如图4所示。由图可见,铰接式带后扫掠角的旋翼的一阶挥舞和一阶摆振振型呈直线式分布,由于该旋翼考虑了旋翼根部倾转系统的扭转刚度对其固有特性的影响,其扭转振型在根部的值具有一定的大小,其它各阶的计算结果也分别符合铰接式旋翼的分布规律。图标F、L、T分别表示挥舞、摆振和扭转模态。
图4 UH-60A旋翼的模态振型Fig.4 The Mode Shapes of UH-60A
计算了旋翼在(0~2)倍工作转速时的固有频率,以研究转速对旋翼固有频率的影响,计算结果,如图5所示。图中转速采用无量纲转速表示。由图可知,转速对旋翼的固有频率具有明显的“刚化效应”,各阶频率随着转速的增大而增大,其对挥舞和摆振频率的刚化效应比较明显,而对扭转频率的刚化效应不是很明显。
图5 转速对旋翼频率的影响Fig.5 The Effect of Rotor Speed on Rotor Frequency
为了进一步分析边界条件对旋翼固有频率的影响,对模型Ⅲ进行了固有特性计算,并与模型Ⅰ进行了对比分析,如表4所示。由表可知,当桨叶根部的扭转刚度很大时,采用两种模型计算结果很接近,当桨叶根部的扭转刚度不是特别大时,就必须考虑旋翼变距系统的扭转刚度对桨叶固有频率的影响。
图6 后扫掠结构对旋翼频率的影响Fig.6 The Effect of Back Swept Structure on the Frequency of Blade
表4 模型Ⅰ和模型Ⅱ的扭转频率计算结果的对比Tab.4 The Comparison of Torsion Frequency of ModalⅠand ModalⅡ
对于带桨毂弹簧的铰接式旋翼,桨叶的根部受到弹簧的弹性约束。当桨叶根部的摆振刚度K2=0~212×103N·m/rad时,频率计算结果如所示,横坐标为2的指数值。由图7(b)可知,与扭转频率的变化规律类似,前二阶摆振频率随着摆振刚度的增加而较快增加,且增加的幅度越来越小,最后趋于水平。当桨叶根部的挥舞刚度K3=0~212×103N·m/rad时,其频率计算结果,如图7(c)所示。横坐标为2的指数值。由图可知,挥舞频率的变化规律与扭转和摆振频率类似,前五阶挥舞频率随着根部挥舞刚度的增加而较快增加,且增加的幅度越来越慢,最后趋于水平。
图7 扭转、摆振和挥舞刚度对旋翼固有频率的影响Fig.7 The Effect of Torsion and Level and Flap Stiffness on Frequency of Blade
基于有限体积法,采用有限体积梁单元建立了旋翼的多体模型,基于非对称特征值QZ分解法计算了高速旋转运动条件下该旋翼的固有频率和振型,并分别分析了转速、后扫掠结构和桨叶根部约束刚度对其固有特性的影响。得出以下结论:
(1)基于有限体积方法能够准确计算旋翼的固有频率和振型,并可方便地计及桨叶根部变距系统的扭转刚度,误差在可接受范围之内;
(2)转速对旋翼的固有频率具有明显的刚化效应,尤其是对摆振和挥舞频率的刚化效应大,而对扭转频率的刚化效应小;
(3)桨叶的后扫掠结构对其摆振、挥舞频率影响很小,而后扫掠结构降低了桨叶的扭转频率,在旋翼的设计中应该适当加以考虑;
(4)桨叶根部的扭转、摆振和挥舞刚度分别对其扭转、摆振和挥舞振型的频率都有较大的影响,在计算其固有特性时,应根据实际结构适当加以考虑,不宜一概完全忽略不计,以免降低计算精度。