民机空气动力设计中的数值优化方法

民机空气动力设计中的数值优化方法

反设计的数值优化方法

Lighthill利用保角变换的方法首先提出了二维翼型的反设计方法,Hicks,Murman和Henne等人将此方法发展为可应用于飞机设计的工程设计方法。后Campbell等提出过一种带约束的直接迭代的表面曲率(CDISC)方法,Yu将其与N-S解算器耦合形成了一种翼型和机翼的设计方法。波音公司则将此方法发展成工程应用的设计方法,并广泛地应用于波音的B757,B777和B737NG等型号的设计过程,取得了很好的效果。例如在B777研制中由于使用了反设计方法,仅经过三轮机翼的设计便取得了满意的结果,使风洞实验的机翼模型大大少于过去B757和B767设计时的数目,充分表明了该设计工具的作用。可以说,反设计方法曾对民机设计起过革新性的推动作用;但反设计方法也有其固有的弱点(参见文献[13]的附录D):首先,对于高度三维的流动要找到“好”的压强分布很困难;其次,不能保证所得结果为最优,即既具有高速巡航低阻的特性又在非设计条件下具有可接受的性能;最后,其他学科的约束会导致反复迭代。

低可信度CFD模型的数值优化方法

随着计算能力和数值优化方法的快速发展,应用基于CFD的数值优化方法于民机设计得到了很大的发展。这一方法的应用也从低可信度CFD模型开始,逐渐发展到采用先进的N-S方程解算器。波音公司发展了一种耦合TRANAIR[16](一种全速势方程的有限元方法,可参见文献[13]附录B)和梯度优化方法的数值优化气动力设计方法,并在1992年形成了TRANAIR优化器的雏形[17]。经过近十年的改进,得到了一个适用于位势流/边界层耦合飞行条件的气动力优化设计工具[18-20],具有多点优化设计能力,可处理高达600个几何自由度和45000个非线性不等式的约束条件(图1表示了TRANAIR优化过程示意图)。作为一个例子,图2给出了采用该软件对机翼/发动机短舱设计计算前后压强分布的对比,图a和图b分别表示了设计前后等马赫数线的分布。可以看出图a中挂架处出现激波;图b中短舱附近的机翼表面上消除了由于短舱干扰形成的激波。算例结果表明该设计软件可以处理很复杂的飞机/发动机综合设计问题。

高可信度CFD模型的数值优化方法现代优化算法可以分为依赖和不依赖梯度的方法两大类。

1.依赖梯度的优化算法

目前可用的大多数依赖梯度的数值优化方法都是从控制理论出发的,Jameson是此类方法的先驱者之一。尽管最初是由Pironneau提出利用控制理论进行椭圆方程系主控的外形优化的[21-22],但Jameson首先提出了通过控制理论自动进行外形优化的伴随方程方法[23]并应用于跨声速流动。后来,Jameson和他的合作者,还有其他研究者,大力发展此方法,从全位势方程到Euler/N-S方程,从无粘设计到有粘设计,甚至从气动设计到气动/结构的耦合设计,形成了大量文献[24-36]。此方法不同于一般梯度优化方法之处在于它将外形作为一个自由表面,促使流动解和最终优化的外形同时趋于收敛,因而使优化方法具有很高的效率(其基本思想可参见文献[13]附录D)。

2.不依赖梯度的优化算法

最早无需梯度的优化算法有Powell(共轭方向法)[37]和Nolder-Mead的单纯形法[38]。最近Sturdza还应用后者于空气动力的设计[39]。近二十多年来人们更多地使用诸如模拟退火法[40]和遗传算法(GeneticAlgorithm-GA)等的搜索方法,特别后者更为人们所关注。Holland利用进化理论创造了遗传算法[41](可参阅文献[13]附录D),即模仿生物的自然选择进行搜索以寻求最优解。与传统的搜索和优化方法相比,遗传算法具有下述4个特点[42-45]:1)不是直接作用于参变量集本身,而是对参变量集的某种编码运算。2)不是对单个点而是对多个点构成的群体进行搜索。3)直接计算适应值(函数),无需导数和其他辅助信息。4)利用概率转移原则,而非传统优化方法中的确定性原则。已有愈来愈多的研究和民机研制机构表现出了对这种随机寻优方法的浓厚兴趣,也已出现了不少利用遗传算法进行翼型或机翼优化计算的文献[46-56]。

3.对高可信度CFD模型数值优化方法的要求

