仿鸟柔性扑翼气动特性与能耗的数值研究

仿鸟柔性扑翼气动特性与能耗的数值研究肖天航，段文博，昂海松( 南京航空航天大学航空宇航学院，江苏 南京 210016)摘 要:建立了适当的三维仿鸟柔性扑翼模型，并以配平重力和平衡阻力为条件，数值计算了它的低雷诺数非定常流场。

研究揭示了翼面初始扭转角度、动态俯仰幅度等重要设计参数与飞行性能的关系，表明扑翼平面的初始扭转 程度、扑翼柔性材料的选择以及两者之间的合理搭配对扑翼机的成功飞行至关重要。

研究分析了仿鸟扑翼的流场 涡结构、升力推力产生原理，下扑过程附着上翼面的前缘涡是升力产生的重要机制。

对扑翼气动功率的比较分析也 发现，人造扑翼机需要的气动功率明显高出同等大小的鸟类，在效率方面尚不及扑翼飞行生物。

关键词:微型飞行器; 柔性扑翼; 非定常流; 气动特性; 数值模拟 中图分类号:V 211． 1 + 5文献标识码:A动机理方面发挥着重要作用。

孙茂教授团队［8 － 9］、R a - m a m u r t h i ［10］、W a n g ［11］等在昆虫扑翼高升力机理的数值 研究方面做了大量的工作，并有不少重要的发现。

仿鸟 扑翼由于尺寸更大，速度更快，其雷诺数比昆虫扑翼更 高些。

对仿鸟扑翼气动特性的数值模拟研究，比较有代表性的工作有，R og e t 等［12］ 通过建立适当的动网格模 型，数值研究了弯折变形扑翼的流场; 左德参等［13］针对仿鸟刚性薄板扑翼、宋书恒［14］、谢辉［15］等针对刚性 NA- CA0014 矩形扑翼，考察了扑动参数对气动特性的影响。

这些研究都有力地推动了扑翼非定常空气动力研究的 进展。

但目前针对仿鸟扑翼的数值研究大多将扑翼假 设为刚性不变形模型，不能反映柔性扑翼流场的真实情 况; 研究内容也主要着重于考察各扑动参数对总气动力 的影响，较少关心柔性参数的协调搭配和气动力能否配 平扑翼机重力、前飞阻力问题; 扑翼机的功率消耗和效 率问题也较少涉及。

本文以南京航空航天大学研制的某型扑翼微型飞 行器为基础，建立适当的柔性扑翼模型，在非结构可变 形动态嵌套网格上，数值模拟该仿鸟扑翼的低雷诺数非 定常流场，并以升力平衡重力、推力平衡阻力为条件，构 建重要扑动参数之间的协调搭配关系，研究仿鸟柔性扑 翼的升力推力产生机理和气动功率消耗，分析流场涡结 构的分布和发展情况，同时尽可能地将数值计算结果与0 引 言随着微型飞行器研究的深入，扑翼越来越受到重 视，根据仿生学和空气动力研究成果预见，翼展小于 15cm 的微型飞行器，采用扑翼方式比固定翼或旋翼更具优势［1 － 2］。

近年来，有关扑翼空气动力的研究，在理论、实验和数值计算方面都已逐步开展。

在实体模型方 面，能实现扑翼飞行的微型扑翼机也相继研制出来，较成功的有“M i c r o B a t ”［3］，南京航空航天大学［4］、西北工业大学的扑翼机［5］等。

目前研制成功的扑翼机大多都有一个共同的特 点［6］: 翅膀由前缘骨架辅以柔性薄膜构成，体型尺寸和 中等偏小的鸟类相差不多，扑动方式也与该体型的鸟类 类似，是仿鸟单对扑翼飞行。

以仿鸟飞行为基础的单扑 翼是研究扑翼飞行器的基础，但即使是单扑翼，仿鸟扑翼的研究也仍然面临相当大的挑战［7］: 其一，很难研制 出像鸟翅膀那样的能弯折、变形、张开和收拢的符合空 气动力学规律的羽翼; 其二，鸟翼复杂的扭转扑动和柔 性变形; 其三，目前的仿鸟扑翼机在稳定性、机动性和效 率等方面仍与鸟类相差很远。

