旋翼翼型气动设计与评估软件HRADesign

2021-09-17孙俊峰卢风顺

空气动力学学报 2021年4期

孙俊峰，卢风顺，黄勇，江雄，牟斌，许勇

(中国空气动力研究与发展中心，绵阳 621000)

0 引言

直升机由于其独特的飞行能力受到各航空大国的普遍重视。旋翼系统是直升机的关键部件，是直升机主要的升力面、推力面和操纵面。旋翼桨叶是由翼型构成的，翼型和桨叶气动外形对旋翼性能有重要影响，高性能旋翼翼型可以提高旋翼悬停效率3%－5%，对直升机的前飞速度、等效升阻比和机动飞行能力、噪声水平等都有很大影响。高性能旋翼翼型设计技术也一直是直升机设计的核心技术，是衡量直升机技术水平的重要标志。

在不同的飞行状态和不同的桨叶半径位置上，旋翼桨叶翼型的运行环境是迥然不同的[1]。翼型的设计需要在较宽的马赫数范围内，有较高的静态和动态最大升力系数，以适应机动过载状态；在较高的马赫数及小迎角时，有较大的阻力发散马赫数，以推迟前行桨叶激波失速；在中等马赫数及中等迎角时，有较高的升阻比，以提高旋翼的悬停效率；在较低的马赫数及大迎角时，有较好的失速特性，以延缓后行桨叶的气流分离；在整个飞行包线内，有较小的俯仰力矩系数，以降低桨叶的操纵载荷，旋翼翼型设计需要综合考虑不同的飞行环境，具有明显的多点、多目标、强约束的特点。

旋翼翼型的设计方法主要分为反设计方法和数值优化设计方法两大类。反设计方法主要是根据给定的速度分布或压力分布设计翼型，缺点是理想的翼型速度分布、压力分布特性难以预先给定。随着CFD技术和计算机技术的不断进步，数值优化设计类方法逐渐成为翼型设计的主流方法，该类方法可以直接实现多目标气动性能的优化，降低了对人经验的依赖，符合旋翼翼型设计多点/多目标的应用环境。数值类方法又包含控制论方法[2-3]和进化算法[4-5]两大类，分别适用于不同的优化环境。

从20世纪70年代开始，直升机技术发达国家都开始致力于新型旋翼高性能翼型的研究，发展了一系列旋翼专用翼型，如美国Boeing-Vertol公司的VR翼型族、Sikorsky公司的SC翼型族，法国宇航院（ONERA）的OA翼型族，俄罗斯的TsAGI翼型族等[6-7]。这些翼型的成功应用，对改善直升机性能起到了重要作用。目前这些翼型族仍在继续发展中，但鲜有相关文献报道。国内也从20世纪90年代开始了对自主旋翼翼型的相关设计研究工作。尚克明等[8]基于Euler方程采用反设计的方法进行了旋翼翼型的设计研究。刘刚等[9]基于进化算法和Kriging模型开展了旋翼翼型多目标优化设计研究，优化结果经过风洞试验验证，满足了设计要求。杨旭东等[10]基于梯度信息改进的响应面方法建立了旋翼翼型多点多约束气动优化策略。韩忠华等[11]发展了基于Kriging模型与遗传算法的旋翼翼型多目标多约束气动优化设计方法，应用于OA209的设计取得了较好的效果。招启军等[12]开展了考虑旋翼翼型定常-非定常影响的综合优化设计研究，优化结果在设计状态下明显改善了翼型的动态失速特性。

近年来直升机工业得到了长足的发展，先进直升机对高性能旋翼翼型的需求愈发迫切。目前国内旋翼翼型的研究正处于探索性研究向旋翼翼型系列研制的起步发展阶段，自主翼型系列的研制还缺乏通用、高效、鲁棒的设计工具，还没有自主翼型应用于型号研究的先例。高性能旋翼翼型的设计，需要有先进的设计和优化框架作为支撑，以框架体系为基础消化吸纳各类先进技术，促进优化系统的不断迭代更新完善，才能更有效地促进旋翼翼型工程实际设计能力的不断提升。

本文通过综合采用基于进化算法的多目标优化设计方法、基于主成分分析的多目标降维技术、翼型参数化技术以及高精度CFD性能分析工具，自主开发了旋翼翼型气动设计与评估软件系统——HRADesign。本文主要介绍了系统的架构设计和工作流程，详细介绍了系统的多目标优化方法模块、多目标降维技术、几何管理模块以及CFD性能分析模块等功能模块。利用该系统对ADODG标准翼型优化算例、旋翼翼型的常规多目标优化设计以及考虑多目标降维的旋翼翼型优化设计进行了研究，验证了优化设计软件系统的有效性。最后部分对本文进行了总结。

1 软件架构设计

