正弦突风下仿生机翼颤振分析

张宴嘉，王新阁，贾文隽

(1.空军航空大学 航空作战勤务学院， 长春 130022； 2.中国人民解放军93286部队， 沈阳 110000)

1 引言

高空长航时无人机多采用大展弦比柔性机翼，来提高机翼的气动性能，满足其长航时的作战要求。但是大展弦比柔性机翼在气动力的作用下会产生较大的弯曲扭转变形[1-2]，当速度达到一定条件时，会发生颤振。经研究发现[3]，仿生机翼不但会提高气动效率，还能有效地抑制颤振的发生。

与以往仿真实验所使用的稳定来流不同，实际大气环境中存在着湍流气体，其分布和运动情况复杂多变，从而现实中相对于飞行器的来流速度不再稳定[4]。在风速不稳定的情况下，要引入突风模型来模拟流场的状态。突风是随时间变化的函数，这里采用工程上最常用的正弦sin突风模型来模拟大气流场，正弦突风模型曲线如图1。

图1 正弦突风模型曲线Fig.1 Sinusoidal gust model

对正弦突风模型进行数值设置，突风的表达式为V=0.01*sin(6.28*2.463 5*t)+0.1，突风的速度值在0.1马赫上下浮动，上下浮动大小为0.01马赫，正弦函数的频率为2.463 5 Hz，对应着仿生机翼的一阶频率。将速度编写为程序如下：

将程序文件以udf的形式导入到fluent中压力远场边界条件中，利用正弦突风模型进行双向流固耦合计算。

2 流固耦合基本理论

按照数据的传递方式，可分为单向流固耦合和双向流固耦合。单向流固耦合是经过流体分析模块计算后把结果数据传递给固体分析模块，然后经固体模块结构分析后得到结果，耦合即结束。用于研究流体载荷对结构体的单向作用问题。双向流固耦合，是在单向流固耦合的基础上，增加了固体变形量数据返回到流体分析模块，再进行流场分析，之后重复该过程，经多次循环后完成耦合分析。研究流体与结构体之间双向作用的效应及影响[5]。

流体分析过程中要遵循质量守恒定律、动量守恒定律和能量守恒定律，三大物理守恒定律：

连续方程为：

(1)

式中：μ、ν和ω代表速度矢量在x、y和z方向的分量。

动量方程为：

(2)

其中：ρ代表密度；u为流体运动速度;F为外力;p为压力;μ0为动力黏性系数。

能量方程为：

(3)

依据弹性力学基本原理，机翼结构力学建立在 3个基本方程上：平衡微分方程、几何变形方程以及应力-应变方程(本构关系)。受到流场与自重作用，机翼发生位移变化，基于有限元法构建结构动力学控制方程为：

(4)

3 计算模型的建立

3.1 CSD模型

1) 模型的建立与导入

大展弦比仿生海鸥机翼设计分为两段[6]。从靠近机体处到半展长2/3的部分，以海鸥翅翼模型50%翼展处的翼型为基础，沿着海鸥翅翼模型的前缘线按比例放样而得。余下的1/3半展长部分为一段平直机翼。目的是为了增大展弦比和减小翼尖涡。仿生机翼模型如图2所示。仿生机翼半展长7.05 m，翼根弦长1.2 m，翼尖弦长0.54 m，展弦比λ=14.8，参考面积6.72 m2，后掠角0°。

图2 仿生机翼模型示意图Fig.2 Bionic wing model

2) 机翼网格划分

使用mesh模块进行网格划分。机翼表面采用六面体网格进行划分，插入网格尺寸为体尺寸，大小为0.04 m，在翼根和翼尖表面插入面尺寸，大小为0.02 m，进一步进行加密。网格在机翼的圆弧处进行加密，设定圆弧角5°。机翼网格划分结果如图3所示，共62 826个网格单元，平均扭曲度为0.22。

图3 仿生机翼网格模型示意图Fig.3 Mesh model of bionic wing

3.2 CFD模型

1) 建立计算域

为了完成对机翼周围流场的分析，需要先建立计算域，建立的流体计算域如图4所示。

图4 仿生机翼计算域模型示意图Fig.4 Calculation domain model of bionic wing