因此，扑翼的相关研究任 重道远，扑翼的空气动力原理及气动特性更是首先要探 明的基础问题。

随着计算流体力学的发展，数值模拟在揭示扑翼气收稿日期:2010-09-03; 修订日期:2011-03-02基金项目:中国博士后科学基金( 20090461113) ; 江苏省博士后科研资助计划基金( 0902086C) ; 南航大基本科研业务费专项( N S 2010025) 作者简介:肖天航( 1979 － ) ，男，湖南衡阳人，讲师，研究方向: 计算流体力学、飞行器设计． E -m a i l : xthang@ n u aa ． e d u ． c n*710 空气动力学学报第29 卷相关实验比较，以期为柔性扑翼微型飞行器的研制和改进提供必要的理论和技术支持。

柔性扑翼模型研究的仿鸟扑翼微型飞行器由一对可上下扑动的主翼和置于其后的固定的平尾、垂尾组成，以电机和机械组成的驱动机构实现主翼的上下扑动( 如图1 所示)。

为使数值计算和分析有针对性，这里以该扑翼微型飞行器为基础建立计算模型，如图2( a) 所示。

由于主翼扑动引起的非定常流场才是我们关注的重点，所以该计算模型只包括主翼和简化的机身，暂不考虑主翼后面的平尾和垂尾，扑翼也被简化成等厚度的薄板。

1( a)计算模型( b)扑动参数定义图2 仿鸟扑翼的计算模型和扑动参数定义F i g．2 Computational m o d e l of the b i r d-l i k e f l a pp i n gwing and d e f i n i t i o n of i t s motion p a r a m e t e r s度的材料用于扑翼非常重要。

基于以上分析，可设定仿鸟扑翼模型的运动规律如下:绕体轴线的扑动角度，Φ( t) = Φm + Φ0 cos( 2πf t)( 1)式中，Φm、Φ分别为平均扑动角度和扑动幅度，f为扑动频率。

各个截面绕前缘的弦向俯仰角度，α( y，t)= f( y) ［αm，t i p+ α0，t i pcos( 2πf t +90°)］( 2)式中，αm，t i p、α0，t i p分别为最大展向截面(翼梢处)的平均俯仰角和俯仰幅度，这两个参数也分别反映了扑翼的初始静态扭转程度和扑翼材料的柔性度。

f( y) 为弦向俯仰角度沿展向的分布规律，这里假定从翼根到翼梢俯仰角线性增大，即f( y) = y/b( 其中y 为展向位置，b为展长)。

上述部分扑动参数可根据扑翼机的传动机构和巡航飞行状态一一设定:根据传动机构可确定平均扑动角度Φm=0°、扑动幅度Φ=36°; 在8m/s 的巡航速度下，基于平均气动弦长和巡航速度的雷诺数为54900，扑动频率为8Hz; 另外两个重要参数αm，t i p、α0，t i p尚不确定，这里以扑翼产生的周期平均升力配平重力和平均推力配平机身阻力为条件来合理选择这两个参数的大小。

