基于传递函数的机翼颤振速度确定方法技术

本发明专利技术涉及基于传递函数的机翼颤振速度确定方法，所述方法利用机翼弯扭振动微分方程和Theodrosen非定常气动力模型，得出长直机翼颤振微分方程，然后对所述机翼颤振微分方程进行Fourier变换，再运用传递函数法求解机翼颤振速度。本发明专利技术的有益效果如下：（1）由于本发明专利技术利用机翼弯扭振动微分方程来准确描述机翼振动，而没有采用机翼振动的低阶模态来近似描述机翼振动，所以计算的机翼颤振速度更加准确；（2）本发明专利技术求解机翼颤振速度的方法与现有技术相比更加简捷。

【技术实现步骤摘要】
基于传递函数的机翼颤振速度确定方法
本专利技术涉及一种基于传递函数的机翼颤振速度确定方法，属于飞行器气动弹性

技术介绍
机翼颤振是指发生在飞机飞行中的动不稳定性，此时的飞行速度称为颤振速度。在飞行达到颤振速度时所发生的自激振动，大多数都会造成灾难性的后果。颤振分析就是求出颤振发生的条件，亦即求出颤振速度，并寻求在飞行器飞行速度范围内避免颤振发生的措施，进而寻找提高颤振临界速度的方法。进行机翼颤振分析和计算时，一方面需要知道机翼结构的动力学特性，另一方面需要知道机翼结构附近非定常流动的空气动力特性。目前，工程中为了避免分析计算过于复杂，在机翼结构动力学计算时尽量减小机翼振动自由度，通常忽略高阶模态，仅利用机翼振动的低阶模态来近似表示机翼弯曲振动位移h和扭转振动转角α，通常选取前一阶、二阶弯曲模态和前一阶扭转模态来描述机翼的动力学特性，即然后，再结合Theodrosen非定常气动理论，求解机翼颤振速度。实际上，机翼是无限多自由度的弹性体，在颤振分析中也没有对自由度的限制。因此，上述计算方法的缺点是：为了避免过于复杂的计算而减小了机翼振动的自由度，对机翼实际振动作出了近似，使计算结果的准确性降低了。
技术实现思路
本专利技术所要解决的技术问题是提供一种计算方法简捷、计算结果较准确的基于传递函数的机翼颤振速度确定方法。本专利技术解决其技术问题所需采用的技术方案：一种基于传递函数的机翼颤振速度确定方法，所述方法利用机翼弯扭振动微分方程和Theodrosen非定常气动力模型，得出长直机翼颤振微分方程，然后对所述机翼颤振微分方程进行Fourier变换，再运用传递函数法求解机翼颤振速度；一、所述方法需要依据的计算公式：a.利用机翼弯扭振动微分方程建立机翼颤振微分方程：利用机翼弯扭振动微分方程式(1)来描述机翼弯扭振动：式中，h为机翼弯曲振动位移，单位为米，α为机翼扭转振动转角，单位为弧度，EI为机翼抗弯刚度，单位牛顿·米2，GJ为机翼抗扭刚度，单位为牛顿·米2，m为机翼单位长度质量，单位为千克，Iα为单位长度机翼绕弹性轴的转动惯量，单位为千克·米2，xα为机翼弹性轴到机翼截面重心的距离，单位为米，Lk为机翼单位长度的升力，单位为牛顿/米，Tα机翼单位长度的扭矩，单位为牛顿，y为机翼展向坐标，单位为米，t为时间，单位为秒，为机翼弯曲振动位移h对机翼展向坐标y的四阶偏导数，为机翼扭转振动转角α对机翼展向坐标y的二阶偏导数，为机翼弯曲振动位移h对时间t的二阶偏导数，为机翼扭转振动转角α对时间t的二阶偏导数；b.根据Theodrosen非定常气动力模型，得到机翼单位长度的升力Lk和机翼单位长度的扭矩Tα为下式(2)：式中，V为空速，单位为米/秒，ρ为空气密度，单位为千克/立方米，b为机翼的半弦长，单位为米，C(k)为Theodrosen函数，为减缩频率，ω为圆频率，单位为弧度/秒，为机翼弹性轴到机翼弦长中点的距离占半弦长的百分比；由于减缩频率k是圆频率ω和空速V的函数，将C(k)写为C(ω，V)，将(2)式代入(1)式得到机翼颤振微分方程式(3)：c.利用传递函数法求解机翼颤振速度①对(3)式作Fourier变换，并经过整理可得到下式(4)：式中，A1(ω，V)、A2(ω，V)、B1(ω，V)、B2(ω，V)的具体表达式如下：其中，虚数②采用传递函数法求解上述(4)式，确定机翼的颤振速度；为了便于应用传递函数理论，定义状态变量向量如下式(6)：式中，T表示向量转置；从而，将(4)式写成状态方程的形式如下式(7)：式中，g(y，ω)＝0边界条件为：式中，η(0，ω)为η(y，ω)在机翼根端y＝0处的值，为η(y，ω)在机翼梢端处的值，为机翼的半展长，单位为米，Mb为机翼根端边界条件选择矩阵，Nb为机翼梢端边界条件选择矩阵，γ(ω)为由边界条件给定的位移或力组成的向量，其表达式分别为：根据传递函数理论，方程(7)式的解为：式中，G(y，ξ，ω，V)为状态空间方程的域内传递函数，H(y，ω，V)为状态空间方程的边界传递函数，其表达式分别为：式中，变量为机翼展向的坐标，单位为米：机翼颤振时弯曲振动位移h和扭转振动转角α的振幅为非零常数，即(10)式有非零解，根据传递函数理论可知，(10)式有非零解的充分必要条件为：令由于A为复矩阵，其行列式值等于零的必要条件为矩阵行列式值的实部与虚部均为零，即求解上述(13)式，可以得到满足方程组的空速V和圆频率ω，分别记为Vcz和ωcz。其中，Vcz即为机翼的颤振速度，ωcz即为机翼的颤振圆频率。二、所述方法的具体步骤如下：步骤(一)：测量长直机翼的下列物理参数：半弦长b，单位为米；半展长单位为米；单位长度机翼质量m，单位为千克/米；机翼弹性轴(2)到机翼弦长中点的距离占半弦长的百分比机翼弹性轴(2)到机翼截面重心的距离单位为米；单位长度机翼绕弹性轴的转动惯量单位为千克·米2；机翼的抗弯刚度EI，单位为牛顿·米2；机翼的抗扭刚度GJ，单位为牛顿·米2；空气密度ρ，单位为千克/立方米；步骤(二)：确定飞机机翼的空速V和圆频率ω的大致范围，假设为步骤(三)：在范围内划分合适步长ΔV和Δω，并进行离散，空速V和圆频率ω的取值为：步骤(四)：取空速V＝V0，圆频率ω依次取ω0+jΔω，j＝0，1，2，3…，将空速V和圆频率ω的取值与步骤(一)中的机翼各物理参数代入(5)式，得到系数A1(ω，V)、A2(ω，V)、B1(ω，V)、B2(ω，V)的值；步骤(五)：将系数A1(ω，V)、A2(ω，V)、B1(ω，V)、B2(ω，V)的值代入(8)式，得到矩阵F(ω，V)的值；步骤(六)：将矩阵F(ω，V)的值代入(13)式，计算Re[detA]和Im[detA]的值；步骤(七)：再依次取空速V＝V0+jΔV，j＝1，2，3…，重复步骤(四)至步骤(六)，计算Re[detA]和Im[detA]的值；步骤(八)：确定同时满足(13)式的空速V的值，即为机翼的颤振速度Vcz。本专利技术的有益效果如下：(1)由于本专利技术利用机翼弯扭振动微分方程来准确描述机翼振动，而没有采用机翼振动的低阶模态来近似描述机翼振动，所以计算的机翼颤振速度更加准确。(2)本专利技术求解机翼颤振速度的方法与现有技术相比更加简捷。附图说明图1为长直机翼示意图；图2为长直机翼弦向剖面示意图；图3为实施例1的Re[detA]的等值线图(V∈(0，50)，ω∈(0，200π))；图4为实施例1的Im[detA]的等值线图(V∈(0，50)，ω∈(0，200π))；图5为图3和图4的等值线合成图(V∈(0，50)，ω∈(0，200π))；图6为实施例1的Re[detA]的等值线图(V∈(30，40)，ω∈(30π，60π))；图7为实施例1的Im[detA]的等值线图(V∈(30，40)，ω∈(30π，60π))；图8为图6和图7的等值线合成图(V∈(30，40)，ω∈(30π，60π))。在图1、2中，1——重心轴，2——弹性轴，3——弹性轴位置，4——重心轴位置。具体实施方式实施例1：为了进一步说明本专利技术所述方法，本实施例1与参考文献(赵永辉的专著《气动弹性力学与控制》，北京：科学出版社，2006)的计算结果(Vcz＝36.6米/秒)相比较。本实施例1本文档来自技高网...