2) 网格划分

划分网格时，对于远场部分，流场变化较为平稳，适当放宽网格尺寸；机翼表面附近流场，变化较为剧烈，划分网格进行加密处理，重点关注近壁面处网格的尺寸和质量。网格划分如图5，共有网格单元数量23 007，节点数量963 742。

图5 网格划分示意图Fig.5 grid generation

4 仿生机翼颤振的发生

4.1 计算结果分析

材料属性选择Al，攻角为0°，计算得翼尖沿y方向的位移如图6所示。

图6 攻角0°时，Al的翼尖位移曲线Fig.6 Angle of attack 0 ° Tip displacement response of Al

将程序文件以udf的形式导入到fluent中压力远场边界条件中，利用正弦突风模型进行双向流固耦合计算。

由图6所示，正弦突风环境下，机翼翼尖位移振荡的幅值提升，机翼的振动呈发散趋势，翼尖最大位移达到了0.053 2 m，机翼发生了颤振。对翼尖位移进行频谱分析[7]，分析其频域下的特征如图7所示。在频谱图中共出现3个峰值，对应频率分别为：2.506 Hz、47.62 Hz和52.63 Hz，其中2.506 Hz为主峰，与仿生机翼的一阶频率接近，47.62 Hz与6阶频率接近。说明在突风作用下，仿生机翼的主要为一阶弯曲颤振，但也夹杂着少量的高阶振型，机翼的振荡和形变具有多样化[8]。

图7 攻角0°时，Al的翼尖位移响应频域分析曲线Fig.7 Angle of attack 0 ° Frequency domain analysis of tip displacement response of Al

通过能量法对仿生机翼的颤振问题进行解释[9]。根据图6的位移响应曲线，取第3个振荡周期进行分析，时间为0.88～1.26 s，总共20个时间步，故t=20，所要分析的时间点在整个计算时长所处的位置如图8蓝点所示。机翼仅在表面上受到气动力的作用，因此仅考虑机翼表面节点即可，表面节点n=332。对0.88～1.26 s上每个时间节点的气动力做功进行计算，得到第3个周期中对机翼做的总功为Wair=3.668 7e5>0，说明机翼是从大气中吸收能量的，因此振动趋于不稳定，位移曲线发散。

图8 周期内时间点的采集曲线Fig.8 Collection of time points in a cycle

综上所述，仿生机翼在正弦突风V=0.01*sin(6.28*2.463 5*t)+0.1马赫的作用下振动发散，出现了颤振现象[10]。

4.2 与稳定风速结果的对比

在正弦突风的环境下，仿生机翼发生了颤振，而稳定风速的条件下振荡却是收敛的[11]，下面将稳定风速下的单向耦合、双向耦合结果与正弦突风下的结果进行对比，探究之间具体的差异。图9为位移响应、升力系数和阻力系数曲线。

图9 稳定风速与正弦突风曲线Fig.9 Comparison of steady wind speed and sinusoidal gust

由图9(a)～图9(c)可以看出，稳定风速的结果与正弦突风的计算结果有着很大的不同。在正弦突风的作用下，仿生机翼的升力系数和阻力系数都保持着较大的振幅，振幅随时间不衰减，升力系数和阻力系数随时间的变化很大，各周期下其最大峰值始终高于稳速下双向耦合的结果，最小峰值低于稳速下双向耦合的结果。说明仿生机翼在正弦突风的作用下处于十分不稳定的状态，正弦突风使机翼出现不稳定发散的位移，因此会造成机翼气动力系数的急剧不稳定。

稳定风速下，双向耦合的计算结果在达到设定的终止计算时间时基本趋于稳定，与单向耦合结果差距不大[12]。而正弦突风的作用下，机翼的气动力和位移结果波动极大，没有趋于收敛，而是保持大幅度振荡状态。因此在计算中要考虑到实际大气情况，不能仅仅将流场考虑为稳定风速加以计算。

4.3 影响仿生机翼颤振特性的因素

为了进一步对正弦突风下仿生机翼的颤振特性进行分析，在计算方法不变的前提下，改变某一特征量来分析该量对机翼颤振的影响。分析材料属性、攻角、速度幅值以及正弦突风频率对颤振的影响。

