刘文君 蔡 毓

（广西大学计算机与电子信息学院，广西 南宁 530004）

自然界中的鸟类能够在复杂和高度动态的环境中飞行，具有稳定性、机动性和高效性，对环境干扰有很强的抵抗力，从而可以进行长距离迁徙。然而，由生物飞行灵感启发研制的各种规模结构、传感、稳定性控制的仿生扑翼飞行器还没有自主地表现出这些特性[1]。仿鸟型扑翼飞行器和仿昆虫型扑翼飞行器的飞行机理有较大差别。仿昆虫型扑翼飞行器主要靠扑翼扑动产生升力和推力，由于其扑动频率高，且机翼质量占比小，所以气动、结构和飞行力学之间的耦合问题影响较小。对于仿鸟型扑翼飞行器而言，其具有低扑动频率、高机翼质量占比和耦合关系性强等特点。

气动效率是仿生扑翼飞行器研究中备受关注的领域，是该飞行器在飞行状态下受到的升力、阻力、推力、力的方向、力的大小以及稳定性等影响的客观反应，也是检验飞行器外形布局及气动关系设计是否合理的重要依据[2]。目前，关于扑翼的气动效率，在生物观测、理论分析、数值模拟、飞行试验4 个层面都在逐步开展研究。在数值模拟方面，Young J[3]采用可压缩二维Navier-Stokes 求解器，对斜向振荡的NACA0012 翼型的流场 进行了数值模拟，用数值粒子追踪法对翼型的尾迹进行了可视化。Lai JCS[4]研究了Strouhal 数推进效率的峰值，得出仅用Strouhal 数不足以表征扑翼推进效率的结论。Trizila P[5]使用代理建模技术对前缘涡、尾流及射流进行了数值研究。肖天航等人[6]提出一种Delaunary 图映射网格变形技术和非结构嵌套网格的数值计算方法。

虽然针对扑翼运动的研究已经初步揭示了产生升推力的飞行机制，然而定性的认识还不能完全满足扑翼飞行器设计层面的需求。一个重要的俯仰运动学案例[7]在中冲程、中下冲程不同的迎角下，分析二维、三维CFD 模型流场非定常瞬态状态，二者经验系数的变化曲线基本一致。该案例表明二维模型为微型飞行器设计中的实际应用提供了重要的仿真依据。该文在进行二维数值模拟的同时，对结果进行了积分计算，得到整个机翼平面的气动力和气动效率。

1 实验模型

1.1 扑动机构模型

在研制仿生扑翼飞行器时，可借助仿生学原理，参考鸟类的翼型构造、形状[8]，来研制机翼结构，如图1 所示。飞行时产生升力克服其自身重量，产生推力克服飞行时的阻力。

该文在实验中的飞行器机构简图，如图2（a）所示，并将其机翼扑动基于叶素理论将机翼平面按展向分成若干翼带，进行三维数学建模，如图2（b）所示。

图1 鸟翼的构造和不同阶段翅膀的形状变化

1.2 机翼运动分析

机构简图的运动角度位移是分析机构运动性能、优化机构的依据。机翼的运动分析，是根据原动件的运动规律，求解从动件的运动规律。机构运动如图3 所示。

图2 扑翼飞行器模型机构简图及机翼建模

在机构简图中的支链上引入闭环结构，利用闭环矢量方程建立连杆机构的约束方程。闭环方程为

式中：R1是曲柄，R2、R3、R4为机械连杆，R5为机构中的辅助杆。R1绕固定铰链点O 转动，即A 点是以O 为原点、R1为半径的圆上运动的动点。α1是连杆机构的设定角度，α2是R1绕定点O 转过的瞬时角，α3是翼杆AB 绕动的瞬时扑动角，α4是平板翼型中的机翼机构的位置角。沉浮运动参数示意图如图4 所示。冲程平面角为0。全局坐标系（X，Y，Z）固定在冲程平面的中心，Z 方向垂直于冲程平面，Y 方向垂直于机翼轴，X 方向与冲程平面平行。扑翼运动规律接近于简谐运动，扑动飞行可以同时产生升力和推力。平板翼型沉浮运动位移被定义为

式中：hm是拍打振幅，f 是拍打频率，t 是拍打时间，h（t）是t 时间内的位移。

2 数值计算方法

该文的数值模拟基于ANSYS 软件平台进行。该平台集成了N-S 方程流场解算器、动态网格生成技术、空气动力学算法等多方面模块，可用于理想和多种设置气体/流体介质的定常/非定常数值模拟仿真计算。

图3 扑翼机构运动简图

流场计算的基本过程是在空间上用有限体积法，将计算域离散成许多小的体积单元，在每个体积单元上对离散后的控制方程组进行求解。这也是流场计算方法的本质。采用分离式求解器（Segregated solver），逐一、顺序地求解各变量。

假定初始压力和速度，确定离散方程的系数及常数项；求解连续方程、动量方程、能量方程；求解湍流方程及其标量方程；判断当前时间步上的计算是否收敛。

图4 沉浮运动悬停参数示意图

2.1 动态网格技术

在非定常流场数值模拟的动态边界运动中，动态网格生成是非常重要的组成部分。动网格计算中，网格的动态变化过程可以用3 种模型进行计算，即弹簧近似光滑模型（spring-based smoothing）、动态分层模型（dynamic layering）和局部重划模型（local remeshing）。该文采用弹簧近似光滑法和局部重划法。