分析最近十余年中出现的大量基于Euler/N-S方程的数值优化方法和文献,可以看出多数仍表现为学院式的探讨,提供可直接用于工程设计的方法和工具显得尚很有限,尽管已开始向这方面努力。这可能是因为:1)只是近几年来随DPW研讨会等的进行,数值模拟才可以比过去更正确地估算阻力值。2)工程界的空气动力外形优化需要在高维搜索空间中进行并存在大量的非线性约束,使优化问题十分复杂且计算开销巨大;3)巨大的计算量要求很丰富的计算资源和很长的计算时间,这与工程问题要求的迅速反馈相悖。

因此要使基于CFD的空气动力优化方法和软件成为日常的工程设计手段和工具需解决如下技术关键:1)具有建立准确计算诸如升力、阻力、力矩等敏感气动特性的正确流动模型的能力。比较现有的气动力优化方法可知,大多数方法还在使用不完善的流动模型,如基于Euler方程,甚至全位势方程等。虽然它们在一定条件下,如巡航小迎角飞行状态,可以提供合理的结果,但工程应用常要求准确地估算出阻力、俯仰力矩等敏感的气动特性,要求可计算整个飞行包线的飞行状态以及不同的复杂的几何外形等,这只能通过求解N-S方程来实现。顺便指出,有些文献(如文献[28])虽以N-S方程为主控方程,但优化时的伴随运算子却是在没有考虑粘性流动的假设下得出的(参见文献[28]第6节)。为了提高计算准确度,最好在离散N-S方程时使用高阶的差分算子[53-54]。2)具有寻求全局最优的能力。通常基于梯度的算法容易陷入局部最优,而遗传算法等随机搜索的方法则具有取得总体最优的优点。3)能有效地处理大量几何和气动力的非线性约束。优化问题的最优解常常是位于不同维超曲面(hyper-surface)的交汇处,遗传算法不同于基于梯度的方法,不限于目标函数的光滑扩展,可应用于多重约束的情况[53-54]。4)可应用于不同的几何外形和设计条件。5)扫描高维搜索空间的计算有效性高,以满足设计周期和研制成本的要求。遗憾的是这正是遗传算法的主要缺点,即估算适应函数的高代价。可以采用多处理器上的有效并行计算来大大减少计算时间[57],或在估算适应函数值时采用近似模型,如降阶模型[54,58]或响应面模型[50]等。

数值优化方法的发展现状和验证研究#p#分页标题#e#

1.空气动力优化设计计算的系列研讨会

近年来CFD学术界和航空业界都十分关注计算阻力的精度问题,这也是CFD应用于工程设计时所面临的第一个具有挑战性的计算。AIAA的应用空气动力学专业委员会在各方支持下,自2001年开始举行了DPW(DragPredictionWorkshop)系列会议[59],参与者都用N-S方程求解相同的几何外形(翼/身组合体,翼/身/短舱/挂架的复杂组合体等),得到了一个巨大的计算结果数据集,可与已有的已经过修正的风洞试验值比较。由于参与的计算者所采用的数值方法、湍流模型、计算网格形式及数目等各不相同,此数据集可用作分析和讨论各种因素对CFD计算结果的影响。该系列会议至今已举行了5届,对推动和提高CFD计算阻力的精度很有意义。文献[13]的附录C中给出了前3届结果的分析和讨论。鉴于DPW系列会议的成功,AIAA应用空气动力学专业委员会针对CFD面临的第二个挑战---计算三维高升力外形的最大升力CLmax,于2009年发起并组织了类似的高升力计算研讨会,其第一次会议(HiLiftPW-I)已于2010年6月在美国举行,文献[60]是该次会议的总结。在上述工作的基础上,2013年1月AIAA又进一步在其ASM会议过程中形成了以加拿大McHill大学Nadarajah教授为首的空气动力优化设计讨论组,作为空气动力优化设计计算系列研讨会实际的组委会。讨论组讨论了:1)建立可供在一个有约束的设计空间中测试气动优化方法的一组标准算例。2)举行研讨会的时间。与会者一致认为,由于工业界对基于CFD的气动外形数值优化方法有强烈的需求,优化方法和工具的研制也已有了相当的发展,可以以类似于DPW的研讨会形式,通过对一系列复杂气动外形的优化,来评估现有的寻求最小阻力外形的各种优化方法的能力,并将结果向工业界/研究机构公布。与会者还认为第一次会议从二维和三维机翼外形开始是合适的,并请加拿大的与会者准备标准算例。第一次会议拟于2013或2014年的AIAA应用空气动力会议期间举行。

2.先导性的研究

事实上为准备此研讨会,波音的Vassberg,斯坦福的Jameson,以色列的Epstein及Peigin等三方从2007年起即开始了先导性的研究(pilotproject),以积累经验和发现问题。三方用各自己开发的优化软件(MDOPT,SYN107,OPTIMAS)对第三届DPW会议的测试机翼DPW-W1独立地作优化计算[61,62]。波音研制的MDOPT[63](也可参见文献[13]的1.7.3节)可使用响应面模型(InterpolatedRe-sponseSurface—IRS)的数值优化格式[64],也可直接从流场解计算设计变量的灵敏度代替IRS模型完成优化。其流场解软件为TLNS3D[65],计算网格点为3582225。Jameson开发的SYN107采用基于梯度的“连续”伴随方程方法[23,31],其流场解软件为FLO107,计算网格点为818,547。