4.3.1材料属性对颤振的影响

选择ANSYS材料库中常用的3种材料：Al、Fe、Ti赋予给仿生机翼，分别在正弦突风下进行流固耦合计算，攻角均为0°，正弦突风V=0.01*sin(6.28*2.463 5*t)+0.1马赫。

图10(a)～图10(c)为3种不同材料属性在算后得到的y方向位移响应曲线。变形关系为：Al>Ti >Fe，Al的最大峰值为0.053 2 m，Fe的最大峰值为0.024 8 m，Ti的最大峰值为0.047 3 m。从图中可以发现，3种材料属性的机翼在正弦突风的作用下位移响应均为发散，由此判定都发生了颤振。

图10 不同材料在突风下的位移响应曲线Fig. 10 Displacement response of different materials under gust

图11(a)～图11(c)为Al、Fe、Ti的位移响应曲线频域下的分析，3种材料在频域下均只出现一个主峰，分别为：2.506、2.339、2.172 Hz，通过对模态分析结果的分析比对，发现它们的主频都在各自的一阶频率附近，说明3种材料的机翼在突风的作用下都以一阶弯曲变形为主。

图11 3种材料频谱分析Fig.11 Spectrum analysis of three materials

综上所述，在正弦突风的作用下，不同材料属性的机翼由于其弹性模量的差异从而造成位移响应不同；但是3种材料机翼的位移响应都呈发散趋势，位移保持大幅度振荡状态，根据频谱分析，3种材料机翼在颤振时均以一阶弯曲变形为主。因此，材料属性的不同对颤振发生与否不起决定性作用。

4.3.2攻角对颤振的影响

探究攻角对仿生机翼颤振的影响选择材料为Al，来流速度为V=0.01*sin(6.28*2.463 5*t)+0.1马赫，攻角取0°、4°、8°、12°、16°进行流固耦合计算，计算得不同攻角下翼尖的y方向位移如图12所示，随后对各攻角下的位移响应进行频谱分析，得到频域图如图13所示。

图12 不同攻角下位移曲线Fig.12 Displacement response at different angles of attack

由图12可知，随着攻角的增大，翼尖的位移响应不断增大，在所有攻角下翼尖 位移都呈现发散状态，说明0°～16°攻角下机翼都发生了颤振。由图13的频谱分析可知，伴随着攻角的增加，位移响应在频域下主峰的幅度也逐步提高，但它们的主峰所对应的频率都为2.506 Hz，说明在0°～12°攻角下，机翼的形变均为一阶弯曲变形。综合看来，在一定攻角范围内，正弦突风不变的情况下，攻角的增加会使机翼的振动位移增大，但对机翼是否发生颤振和颤振的形变姿态影响不大。

图13 不同攻角下位移频谱Fig.13 Spectrum analysis of displacement at different angles of attack

4.3.3速度幅值对颤振的影响

探究突风的幅值对仿生机翼颤振特性的影响，分别取3种不同的突风速度分别为：

V1=0.005*sin(6.28*2.463 5*t)+0.05

V2=0.01*sin(6.28*2.463 5*t)+0.1

V3=0.015*sin(6.28*2.463 5*t)+0.15

其中速度的幅值关系为：V1

由图14可知，随着突风速度幅值的增大，翼尖位移响应也随之增大，对于3种不同的突风速度翼尖位移都呈现发散状态，3种工况下机翼均发生了颤振。由图15的频谱分析可知，伴随着速度幅值的增加，位移响应在频域下主峰的幅度也逐步提高，但它们的主峰所对应的频率都为2.506 Hz，说明在V1、V2、V3三种突风环境下，机翼的形变均为一阶弯曲变形。综上所述，正弦突风速度幅值的提高会使机翼的位移增加，但是速度的改变不会影响到机翼形变的姿态，在一定速度幅值范围内对其是否发生颤振影响不大。

图14 不同速度突风下位移响应曲线Fig.14 Under different speed gust displacement response

图15 不同速度突风下位移频谱Fig.15 Frequency spectrum analysis of displacement under different speed gusts

4.3.4突风频率对颤振的影响