【技术保护点】
一种基于传递函数的机翼颤振速度确定方法，其特征在于所述方法利用机翼弯扭振动微分方程和Theodrosen非定常气动力模型，得出长直机翼颤振微分方程，然后对所述机翼颤振微分方程进行Fourier变换，再运用传递函数法求解机翼颤振速度；一、所述方法需要依据的计算公式：a.利用机翼弯扭振动微分方程建立机翼颤振微分方程：本专利技术利用机翼弯扭振动微分方程式(1)来描述机翼弯扭振动：EI∂4h∂y4+m∂2h∂t2-mxα∂2α∂t2-Lh=0GJ∂2α∂y2-Iα∂2α∂t2+mxα∂2h∂t2+Tα=0---(1)]]>式中，h为机翼弯曲振动位移，单位为米，α为机翼扭转振动转角，单位为弧度，EI为机翼抗弯刚度，单位牛顿·米2，GJ为机翼抗扭刚度，单位为牛顿·米2，m为机翼单位长度质量，单位为千克，Iα为单位长度机翼绕弹性轴的转动惯量，单位为千克·米2，xα为机翼弹性轴到机翼截面重心的距离，单位为米，Lh为机翼单位长度的升力，单位为牛顿/米，Tα机翼单位长度的扭矩，单位为牛顿，y为机翼展向坐标，单位为米，t为时间，单位为秒，为机翼弯曲振动位移h对机翼展向坐标y的四阶偏导数，为机翼扭转振动转角α对机翼展向坐标y的二阶偏导数，为机翼弯曲振动位移h对时间t的二阶偏导数，为机翼扭转振动转角α对时间t的二阶偏导数；b.根据Theodrosen非定常气动力模型，得到机翼单位长度的升力Lh和机翼单位长度的扭矩Tα为下式(2)：Lh=πρb2(-∂2h∂t2+V∂α∂t-ba‾∂2α∂t2)+2πρVbC(k)(Vα-∂h∂t+b(12-a‾)∂α∂t)Tα=πρb2(-ba‾∂2h∂t2-Vb(12-a‾)∂α∂t-b2(18+a‾2)∂2α∂t2)+2πρVb2(12+a‾)C(k)(Vα-∂h∂t+b(12-a‾)∂α∂t)---(2)]]>式中，V为空速，单位为米/秒，ρ为空气密度，单位为千克/立方米，b为机翼的半弦长，单位为米，C(k)为Theodrosen函数，为减缩频率，ω为圆频率，单位为弧度/秒，为机翼弹性轴到机翼弦长中点的距离占半弦长的百分比；由于减缩频率k是圆频率ω和空速V的函数，将C(k)写为C(ω,V)，将(2)式代入(1)式得到机翼颤振微分方程式(3)：EI∂4h∂y4+m∂2h∂t2-mxα∂2α∂t2-πρb2(-∂2h∂t2+V∂α∂t-ba‾∂2α∂t2)-2πρVbC(ω,V)(Vα-∂h∂t+b(12-a‾)∂α∂t)=0GJ∂2h∂y2-Iα∂2α∂t2+mxα∂2h∂t2+πρb2(-ba‾∂2h∂t2-Vb(12-a‾)∂α∂t-b2(18+a‾2)∂2α∂t2)+2πρVb2(12+a‾)C(ω,V)(Vα-∂h∂t+b(12-a‾)∂α∂t)=0---(3)]]>c.利用传递函数法求解机翼颤振速度①对(3)式作Fourier变换，并经过整理可得到下式(4)：∂4α∂y4=A1(ω,V)h+B1(ω,V)α∂2α∂y2=A2(ω,V)h+B2(ω,V)α...