HRADesign作为一个综合的旋翼翼型设计和评估的软件系统，系统开发基于Eclipse集成环境，选用Python和C++语言作为开发工具，XML文件用于交换信息，采用wxWidgets软件包提供可视化支持。用户可以在PC机上完成优化问题的描述和参数的输入后，通过集群系统完成优化设计的流程。系统根据高内聚、低耦合的原则分成层次结构，包含用户界面层、功能模块管理层、数据传输层、以及基础服务层。图1给出了软件的体系结构示意图，各层的功能定位如下：

图1 软件体系结构Fig.1 A diagram of the software architecture

1）用户界面层。界面层给用户提供了友好的人机交互接口来操作整个软件，借助该层功能，用户可实施优化设计问题建立、优化进度监控、优化结果查验等操作。

2）功能模块管理层。优化设计层集成了优化方法、代理模型、翼型参数化、优化目标处理、PCA降维技术、专家系统等模块，借助CFD或试验数据库的数据支撑，完成各种气动外形的优化设计功能。

3）数据传输层。该层主要负责为自适应优化设计层生产、收集和存储数据，提供目录管理、文件传输等基础功能。

4）基础服务层。该层提供整个软件运行所需要的基础服务，包括软件服务（CFD解算器、网格程序、数据库系统）和硬件服务（本地计算资源和大规模集群资源）。

图2给出了HRADesign系统进行优化设计流程的示意图。初始翼型进入系统，进行参数化并提取设计变量；选择试验设计方法，对每个样本点进行网格生成和气动特性分析；用气动分析得到的数据来构建代理模型，利用该模型分析优化问题的优化目标和约束；最后进入优化迭代流程，直至获得最终的优化翼型外形。

图2 系统优化设计流程示意图Fig.2 A diagram of the optimization process

2 系统功能模块

HRADesign系统由多目标优化方法模块、翼型参数化管理模块、代理模型模块、CFD性能分析模块等主要的功能模块构成，此外还包含问题定义，显示监控等辅助模块。系统提供了基于Windows的图形操作界面，用户在前台通过操作界面完成优化问题的定义以及各项参数的设置，系统底层封装了与后台集群系统的连接以及信息的交换，可以实现本地计算与集群系统计算的无缝切换，用户可以通过输入/输出系统监视优化进程，方便处理各类操作问题。图3给出了系统的主要应用界面示意图。

图3 软件应用界面示意图Fig.3 A diagram of the graphical user interface

2.1 多目标优化方法

多目标优化方法模块主要完成各种全局和局部数值优化方法的封装，包含有进化算法和伴随方法等，用户通过应用界面可以选择不同的优化方法完成翼型的鲁棒设计、PCA降维分析等多目标设计功能。

HRADesign系统是以优化算法为核心，驱动整个优化流程的发展来进行翼型的优化设计。系统采用的多目标优化方法大多是以进化算法为基础，结合多目标Pareto解的概念以及约束处理机制发展起来的，进化类算法属于全局类优化方法，优化过程不依赖目标函数与设计变量的梯度信息，适合处理旋翼翼型复杂流动中的各类非线性问题。进化算法通过模拟生物种群的进化过程，利用选择、交叉、变异等进化算子来找到多目标问题的Pareto解，解决了传统优化设计中多目标加权的权重系数难以给定的难题，得到广泛应用。图4给出了利用该方法进行优化设计的流程图。

图4 进化算法流程图Fig.4 The flow chart of the optimization process based on the evolutionary algorithm

2.2 主成分分析

主成分分析PCA的方法[13]是Deb在2005年提出的用于多目标降维的算法。算法首先对数据集进行标准化处理，求解目标函数的相关矩阵，得到相关矩阵的特征值及其特征向量，将特征值按从大到小的顺序排序并求出其贡献率，当贡献率的累积大于初始给定的参数阈值时停止累积。阈值的选取在很大程度上影响实验结果，如果阈值取太大，可能所有目标都会被选取，如果太小又容易丢失非冗余目标的信息。文献一般建议阈值取值为95%。算法简介如下[14]：

1）设置阈值，非冗余目标集合Ⅰ为空。

2）计算数据集。随机初始化种群，对种群中的个体计算目标集中的所有目标，得到个体目标的数据集合P。

3）对数据集P进行PCA处理，选取非冗余目标Ⅰ：

（a）对目标集进行归一化处理，计算相关矩阵，对相关矩阵计算特征值及对应的特征向量，依据特征向量来选取目标。

（b）将特征值按照从大到小的顺序依次排序，计算每个特征值占总特征值的比率，则每个特征值对应的特征向量依次称为第一、第二、···、第n主成分。

（c）依次分析每个主成分。如果主成分里的元素有正有负，选取最大和最小元素对应的目标加入非冗余目标集Ⅰ；如果所有元素均为正，选取最大元素值对应的目标加入非冗余目标集Ⅰ；如果所有元素均为负，将所有目标均加入到非冗余目标集Ⅰ。