以色列航空公司开发的OPTIMAS采用降阶模型的GA算法,流场解软件为NES[66-68],计算网格点为250,000。对三方独立优化后所得的外形再用不参与优化的流场解软件OVERFLOW[69]作评估计算,计算网格点数为4,000,000,以便能准确地计算阻力。结果表明,4个分析软件计算得到的阻力增量值的分散度在Ma=0.76时为5counts(1count=0.0001),Ma=0.78时为10counts,因此很难确定哪个优化后外形最优。但从Ma=0.76,C=0.5单设计点的阻力改进结果(表1)[61]看,OPTIMAS优化后的04外形明显优于MDOPT优化后的M5和SYN107优化后的S4。文献[61]还讨论了从比较中可吸取的经验和教训。

一种基于高可信度CFD模型的数值

优化方法的构造本节将以OPTIMAS为例对如何满足可应用于工程实践的高可信度CFD模型数值优化方法的要求做一说明。

1.优化方法的构造及其特点

OPTIMAS是将遗传优化算法和求解全N-S方程的分析算法相结合的一种有效并鲁棒的三维机翼优化方法。1)其全N-S方程的流场并行解算器NES[66-67]基于高阶低耗散的ENO概念(适用于在多区点对点对接网格中的多重网格计算)[66,71]和通量插值技术相结合的数值格式,采用SA湍流模型,可快速准确地完成气动力计算,因此具有计算大量不同流动和几何条件的鲁棒性。作为例子图3和4给出了ARA翼身组合体Ma=0.80,Re=13.110时的升阻极线和CL=0.40时的阻力发散曲线[68],使用的网格点数分别为,细网格(3lev):900,000,中等网格(2lev):115,000。可见升阻极线直到大升力状态的计算与实验都很一致。对比图中还给出的TLNS3D在细网格(2,000,000)中的计算值可见,无论升阻极线或阻力发散曲线NES的都更优。作为数值优化软件的特点之一是其在流场解算器中首次使用了高精度格式。2)优化计算的遗传算法中采用了十进制编码、联赛选择算子[42]、算术交叉算子、非均匀实数编码变异算子[72]和最佳保留机制。为解决搜索时总体寻优耗时大和求解N-S方程估算适应函数代价高的问题,在寻优过程中估算适应函数时采用当地数据库中的降阶模型[54,58]获取流场解(当地数据库是在搜索空间中离散的基本点处求解全N-S方程建立的),并以多区预测-修正方法来弥补这种近似带来的误差。多区预测-修正方法即在搜索空间的多个区域并行搜索得到各区的优化点,再通过求解全N-S方程的验证取得最优点。为保证优化的收敛,寻优过程采用了迭代方法。3)在整个空间构筑寻优路径(图5),扩大了搜索空间和估算适应函数的区间[54]。4)为提高计算效率,OPTIMAS包含了五重并行计算:Level1并行地求解N-S方程Level2并行地扫描多个几何区域,提供多个外形的适应函数的计算(level1隐于level2中)。Level3并行的GA优化过程(level3隐于level4中)。Level4并行地GA搜索多个空间。Level5并行地生成网格。5)采用单参数或双参数的BezierSpline函数对搜索空间参数化;并基于优化外形与原始外形的拓扑相似自动地实现空间网格的快速变换。

2.优化设计的典型结果