【技术特征摘要】
1.一种基于传递函数的机翼颤振速度确定方法，其特征在于所述方法利用机翼弯扭振动微分方程和Theodrosen非定常气动力模型，得出长直机翼颤振微分方程，然后对所述机翼颤振微分方程进行Fourier变换，再运用传递函数法求解机翼颤振速度；一、所述方法需要依据的计算公式：a.利用机翼弯扭振动微分方程建立机翼颤振微分方程：利用机翼弯扭振动微分方程式(1)来描述机翼弯扭振动：式中，h为机翼弯曲振动位移，单位为米，α为机翼扭转振动转角，单位为弧度，EI为机翼抗弯刚度，单位牛顿·米2，GJ为机翼抗扭刚度，单位为牛顿·米2，m为机翼单位长度质量，单位为千克，Iα为单位长度机翼绕弹性轴的转动惯量，单位为千克·米2，xα为机翼弹性轴到机翼截面重心的距离，单位为米，Lh为机翼单位长度的升力，单位为牛顿/米，Tα机翼单位长度的扭矩，单位为牛顿，y为机翼展向坐标，单位为米，t为时间，单位为秒，为机翼弯曲振动位移h对机翼展向坐标y的四阶偏导数，为机翼扭转振动转角α对机翼展向坐标y的二阶偏导数，为机翼弯曲振动位移h对时间t的二阶偏导数，为机翼扭转振动转角α对时间t的二阶偏导数；b.根据Theodrosen非定常气动力模型，得到机翼单位长度的升力Lh和机翼单位长度的扭矩Tα为下式(2)：式中，V为空速，单位为米/秒，ρ为空气密度，单位为千克/立方米，b为机翼的半弦长，单位为米，C(k)为Theodrosen函数，为减缩频率，ω为圆频率，单位为弧度/秒，为机翼弹性轴到机翼弦长中点的距离占半弦长的百分比；由于减缩频率k是圆频率ω和空速V的函数，将C(k)写为C(ω，V)，将(2)式代入(1)式得到机翼颤振微分方程式(3)：c.利用传递函数法求解机翼颤振速度①对(3)式作Fourier变换，并经过整理可得到下式(4)：式中，A1(ω，V)、A2(ω，V)、B1(ω，V)、B2(ω，V)的具体表达式如下：

【专利技术属性】
技术研发人员：段静波，江涛，路平，蔚建斌，邱金刚，左军，
申请(专利权)人：中国人民解放军军械工程学院，
类型：发明
国别省市：河北;13