在之前的计算中，突风速度的频率设置的是仿生机翼一阶固有频率，即2.463 5 Hz，发现在改进机翼材料属性、攻角以及速度幅值的情况下均发生了颤振，对应的形变为一阶弯曲。机翼的材料选择Al，计算攻角为0°，对突风的频率进行修改，首先将仿生机翼的一阶固有频率作为基准，以±2%的幅度修改突风频率，改变后的突风速度分别为：

Va=0.01*sin(6.28*2.41*t)+0.1

Vb=0.01*sin(6.28*2.51*t)+0.1

其中：Va的输入频率的2.41 Hz，为98%一阶固有频率大小；Vb的输入频率为2.51 Hz，为102%一阶固有频率大小。在Va和Vb的速度下，翼尖沿y方向的位移响应曲线和对应的频谱如图16～图17所示。

由图16～图17可知，在2.41 Hz和2.51 Hz的速度输入频率下，翼尖位移都呈发散状态，在经历了8个周期的振荡后，位移曲线基本保持等幅振荡，机翼发生了颤振。在Va的速度下，翼尖的最大位移为0.050 505 m，在Vb的速度下，翼尖的最大位移为0.049 917 m，两者的差距很小。根据5.13～18的频谱分析可知，位移响应在频域下的峰值分别为2.425 Hz和2.506 Hz，分别接近各自的速度输入频率，可见突风速度的输入频率影响着机翼振动的频率。

图16 输入频率为2.41 Hz时的突风速度(Va)时的位移曲线和频谱Fig.16 Gust velocity(Va) at 2.41 Hz input frequency

图17 输入频率为2.51 Hz时的突风速度(Vb)时的位移曲线和频谱Fig.17 Gust velocity(Vb) at 2.51 Hz input frequency

图18～图19为以±5%的幅度修改突风频率，计算所得的翼尖位移响应和频谱分析。改变后的突风速度分别为：

图18 输入频率为2.34 Hz时的突风速度(Vc )时的位移曲线和频谱Fig.18 Gust velocity(Vc) at 2.34 Hz input frequency

图19 输入频率为2.58 Hz时的突风速度(Vd )时的位移曲线和频谱Fig.19 Gust velocity(Vd) at 2.58 Hz input frequency

Vc=0.01*sin(6.28*2.34*t)+0.1

Vd=0.01*sin(6.28*2.58*t)+0.1

由图18～图19可知，在2.34 Hz和2.58 Hz的速度输入频率下，翼尖位移曲线先发散再收敛最终保持等幅振荡状态，机翼发生了颤振。在Vc的速度下，翼尖的最大位移为0.045 392 m，在Vd的速度下，翼尖的最大位移为0.046 073 m，两者的差距很小。根据图18～图19的频谱分析可知，位移响应在频域下的峰值分别为2.339 Hz和2.597 Hz，分别接近各自的速度输入频率，可见突风速度的输入频率影响着机翼振动的频率。

图20～图21为以±10%的幅度修改突风频率，计算所得的翼尖位移响应和频谱分析。改变后的突风速度分别为：

Ve=0.01*sin(6.28*2.21*t)+0.1

Vf=0.01*sin(6.28*2.71*t)+0.1

由图20～图21可知，在2.21 Hz和2.71 Hz的速度输入频率下，翼尖位移曲线先收敛再发散最终保持等幅振荡状态，机翼发生了颤振。在Ve的速度下，翼尖的最大位移为0.043 244 m，在Vf的速度下，翼尖的最大位移为0.043 067 m，两者的差距很小。根据图20～图21的频谱分析可知，位移响应在频域下的峰值分别为2.172 Hz和2.673 Hz，分别接近各自的速度输入频率，可见突风速度的输入频率影响着机翼振动的频率。

图20 输入频率为2.21 Hz时的突风速度(Ve )时的位移曲线和频谱Fig.20 Gust velocity(Ve) at 2.21 Hz input frequency

图21 输入频率为2.71 Hz时的突风速度(Vf )时的位移曲线和频谱Fig.21 Gust velocity(Vf) at 2.71 Hz input frequency

