石胜庆，张祥谋，程钰锋

(1.江西上饶海军91551部队，江西 上饶 334000;2.北京装备学院研究生院，北京 101416)

0 引言

自1903年莱特兄弟发明飞机以来，螺旋桨动力系统便作为飞行器的主要动力装置之一而被人们广泛认识和深入研究［1］。

近年来，国内对于螺旋桨气动方面的各种研究也逐渐增多。2008 年，中科院的聂营、王生等人［2，3］采用Gambit软件对螺旋桨进行几何建模，再用Fluent软件基于滑移网格，仿真研究了临近空间螺旋桨的气动性能。2008年，西工大的许建华、宋文萍等人［4］采用雷诺平均N-S方程和嵌套网格技术对美国国家可再生能源实验室的CER实验型风力机叶片进行了仿真研究，湍流方程是B-L代数模型;计算中采用了有限体积空间离散法和改进型隐式LU-SGS时间推进格式，所用的嵌套网格技术有效地捕捉了螺旋桨的尾涡;在此基础上，2009年［5］，他们采用在不同粗细网格上消除不同频率误差加速解的收敛的多重网格技术，仿真研究了螺旋桨侧流粘性流场，其中多重网格采用非线性方程的全近似格式(FAS)。2010年，西工大的罗淞、杨永等人［6］采用多块点对点网格生成技术，生成了螺旋桨的空间计算网格，分别研究了欧拉方程组和N-S方程组对不同进距比下螺旋桨运动过程的仿真结果，与实验进行对比后得出欧拉方程比较适合于螺旋桨的工程设计。

以螺旋桨动力系统作为动力装置的飞行器在上升、下降或转弯的过程中，合成气流的实际速度方向与桨盘所在面并不垂直，这时螺旋桨处于有偏角的侧流状态，其运动特点会发生改变。本文基于滑移网格模型，通过求解三维非定常N-S方程，数值模拟了螺旋桨侧飞状态时的运动过程，仿真研究了侧偏角对螺旋桨气动性能的影响。

1 数学模型

对于N-S方程，连续方程、动力方程和能量方程的通用形式可以写成如下形式。

其中:ρ是气体密度，U是速度矢量，φ是通用变量，Γ是广义扩散系数，S是广义源项。对于连续方程、动力方程和能量方程，φ分别为1、ui和T;Γ分别为0、μ 和k/cp;S分别为0、-∂p/∂xi和ST。ui是速度分量，T是温度，μ是粘性，k是流体的传热系数，cp是比热容，ST是粘性耗散项，即流体的内热源及由于粘性作用流体机械能转换为热能的部分。

理想气体状态方程为:

式中R是气体常数。

湍流模型采用RNGk-ε两方程模型，它是基于k-ε标准两方程的湍流模型，采用一种叫做重正规化群的数学方法对N-S方程进行暂态推理得到的改进型k-ε两方程湍流模型，它是由V.Yakhot和 S.A.Orszag于1986 年提出并逐步完善的［7，8］，其基本思想是认为，在流场中小涡是各向同性的，处于统计定常和统计平衡的状态。RNGk-ε湍流模型考虑了低雷诺数流动粘性，改进了标准k-ε模型的高雷诺数性质，并且提供了Prandtl数的解析公式，考虑了湍流漩涡，因此更加适合于雷诺数不是很高和带有强漩涡运动状态的数值仿真。

采用耦合求解器，首先同时求解连续方程、动力方程和能量方程，然后求解湍流方程。在耦合算法中使用隐式格式，即通过求解方程组的形式求解流场变量，它是使用块Gauss-Seidel法与AMG法(Algebraic Multi-Grid，代数多重网格法)联合完成的。

采用二阶精度的有限体积AUSM(Advection Upstream Splitting Method)离散格式对粘性流体的控制方程和湍流方程进行空间离散。另外，本文采用了滑移网格模型来处理螺旋桨的运动过程，详见文献［9］。

采用滑移网格模型处理螺旋桨的旋转运动过程。滑移网格是在动参考系模型和混合面法的基础上发展起来的，常用于风车、转子、螺旋桨等运动的仿真研究。在滑动网格模型计算中，流场中至少存在两个网格区域，每一个区域都必须有一个网格界面与其他区域连接在一起。网格区域之间沿界面做相对运动。滑动网格模型允许相邻网格间发生相对运动，而且网格界面上的节点无需对齐，即网格交界面是非正则的。在使用滑动网格模型时，计算网格界面上的通量需要考虑到相邻网格间的相对运动以及由运动形成的重叠区域的变化过程。详见文献［9］。