（d）如果主成分对应的特征值小于0.1，选取最大元素对应的目标加入到非冗余目标集Ⅰ。

（e）考察特征值比率，如果大于阈值，停止分析过程，输出非冗余目标解集Ⅰ；否则转到（c）继续非冗余目标的选取。

在优化设计中通过PCA主成分分析得到各目标之间的关联关系后，提取决定问题本质的主要目标，将冗余目标剔除、或者转化为约束条件，将高维多目标优化转化为低维优化问题。可以解决多目标优化收敛慢甚至不收敛的问题，提高优化结果的可靠性。图5给出了基于PCA分析的优化设计流程图。

图5 基于PCA分析的多目标优化设计流程图Fig.5 The flow chart of PCA

DTLZ测试函数[15]是由Deb提出的一组测试多目标优化算法性能的测试函数，共有9组函数，DTLZ测试函数有已知的Pareto最优解，其中DTLZ5(I,M)（M表示目标个数，I表示非冗余目标个数），用于测试算法处理包含冗余目标的能力。选取DTLZ5（2,10）函数用于测试，该函数有10个目标，其中2个为非冗余目标，函数公式见式（1）：

对该函数进行优化测试，经过PCA分析，得到目标9和目标10两个非冗余目标。首先对10个目标进行优化，得到优化结果对应的两个非冗余目标的Pareto前沿如图2（a）所示，去除冗余目标后再进行优化，得到优化结果如图2（b）所示。DTLZ5（2,10）函数的Pareto前沿收敛到圆弧曲线[16]，通过图6可以看出，算法在处理冗余目标后，达到了最优Pareto前沿。

图6 DTLZ5（2,10）测试函数收敛比较Fig.6 Convergence of the test function DTLZ5(2,10)

2.3 翼型参数化

翼型参数化管理模块主要完成翼型的输入/输出，约束评估、翼型曲线参数化和网格自动重构、设计变量的选取和确定设计空间范围等功能。目前常用的参数化方法有解析函数线性叠加法，NURBS曲线[17]、CST[18]方法、FFD[19]方法等。本系统主要采用CST技术实现翼型的参数化表示。

CST方法是波音公司B.M.Kulfan等提出的一种通用几何参数化表示方法。可以用统一的解释函数表示钝前缘/尖后缘类翼型和双钝头翼型等新型翼型。用该方法描述翼型，容易控制前缘半径、弯度/厚度分布、后缘角以及后缘厚度等关键参数，而且设计参数数目容易控制，也具有局部修改控制的能力。

翼型上下表面均用以下公式描述：

其中类型函数：

形状函数：

翼型表面形状改变后采用基于双曲方程[20]的网格生成方法完成计算网格重构，可以直接输出给CFD解算器使用。

2.4 Kriging代理模型

优化设计过程中，采用高精度CFD分析工具进行翼型气动性能计算、优化目标函数和约束函数的评估，会带来计算资源成本过高、计算周期过长的问题。为了提高设计的效率，可以采用代理模型的方法构造目标函数和约束函数的近似函数，优化算法作用于近似目标函数和近似约束函数以寻找优化问题的最优解，通过近似函数的不断改进和优化算法的不断迭代，直至最终满足收敛条件。

HRADesign系统中，代理模型技术主要包括试验设计方法和代理模型方法两部分。试验设计方法决定了样本点的个数和样本点的空间分布情况，系统采用拉丁超立方采样[21]和均匀采样[22]的方法选取样本点，以保证样本点在设计空间的均匀分布。Kriging模型[23-24]则用来作为目标函数和约束函数的近似模型。

Kriging代理模型起源于地理空间统计学，是一种估计方差最小的无偏估计插值模型，具有全局近似和局部随机误差估计相结合的特点。通过Kriging模型可以得到未知点的函数值和不确定性，因而在优化过程中，Kriging模型需要根据优化进程自适应地更新，增加样本点在非线性区域的分布，提高模型的预测精度，系统采用EI方法或者最小值加点准则等方法，提高模型的精度和自适应能力。图7给出了基于代理模型的优化流程示意图。

图7 基于代理模型的优化流程示意图Fig.7 The flow chart of a surrogate-based multi-objective optimization

第一步：试验设计。采用均匀设计等试验设计方法在设计空间中选取样本点。

第二步：样本点性能评估。对样本点分别生成计算网格，采用CFD工具进行气动性能评估。

第三步：构建Kriging模型。利用样本点的性能计算结果构建初始Kriging代理模型。

第四步：利用构建的Kriging代理模型，结合基于进化算法的多目标优化方法进行设计优化流程。