综合分析图16～图21可知，当突风速度的输入频率在仿生机翼一阶固有频率±10%内波动时，机翼的位移曲线最终会保持大幅度等幅振动，机翼都会出现颤振现象。对于机翼振荡的幅值，对比发现突风速度的输入频率越接近一阶固有频率，振荡幅值越大，机翼的颤振越剧烈。在突风的作用下，翼尖位移响应的频率接近各自的突风速度输入频率，可见机翼的振动频率由外界的输入频率决定。

图22为以+20%的幅度修改突风频率，计算所得的翼尖位移响应和频谱分析。改变后的突风速度为：

Vg=0.01*sin(6.28*2.95*t)+0.1

由图22可知，在2.95 Hz的速度输入频率下，翼尖位移曲线呈收敛趋势，机翼没有发生颤振。在Vg的速度下，翼尖的最大位移为0.037 711 m。根据图22(b)的频谱分析可知，位移响应在频域下的峰值为2.506 Hz，接近机翼的一阶固有频率。

图22 输入频率为2.95 Hz时的突风速度(Vg )时的位移曲线和频谱Fig.22 Gust velocity(Vg) at 2.95 Hz input frequency

继续增大突风的修改幅度，分别采用+200%、+300%一阶固有频率和二阶固有频率作为突风的速度输入频率，计算所得的翼尖位移响应和频谱分析如图24～图26所示。改变后的突风速度分别为：

图24 输入频率为7.38 Hz时的突风速度(Vi )时的位移曲线和频谱Fig.24 Gust velocity(Vi) at 7.38 Hz input frequency

Vh=0.01*sin(6.28*4.92*t)+0.1

Vi=0.01*sin(6.28*7.38*t)+0.1

Vj=0.01*sin(6.28*10.2*t)+0.1

图23 输入频率为4.92 Hz时的突风速度(Vh )时的位移曲线和频谱Fig.23 Gust velocity(Vh) at 4.92 Hz input frequency

在Vh的速度下，翼尖的最大位移为0.038 544 m，在Vi的速度下，翼尖的最大位移为0.038 179 m，在Vj的速度下，翼尖的最大位移为0.038 032 m。根据图24的频谱分析可知，位移响应在频域下的峰值均为2.506 Hz。

综合对比图22～图25分析可知，当输入频率在120%一阶固有频率及以上时，翼尖位移曲线均呈收敛趋势，机翼都没有发生颤振。不同输入频率下，翼尖的最大位移差距很小，说明机翼的位移在收敛状态下，输入频率的大小对机翼振荡的幅值影响不大。翼尖位移响应的频率均接近一阶固有频率，说明当输入频率在120%一阶固有频率及以上时，输入频率的改变不会影响机翼的振动频率，机翼的形变以一阶弯曲变形为主。

图25 输入频率为10.2 Hz时的突风速度(Vj )时的位移曲线和频谱Fig.25 Gust velocity(Vj) at 10.2 Hz input frequency

图26为以-20%的幅度修改突风频率，计算所得的翼尖位移响应和频谱分析。改变后的突风速度分别为：

Vk=0.01*sin(6.28*1.97*t)+0.1

由图26可知，在1.97 Hz的速度输入频率下，翼尖位移曲经历了数次收敛和发散的过程，呈现不规律振荡，机翼的振动情况复杂。根据图26的频谱分析可知，位移响应在频域下出现2个峰值，主峰为2.005 Hz，次峰为2.506 Hz接近机翼的一阶固有频率，说明在80%一阶固有频率的输入频率下，机翼的振动频率主要为外界输入频率，以一阶固有频率振动的幅度较弱。

图26 输入频率为1.97 Hz时的突风速度(Vk )时的位移曲线和频谱Fig.26 Gust velocity(Vk) at 1.97 Hz input frequency

继续增大突风的修改幅度，分别采用70%、50%、20%一阶固有频率作为突风的速度输入频率，计算所得的翼尖位移响应和频谱分析如图27～图29所示。改变后的突风速度分别为：

图27 输入频率为1.72 Hz时的突风速度(Vl )时的位移曲线和频谱Fig.27 Gust velocity(Vl) at 1.72 Hz input frequency

Vl=0.01*sin(6.28*1.72*t)+0.1

Vm=0.01*sin(6.28*1.23*t)+0.1

Vn=0.01*sin(6.28*0.492*t)+0.1