弹簧近似法适用于该研究中的非结构化网格，网格拓扑不变，通过近似弹簧的压缩或拉伸，实现网格和模拟流场计算区域的变化。网格边界处位移变化后，弹簧系统会经过调整，使已破坏的原有平衡达到新的平衡，生成力的大小由胡克定律计算得到。在控制台窗口掌握区域范围、体积统计以及连通性信息。网格检查没有出现网格体积为负数的状态。该次数值模拟实验需对模型进行不同方向的简谐运动，多次模拟之后得到相应良好的设定尺寸，根据预先设置的最大、最小网格规格，同时采用网格重划法在网格更新计算区域进行网格重划。在数值模拟中，还要注意更新区域的完整性。

2.2 非定常数值模拟方法

数值模拟中的控制方程是具有恒定密度和黏度的非定常Navier-Stokes 方程[9]。基于动态网格，采用二阶迎风格式（second order upwind）有限体积法，考虑量值在物理上的值点分布，即曲线的曲率影响，具有二阶精度截差，同时也会考虑气流的流动方向，所以数值解在物理上是合理的。利用BLU-SGS 隐式算法进行时间推进，求解低速非定常雷诺平均N-S 方程。同一算例应用不同的湍流模型，得到的数值结果会存在很大不同，与此同时，采用湍流模型中试验最成功的一种模型，Spalart-Allmaras 湍流模型[10]，可以合理地处理边界层的黏性影响区域。

文中分别进行了稳态和瞬时模拟计算，对于稳态问题，压力修正法采用SIMPLEC 算法，对于瞬态问题，同时为了更好地满足动量方程、连续性方程，采用PISO 算法。

2.3 算例准确性验证

该算例重点讨论外部空气动力学应用的无量纲系数，对气动进行配置，计算典型的阻力系数。基于该文研究的低速问题，周围流场中的空气密度保持不变，为不可压缩流，并且不会出现马赫数。气动系数是关于雷诺数的函数CD=f（M，Re）。如图5 所示，平板模型垂直于流动气流，这种机构在任何常规配置中都会产生最大的阻力系数。

雷诺数计算公式为：

阻力系数计算公式为：

其中，S 是平板单位跨度的正面面积。文献[7]中的阻力系数为2.0，该文数值模拟方法下，迭代收敛后的阻力系数为1.9986。数值模拟和理论方法对该算例进行比较，数据吻合，所以将利用该数值模拟方法进行进一步研究。

2.4 数据处理方法

对不同参数的二维模型进行非定常流场数值模拟，得到已设置的气流条件下的升力系数和阻力系数，与此同时，将每个翼带的二维气动力，通过积分计算形式，得到整个机翼平面气动力。后续通过如下数据处理方法做进一步研究。

减缩频率定义为

雷诺数表示迁移惯性力与黏性力的比值，定义为

经过无量纲化，拍打运动的升力系数和推力系数定义为

式中：CL、CT分别是升力系数和推力系数定义，L、T 分别是升力和推力，是来流密度，S 是平板跨度的正面面积。

依据文献[11]，平均推力系数和平均气动力输入系数分别定义为

同时平均气动力系数定义为

依据文献[12]，升力效率和推进效率分别定义为

3 数值结果与分析

对于扑翼飞行，机翼产生显著升力来抵消自身重力。气动效率是与给定功率的升力和推力有关的运动效率，能保持在空中的能源供应。该文从上述2 个方面对扑翼气动效率展开研究。参数空间包括扑打频率f，风速，减缩频率k，雷诺数Re。在二维扑翼无粘平板经典理论中，对于低振幅、低减缩频率的情况，扑翼在周期内的平均推力不依赖初始攻角[13]，并且对于三维流动的扑翼，只要运动不发生明显的流动分离现象，同样适用该理论[13]。该文数值研究中平均攻角α 设置为零。

3.1 扑动频率f 和风速对气动效率的影响

低频振荡是导致流动失稳的一个主要因素[13]，由飞行速度和扑打频率共同作用。对于扑翼飞行，机翼产生显著升力来抵消自身重力。在不同扑打频率f 和风速下，总气动力会随扑打频率的提高而增大。在飞行器起飞阶段，提高扑打频率，增加推力，升力会在达到一定扑打频率后有增加的趋势。结合如图5 所示的气动效率实验结果，在扑打频率为7 Hz、8 Hz 时，具有较高速度，达到最大升力效率。进入巡航阶段，降低扑打频率，达到最大推进效率。

3.2 减缩频率k 和雷诺数Re对气动效率的影响

该节实验计算考虑有限翼展的拍打运动，针对k=0.2～1.5、Re=1.7×104～5.2×104的减缩频率、低雷诺参数进行大量数值计算，得到定性以及定量的结果。

随着减缩频率k 的增加，实验中雷诺数下的气动力均增加。升力效率和推进效率分别如图6 所示。计算出的气动效率表明，随着雷诺数Re 的增加，实验中减缩频率k 下的推进效率逐步下降至接近平缓，而升力效率的结果变化更陡峭，出现较高的峰值。

从续航时间考虑，应以较低的扑动频率飞行，较高的扑动频率对扑翼飞行没有本质优化，却会造成机翼架构的损坏，正负升力抵消。

图5 气动效率随扑打频率和风速变化的比较

图6 气动效率随减缩频率和雷诺数变化的比较

4 结论

为了探讨仿生扑翼飞行器产生升力和推力的飞行机理，该文在研究考虑扑翼气动力和气动效率数值计算方法的基础上，利用数值模拟仿真，针对不同扑动频率、风速、减缩频率、雷诺数、扑动幅度等参数变量，研究揭示了该飞行器扑翼机构数值模拟计算的总气动力和气动效率的变化，从而研究分析主要飞行规律。为仿生扑翼飞行器的进一步研究提供了参考依据。