文献[53]~文献[58]给出的大量算例充分表明了OPTIMAS优化软件的优异性能。本文5.2中给出了其优化三维机翼的性能,这里再补充两例。1)翼身组合体整流(fairing)外形的优化文献[73]讨论了某公务机翼身组合体机翼外形优化的单点和多点设计两者性能的比较。结果表明,多点优化设计能同时保证设计的巡航状态时,和高Ma数飞行,起飞等非设计状态时的良好性能。文献[74]进一步讨论了翼身组合体整流外形的优化设计。流动的复杂性(三维粘流/无粘流强相互作用的流动区域)和几何的复杂性(三维非线性表面)使整流外形的设计经历了传统的试凑法,基于Euler解的试凑法等,最后才发展为现代完全N-S解的数值优化方法。文献[74]采用了这种方法,先作机翼外形优化,再作整流外形优化,然后再作机翼优化,整流外形优化,……依次迭代,直至收敛。优化中用双参数的BezierSpline函数将整流外形参数化,所得搜索空间的维数ND=(2N-2)*(M-1)决定的参数化整流外形与实际外形的差别在M=10,N=4,ND=54时可准确到0.3mm(满足工程需求)。计算网格数为90万。表3给出了设计条件和约束,表4给出了设计点的阻力值比较。由表4可知,GBJ2的减阻为16.7,50%DC,GBJFR1的减阻为10.7,32.1%DC,GBJFR2的减阻为5.9,两次优化机翼的减阻总计为22.6,67.9%DC,优化机翼和优化整流外形减阻作用分别约占2/3和1/3,可见整流外形的优化也是十分重要的。约束则使减阻损失4.6(如GBJFR3-GBJFR1)。图6至图9分别为原始外形,GBJ2,GBJFR2和GBJFR4的整流处等压线分布,可见整流外形的优化消除了原始外形和GBJ2中存在的激波。图10和图11分别给出了Ma=0.8时升阻极线和CL=0.4时阻力发散曲线的比较,可见优化设计不仅对设计点,对非设计状态也都有好处。2)翼身融合体飞机气动外形的优化[75]优化设计以英国克朗菲尔德大学设计的BWB外形[76]为出发外形,该外形的主要设计点为,。在数值优化计算中还考虑了,的第二个设计点和,(起飞状态)的第三个设计点。几何约束有剖面相对厚度,前缘半径,后缘角,每个剖面的樑处还附加两个厚度约束,其中上标b表示出发外形,*表示优化外形,下标i表示第i个剖面。附加的空气动力约束为对俯仰力矩的规定。采用Bezier样条描述几何外形,总设计变量为93个。表5给出了设计计算各状态的条件和约束,其中是权系数。表6给出了优化计算结果。#p#分页标题#e#

单点优化的BWB-1结果与文献[77]的结果相比较可见,文献[77]采用Euler方程的无黏优化使阻力降低了26counts;而这里的BWB-1全N-S方程优化使阻力降低了52counts,显示了此黏性优化方法的优点。比较有、无俯仰力矩约束时优化得到的BWB-2和BWB-1表明,尽管BWB-1阻力降低的效果突出,但其值过大,出于稳定性考虑而不能接受;BWB-2的阻力虽比BWB-1大了1.9counts,却满足了力矩的要求。表6中的双点优化设计(BWB-4),使第三设计点(低速状态)的达到1.671(消除了BWB-2达不到设计要求1.63的缺点),且基本保持了主设计点的阻力收益,为196.6。然而BWB-4在时的阻力达216.6,高于BWB-2的213.4,表明需要三个设计点的优化设计(BWB-3)。BWB-3在时,为202.5(比两点设计值减小了14.1),同时满足了其它两个设计点的性能要求。图12至图15给出了所有设计状态和时的极曲线,时的阻力发散曲线和时的随迎角α变化的曲线。由图可见,时所有优化设计的极曲线都非常接近,相比于原始外形的极曲线,性能有了很大改进;时也一样,特别是三点优化设计的BWB-3,优点更明显。阻力发散曲线也都有了很大改进,在前所有的总阻力基本保持常值,单点与两点优化的阻力发散点接近,而三点优化的可达附近。由图15可知,没有考虑低速目标的BWB-1和BWB-2具有较低的,将低速目标计入设计状态的BWB-3和BWB-4所得的皆优于原始外形的。上述结果表明三点优化设计具有最佳的优化效果和总体最好的气动性能。Fig.15LiftcoefficientCLvsangleofattackatMa=0.2最后,上述各优化结果在(主设计点)时的阻力值基本相同,但几何外形却差别不小,由此可见,外形阻力优化问题没有唯一解[75]。上述计算是在具有456GBRAM,114MB二级高速缓存的机群环境下通过“过夜”方式完成单点优化设计,在1.5-2天的计算时间内完成三点优化设计的,计算时间可满足应用于工程设计的要求[75]。

结束语

本文概要地叙述了数值优化方法在民机发展中的应用历史和现状;介绍了即将举行的空气动力优化设计计算系列研讨会;重点讨论了可应用于民机气动设计的基于高可信度CFD模型的数值优化方法,对其要求及其构造方法。以算例表明了这种方法不仅可用于传统圆筒机身+机翼的民机外形,也可用于非常规布局的民机外形(BWB飞机),且从所需计算机资源和计算时间看,可用于民机日常的工程设计。我们如能尽快掌握并应用这种方法无疑将大大缩短我国民机的研制周期和节约研制成本,有利于我们迅速赶超世界先进水平。相信世界范围的空气动力优化设计计算系列研讨会必将进一步推动数值优化方法的发展和应用。(本文图、表略)

本文作者:朱自强 单位:北京航空航天大学 航空科学工程学院