图28 输入频率为1.23 Hz时的突风速度(Vm )时的位移曲线和频谱Fig.28 Gust velocity(Vm) at 1.23 Hz input frequency

图29 输入频率为0.492 Hz时的突风速度(Vn )时的位移曲线和频谱Fig.29 Gind burst speed(Vn) at input frequency of 0.492 Hz

综合对比图26～图29分析可知，当输入频率在80%一阶固有频率及以下时，翼尖位移曲线均经历了数次收敛和发散的过程，呈现不规律振荡，机翼的振动情况复杂，在这种无规则振动下极易造成破坏，影响安全性能。根据频谱分析可知，位移响应在频域下均出现两个峰值，分别对应的是各自的外界输入频率和一阶固有频率，当输入频率较大时，机翼以输入频率振动为主，随着输入频率的减小，机翼的振动逐渐转为以一阶固有频率振动为主导，说明输入频率的大小对机翼的振动周期有一定影响。

综上所述，对正弦突风的输入频率进行了改变，将20%一阶固有频率到二阶固有频率作为突风速度的输入频率进行了流固耦合计算，得到的规律如表1所示。

表1 输入频率与机翼振动的关系

当正弦突风速度的输入频率大于0小于机翼一阶固有频率时，机翼振动时而收敛时而发散，无规则混乱状态。机翼以一阶固有频率和输入频率振动，输入频率越大，以输入频率振动越明显、以一阶固有频率振动较弱，反之以一阶固有频率振动更明显、以输入频率振动较弱。

当正弦突风速度的输入频率大于等于90%一阶固有频率小于等于110%一阶固有频率时，机翼的位移发散，最终保持等幅振动，发生了颤振现象。输入频率越接近一阶固有频率，振动的幅值越大，反之振动的幅值越小。机翼其以外界输入频率发生振动。

当正弦突风速度的输入频率大于110%一阶固有频率小于等于二阶固有频率时，机翼的振动呈收敛状态，没有发生颤振。输入频率的大小对机翼的振动幅值影响不大，机翼其以一阶固有频率发生振动。

因此，仿生机翼在飞行过程中要尽量避免出现在正弦突风的环境下，尤其突风的速度输入频率在0～110%一阶固有频率之间，机翼会随时间发生较大的弯扭变形甚至出现结构破坏。当输入频率超过110%一阶固有频率时，变形逐渐收敛，飞行较为安全。

4.3.5平直机翼颤振分析

本节引入大展弦比平直机翼模型，应用双向流固耦合的计算方法，探究不同的输入频率对其颤振特性的影响。截面翼型为NACA0417翼型，为大展弦比机翼。平直机翼半展长7.05 m，翼根弦长1.2 m，翼尖弦长0.54 m，展弦比λ=16.2，参考面积6.05 m2，后掠角0°。平直机翼的材料选择Al，攻角为0°经过计算后总结出平直机翼的颤振特性随输入频率的变化关系，如表2所示。仿生机翼在外界输入频率为 90%～110%一阶固有频率的情况下发生颤振，颤振的范围小于平直机翼的颤振范围。仿生机翼在外界输入频率为110%一阶固有频率至二阶固有频率的情况下振动趋于收敛，振动收敛的频率范围大于平直机翼的频率范围。综合对比颤振特性可知，仿生机翼优于平直机翼，仿生机翼的构型设计有效地抑制了颤振的发生。

表2 输入频率与平直机翼振动的关系

5 结论

1) 在正弦突风的作用下，仿生机翼的位移、升力系数和阻力系数与稳定风速的结果有着很大的差距，说明仿生机翼在突风作用下处于十分不稳定的状态。

2) 从仿真结果分析可得，材料属性、攻角和速度幅值只会对仿生机翼位移响应的幅值造成影响，而对振动频率和颤振的影响不大，机翼均以一阶弯曲变形为主。而在不同的输入频率下，机翼振动形式和振动频率差异极大，因此，正弦突风的输入频率是颤振分析的重点。

3) 通过分析正弦突风的输入频率对大展弦比平直机翼颤振特性的影响可知，平直机翼的颤振的频率范围更大，仿生机翼在颤振方面的性能优于平直机翼。