第五步：Kriging模型的重构。优化过程中评估种群个体的EI值，根据种群个体的EI值选定附加的样本点，对样本点进行计算网格的生成和CFD性能评估，在N=N+ 1个样本点的基础上进行Kriging模型的重构。

第六步：返回第四步继续迭代循环。

第七步：直到系统收敛或达到设定的优化步数，返回优化结果。整个优化流程结束。

2.5 气动性能分析

旋翼翼型CFD性能计算的精度和效率是进行气动设计优化的关键，对优化设计的结果有直接影响。

HRADesign系统中采用的CFD分析工具是自主开发的二维RANS解算器MBNS2D。MBNS2D采用格心型有限体积方法求解雷诺平均Navier-Stokes方程，空间离散采用Roe格式，湍流模型包含SA一方程模型和SST两方程模型，通过多重网格技术进行流场加速收敛，提高计算效率。

在旋翼翼型的计算中，转捩对翼型前缘流动影响很大，如果不考虑转捩的影响，阻力及零升阻力的计算结果会与试验差别很大，为了增强方法的适应性，HRADesign系统中采用了 γ -Reθ湍流转捩模型[25]来提高翼型阻力系数的计算精度。

3 优化算例及分析

3.1 算例1

为了提高优化设计软件的可信度，AIAA气动优化设计讨论组（ADODG）给出了一套优化设计的标准算例，用于优化设计软件的验证和确认。这些算例包含了翼型和机翼在气动和几何约束条件下阻力最小化优化设计问题。本文以其中跨声速条件下RAE2822翼型阻力最小化算例来考核优化设计系统的性能。

以RAE2822翼型作为初始翼型开始优化，翼型的跨声速设计条件为：Ma∞= 0.734,Re= 6.5×1 06，CL= 0.824，选取阻力系数CD最小化作为优化目标，几何约束条件为保持优化后翼型的面积不减少。优化问题数学描述为：

图8和图9分别给出了优化前后翼型的外形和表面压力分布的比较，图10给出了优化前后翼型压力云图的比较。可以看出，优化翼型前缘吸力峰增强，基本消除了初始翼型中段的强激波，等值线变得平顺光滑。表1给出了翼型在优化前后的气动特性比较，在保持升力和力矩的气动约束以及面积不减的几何约束条件下，优化翼型的阻力系数降低了约87 counts，优化翼型的阻力特性得到明显改善。文献[26]对该算例给出了进一步的分析。

表1 优化前后翼型性能比较Table 1 Performance of the original and optimized RAE2822 airfoils

图8 优化前后翼型的形状比较Fig.8 Configurations of the initial and optimized RAE2822 airfoils

图9 优化前后翼型表面压力分布比较Fig.9 The comparison of airfoil pressure distributions between the initial and optimized airfoils

图10 优化前后翼型压力云图比较Fig.10 The comparison of pressure contours between the initial and optimized airfoils

3.2 算例2

直升机旋翼翼型的设计，需要考虑多种飞行状态下的性能，一般要求：较高的最大升力系数；较高的阻力发散马赫数；在较大的马赫数范围内有较高的升阻比以及较小的俯仰力矩，属于典型的多目标设计优化问题。以某翼型为基本翼型，在保持厚度不减的条件下开展多点设计优化，使得前飞、机动、悬停条件下满足设计指标的要求。优化问题描述为：

其中Maddm为阻力发散马赫数，K为升阻比，下标0表示基本翼型的性能指标。

图11给出了优化翼型和基本翼型外形的比较，优化外形头部半径增大，上表面厚度增加。

图11 优化前后翼型的形状比较Fig.11 Configurations of the initial and optimized airfoils

图12给出了机动状态Ma= 0.4下升力特性的比较，可以看出优化翼型最大升力系数有明显提高。图13给出了悬停状态下极曲线的对比，可以看出优化翼型的阻力系数有所降低。表2给出了优化翼型和基本翼型的性能结果比较，从结果可以看出，优化翼型的机动和悬停性能有明显提升，前飞性能在保持力矩性能的前提下阻力发散马赫数略有降低。该厚度翼型主要位于桨叶的中段，在满足厚度和力矩约束的条件下，需要兼顾低速和高速性能，要提高前飞阻力发散马赫数，就要损失低速升力性能，所以多目标优化设计中往往需要在性能之间进行折衷平衡。

表2 优化前后翼型性能比较Table 2 Performance of the initial and optimized airfoils

图12 翼型升力特性比较（ M a=0.4, R e=2.8×106）Fig.12 The comparison of lift coefficients between the initial and optimized airfoils （ M a=0.4,Re=2.8×106）

图13 优化翼型与初始翼型极曲线特性比较（Ma=0.6, Re=4.2×106）Fig.13 The comparison of polars characteristics between the initial and optimized airfoils （ M a=0.6, R e=4.2×106）