扑翼机重约100g，在8m/s 巡航速度条件下(重量和巡航速度基本符合Tennekes 统计的飞行生物的比例关图1 仿鸟扑翼微型飞行器F i g．1 T h e b i r d-l i k e flapping wing M A V要使数值结果能反映真实的流场情况，首先重要的一步是设定合适的和真实扑翼相符的运动规律。

仿鸟扑翼微型飞行器翅膀的扑动本质上和鸟翼扑动一样是一种周期运动，主要可以分解成两个自由度方向的运动:( 1) 绕体轴线的上下扑动。

( 2) 绕前缘或某一位置弦向的俯仰运动。

扑动参数的定义如图2( b) 所示。

绕体轴线的展向扑动相对简单，其扑动角呈正弦或余弦变化的规律。

弦向俯仰运动则与扑翼的结构设计方式有关。

目前大多数研制成功的仿鸟扑翼机，包括本文的研究对象，扑翼都采用前缘骨架的结构设计方式，同时翼根处固定于机身。

弦向俯仰运动以扑翼前缘为轴，且沿展向各截面因相对气流速度和所受气动力大小不同而不相同，弦向俯仰角翼根处为零，翼梢处最大，使得整个翼面沿展向呈现扭转变形。

各个展向截面处的俯仰角和扑动角一样呈正弦或余弦规律，已有实验结果表明［16］，弦向俯仰运动与展向扑动运动相差90°相位。

需要说明的是，扑动过程中展向各截面的扭转角度很大程度上由扑翼材料的柔性决定，材料越柔软，扭转角将越大，选择合适柔软［17］系) ，配平重力需要的周期平均升力系数为CL= m g/( 0．5ρUr e fSref= 0．52) ( 其中，mg 为扑翼机重量，参考速度Uref为巡航速度，Sr e f为扑翼面积)。

2数值计算方法22．1 非定常N-S 方程及求解A r b i t r a r y-L a g r a n g i a n-E u l e r i a n( A L E)格式的非定常可压缩N-S 方程为:第6 期肖天航等:仿鸟柔性扑翼气动特性与能耗的数值研究711t∫Ω( t)WdV + ∮ Ω(t)( F(W)－ v W)d Sgn= ∮ Ω(t) F v d S( 3)式中W = ( ρ，ρu，ρv，ρw，ρe)T 为守恒量，F(W)、F 分别v为无粘通量和粘性通量，vg n为控制体表面的法向运动速度。

众所周知，可压缩N-S 方程求解低速流场会面临所谓“刚性”问题。

为此对方程(3) 的时间导数项进行预处理，并用双时间步推进求解。

即:( a)初始网格( b) 某时刻的嵌套网格图3 仿鸟扑翼模型流场的计算网格F i g．3 Computational g r i d s f o r t h eb i r d-l i k e flapping w i n gΓ∫QdV + WdVt∫Ω( t)Ω( t)2．3 算例验证选取一个来流条件和流动情况与本文研究对象类( F(W)－vg nW)dS = ∮F v d S( 4)Ω(t)式中:、t 分别为伪时间和物理时间，Γ为伪时间导数项的预处理矩阵，Q 为基于压力、速度和温度的原参变量。

该方程在空间上采用二阶迎风格式的有限体积法离散，物理时间采用二阶后向欧拉格式离散，两物理时间步之间用预处理的伪时间步推进，离散后的方程用隐式似的实验作为算例，来验证上述数值方法用在该类问题上的准确性。

F e j t e k等［20］在风洞中测量了NACA8318 直机翼绕翼根上下扑动的升力推力情况，该机翼弦长76mm，展弦比4．01。

这里取其中一个实验条件进行计算:来流速度20．4m/s、扑动幅度、扑动频率3．3H z，雷诺数107000，与本文要研究的扑翼机模型较为相似。

算例的网格尺寸标准和扑翼机模型保持一致，即机翼面网格大小、附面层厚度、空间网格增长比例都相同。

计算过程中，时间方向上的离散精度也保持一致，每个扑动周期分成200 个物理时间步，每物理时间步的内迭代进行15 步或残差下降2 阶。

图4 是计算得到的周期变化的力系数与实验值和理论分析结果的比较，图4( a) 为升力系数、图4( b) 为推力系数。

就升力系数来说，无论是变化趋势还是幅值，计算结果、实验值和理论分析值三者都吻合得比较好，计算得到的周期平均升力系数为0．392，与实验值0． 390也很接近。

就推力系数来说，因为量值较小，实验值得到的曲线并不光滑，其峰值与理论分析结果也有较大差迭代求解，湍流粘性系数通过求解L U-S G SS p a l a r t-A ll m a r a s湍流模型得到。