2 螺旋桨模型及边界条件

2.1 模型及网格

由叶素理论可知，叶素的气动性能对螺旋桨整体气动性能有很大影响。由于临近空间环境大气密度较低，叶素雷诺数较低，因此本文选择低雷诺数翼型SD-8000PT作为螺旋桨的叶素［10］。螺旋桨的安装角和弦长详见文献［10］。两桨叶螺旋桨直径D为1.5 m，两个桨叶之间的桨毂简化为圆柱面，其半径和长度分别为0.01D和0.05D。

计算区域是一个长10D、直径8D的圆柱体。速度入口距螺旋桨4D，给定气流速度及总温，压力出口距螺旋桨6D，给定总温和总压;远场距螺旋桨转轴4D，给定气流速度、总压及总温。

由于滑移网格模型允许相邻网格之间发生相对运动，而且网格界面上的点无需对齐，即网格是非正则的。利用这一特点，可以更好地分布网格的疏密度，既保证了计算流场所需要网格数，又使网格总数减少，从而节约计算资源。

对螺旋桨的仿真而言，螺旋桨近区域流场变化剧烈，所以应该加密网格，远离螺旋桨的区域流场较为平缓，对网格数目要求不高，所以网格数目较少。基于此，本文将计算区域分为气流入口区域、旋转区域和气流出口区域3个小的计算区域。气流入口和出口区域采用TTM［11］方法生成结构网格，并在各自靠近螺旋桨的一端加密网格。旋转区域采用TGrid网格划分法生成非结构网格，螺旋桨表面网格节点之间相隔1 mm。

图1 计算区域及网格示意图

2.2 边界条件

以螺旋桨动力系统作为动力装置的飞行器在上升、下降或转弯的过程中，合成气流的实际速度方向与桨盘所在面并不垂直，这时螺旋桨处于有偏角的侧流状态，其运动特点会发生改变。本文主要研究螺旋桨处于此类侧流工况下的气动性能。

图2是笛卡尔坐标系下飞行器运动方向即螺旋桨速度方向与桨盘关系的示意图。图中，红色加粗线段CO是飞行器的运动方向即飞行器和螺旋桨的速度方向，桨盘所在面是坐标面y-z，BO是CO在x-z平面上的投影，O是坐标系的原点。假设飞行器的运动速度矢量为V，螺旋桨速度方向与x-y平面的角度α=∠AOB，螺旋桨速度方向与x-z平面的角度β=∠BOC，则螺旋桨速度矢量在x、y、z三个方向上的分量分别如下:

图2 笛卡尔坐标系中螺旋桨速度矢量示意图

本节主要研究在20km的临近空间环境下保持螺旋桨速度|V|为5 m/s和转速为300 rpm不变的条件下，偏角β=45°时，螺旋桨的气动性能随偏角α的变化规律。

3 仿真结果及分析

本文在保持偏角β为45°不变的情况下，数值模拟了不同偏角α的工况下螺旋桨的运动过程。表1是不同偏角条件下螺旋桨的气动参数，其中Tx、Ty、Tz分别表示螺旋桨在x、y、z三个方向上拉力的平均值，Mx、My表示螺旋桨拉力在x、y方向上扭矩的平均值，ηx是螺旋桨在x方向上的效率，效率计算公式中的拉力取螺旋桨在x方向上拉力即轴向拉力的平均值，扭矩取螺旋桨在x轴上扭矩即轴向扭矩的平均值。

由表1可见，螺旋桨在x方向上的拉力分量Tx随着偏角α的增大逐渐增大，由公式(3)可知，当偏角β=45°时，螺旋桨在x方向上的速度分量vx随着偏角α的增大逐渐减小，当螺旋桨的转速不变前进速度减小时，螺旋桨产生的拉力增大［9］，所以Tx随着偏角α的增大逐渐增大。螺旋桨在y方向上的拉力分量Ty有随α增大而逐渐减小的趋势。Ty是螺旋桨侧飞时偏角带来的阻力，Ty越大侧向阻力越大，对螺旋桨的飞行越不利。螺旋桨在z方向上的拉力分量Tz随α的变化规律不是很明显。比较可见，Tx的约为Ty的2倍，约为Tz的100倍，可见，螺旋桨在z方向的拉力分量可能是力学噪声，即可以认为Tz为0。

由扭矩分布可知，桨叶对x轴的扭矩Mx随α的增大逐渐增大，桨叶随y轴的扭矩My随α的增大有逐渐减小的趋势，这可能是由拉力变化而引起的。由效率的分布可见，不同偏角工况中螺旋桨的效率变化不大，最大变化幅度小于0.0051。

图3(a)是偏角β=45°、偏角α=0°时桨叶周围的流线分布图。由图3(a)可见，当桨径大于0.3 m时桨叶迎风面就会发生流动分离现象，且分离涡的尺度随着桨径的增大逐渐增大，虽然r=0.1 m叶素处于负攻角工况，但桨叶的背风面并没有发生流动分离现象。

表1 不同偏角下螺旋桨的气动参数

图3 桨叶周围流线分布比较图

图3(b)是偏角 β=45°、偏角 α=30°时桨叶周围的流线分布图。由图3(b)可见，在该工况下桨叶上发生流动分离的桨径范围是小于0.1r和大于0.3r的桨径范围。在小于0.1r的桨径区域内，叶素处于负攻角工况，流动分离发生在桨叶的背风面，在大于0.3r的桨径区域，流动分离发生在桨叶的迎风面，与前述情况一样，分离现象随着桨径的增大逐渐严重。

图3(c)是偏角 β=45°、偏角 α=60°时桨叶周围的流线分布图。比较图3(b)和图3(c)可见，偏角α=60°时桨叶周围的流动分离现象比偏角α=30°时桨叶周围的流动分离现象更为严重，首先，r=0.1处叶素背风面分离涡的尺度加大，其次，r=0.3处叶素迎风面分离涡的尺度加大。但这两个工况下桨叶上发生分离的区域是相同的，分离涡的尺度随桨径的变化规律也是相同的。

图3(d)是偏角 β=45°、偏角 α=90°时桨叶周围的流线分布图。由图3(d)可见，当α=90°时，螺旋桨桨叶周围的流场非常复杂，且与α为0°、30°和60°的工况有所不同。该工况下，螺旋桨桨叶背风面没有发生流动分离现象，但整个桨叶迎风面都发生了强度很大的流动分离，且分离涡的尺度随着桨径的增大逐渐减小。

4 结论

本文基于滑移网格模型，数值研究了螺旋桨侧飞状态时侧偏角对螺旋桨气动性能的影响。主要结论有:由于螺旋桨轴向速度随着偏角的增大逐渐减小，所以螺旋桨轴向拉力随着偏角的增大逐渐增大。螺旋桨的侧向阻力随着偏角的增大逐渐减小。桨叶上的分离区域随着偏角的增大逐渐增大。当偏角为0°、30°、60°时，桨叶上分离涡的尺度随着桨径的增大逐渐增大。当偏角为90°时，桨叶周围的流场非常复杂，且分离涡的尺度随着桨径的增大逐渐减小，与其余工况有所不同。

［1］Roosenboom Eric W M ，Heider A，Schröder A.Propeller Slipstream Development［C］.25th AIAA Applied Aerodynamics Conference，25-28 June 2007，Miami，FL，AIAA 2007-3810.

［2］聂营，王生，杨燕初.螺旋桨静推力数值模拟与实验对比分析［J］.计算机仿真，2009，26(3):103-107.

［3］聂营.平流层飞艇高效螺旋桨设计与试验研究［D］.北京:中国科学院光电研究所硕士学位论文，2008.

［4］许建华，宋文萍，韩忠华.基于嵌套网格技术的风力机叶片绕流数值模拟［J］.太阳能学报，2010，31(1):91-94.

［5］许建华，宋文萍，韩忠华，等.基于CFD技术的螺旋桨气动特性研究［J］.航空动力学报，2010，25(5):1103-1109.

［6］罗淞，杨永，左岁寒.螺旋桨流场Euler/NS数值解对比分析［J］.航空计算技术，2011，41(1):74-77.

［7］Yakhot V，Thangam S，Gatski T B，et al.Development of Turbulence Models for Shear Flows by a Double Expansion Technique［R］.NASA Contractor Report 187611，ICASE Report No.91-65.

［8］Rubinstein R.Renormallzation Group Theory of Bolgiano Scaling in Boussinesq Turbulence［R］.NASA Technical Memorandum 106602，ICOMP-94-8;CMOTT-94-2.

［9］程钰锋，聂万胜.基于滑移网格的临近空间螺旋桨流场数值仿真［J］.直升机技术，2012(2).

［10］Colozza A J.APEX 3D Propeller Test Preliminary Design［R］.NACA/CR-2002-211866，2002.

［11］Thompson.J F，Thames F C，Martin C W.Automatic Numerical Generation of Body-Fitted Curvilinear Coordinate System for Field Containing Any Number of Arbitrary Two-Dimensional Bodies［J］.J.Comput.Phys.1974(15):299-399.