考虑三维流动效应的自然层流短舱压力分布反设计

刘红阳，宋超，罗骁，周铸，吕广亮

中国空气动力研究与发展中心 计算空气动力研究所，绵阳 621000

经济性和环保性一直是民航客机永恒不变的追求之一，气动阻力与之紧密相关。客机各部件的阻力一般分为压差阻力和摩擦阻力，湍流状态下的摩擦阻力远高于层流状态下的摩擦阻力，因此使客机表面保持大部分层流区域将会大大减小客机的总阻力。早期已有研究表明，如果在机翼和尾翼保持一定的层流区，总阻力将降低30%［1］。然而机翼后掠角较大，且对升力和力矩有一定要求，自然层流机翼的设计较难实现。

随着涡扇发动机技术的发展，大涵道比发动机以其噪声低、燃油效率高、使用寿命长等优势被更多的应用于民航客机。一般尺寸的民用涡扇发动机短舱的阻力占客机总阻力的5%左右，随着发动机涵道比的不断增大，短舱的尺寸随之增大，短舱阻力的占比也不断提升。由于飞机设计中对短舱的升力和力矩等特性不作要求，且没有横流不稳定和附着线转捩等问题，自然层流短舱的设计更容易实现，进而达到减小阻力的目的［2］。

国外学者对层流短舱的研究工作开展的比较早。从20 世纪80年代中期开始，美国NASA兰利研究中心开展了自然层流短舱的风洞和飞行试验研究［3］。1984年，Younghans 等［4］对设计的单独自然层流短舱和翼吊自然层流短舱开展了风洞试验验证研究，获得了较大范围的层流区域。1990年，Radespiel 等［5］研究了层流短舱设计的可行性，使用数值流动方法实现了层流短舱设计，在巡航飞行状态下自然层流长度占短舱长度的60%。1998年，德国宇航院DLR、罗尔斯·罗伊斯公司和Motorenund Turbinen-Union Friendrichshafen（MTU）合作开展了发动机短舱自然层 流 和 混 合 层 流 飞 行 试 验［6］，在VFW614/ATTAS 飞机上获得了短舱长度60%的层流区。2011年，Lin［2］对自然层流短舱开展了气动设计、风洞试验和转捩预测方法的验证研究。2014年，Vermeersch 等［7］通过求解欧拉方程并分析边界层稳定性，对自然层流短舱和混合层流短舱开展了研究，分析了除冰和加热等对转捩的影响。随着自然层流短舱技术的逐步发展和工业制造水平的不断提升，波音公司已成功将其应用于波音787 系列客机上。

国内学者对自然层流短舱的研究起步较晚。2014年，何小龙等［8］基于EFFD 方法实现了轴对称自然层流单独通气短舱和带动力短舱的优化设计，获得了短舱长度40%以上的层流覆盖区。2016—2017年，Zhong 和Li 发展了短舱截面造型设计方法，根据二维型线生成三维短舱并进行评估和优化［9-10］。2019年，杜玺等［11］开展了自然层流短舱的气动设计、数值计算和风洞试验验证研究，获得了短舱长度30%～55%的层流区，Wang和Sun 等［12］使用有限差分算法对跨声速自然层流短舱进行了优化设计研究，层流区面积增加16.64%，总阻力系数减小11.6 counts，孟晓轩等［13］基于线性稳定性分析方法，将双eN方法同RANS 方程求解器耦合，研究了来流马赫数、雷诺数、湍流度以及攻角对短舱转捩的影响。2020—2021年，胡骁等［14-15］通过设计轴对称层流短舱基准面，并二次开发CATIA 工具，实现了三维非轴对称自然层流短舱设计。

实际上短舱属于非轴对称构型，受机翼机身的气流影响，具有复杂的流动状态。传统的设计方法以二维或二维假设的形式开展，没有考虑短舱较强的三维流动效应，这种方式会损失设计过程中获得的层流区面积；如果直接以非轴对称短舱为对象进行优化设计，则面临设计变量数目多、CFD 求解次数多、计算成本高的问题。

本文提出一种基于生成拓扑映射（Generative Topographic Mapping，GTM）［16-17］的压力分布反设计方法，充分考虑复杂三维流动效应，直接对短舱三维外形进行设计，通过对短舱外形及压力分布组成的数据集进行降维处理，有效减少所需样本数目，且设计过程中无需反复计算流场，能够大大提高设计效率。

结合自由曲面变形技术（Free Form Deformation， FFD）［18］、RANS 方程、基于SST（Shear Stress Transport）湍 流 模 型 的γ-转 捩 模型［19］、基于径向基函数插值的网格变形方法［20-21］、GTM 模型和遗传算法，构建了考虑三维流动效应的自然层流（Natural Laminar Flow，NLF）短舱优化设计方法，利用该方法对单独通气NLF 短舱进行了优化设计研究。

1 基于生成拓扑映射的压力分布反设计

1.1 生成拓扑映射模型

生成拓扑映射是一种非线性隐变量模型，能够实现高维数据到低维空间的映射，并保持原有的拓扑相对关系，目前广泛应用于数据分析及可视化分析。

设 有D维 数 据 集T={t1，t2，…，tN}，但 数 据本质是L维的（L＜D）。映射函数y(x，W)将L维空间的K个隐变量点x={x1，x2，…，xK}映射到D维数据集T。采用中心位于y(x，W)且方差为β的高斯函数作为D维数据集T的分布：

式中：p(t|x，W，β)为条件概率密度函数，对x积分可得：

其中：p(x)为隐变量空间维度的先验概率密度函数：

表示均匀分布在隐变量空间内的K个点上。

由此得到：

通过最大对数似然函数可求得参数W和β：

分布p(t|W，β)为中心固定的带约束混合高斯模型。GTM 算法可归纳为：已知隐变量空间中x的分布和数据空间的数据集合T，使用最大似然法求解参数W和β。常采用期望最大化（Expectation Maximization，EM）算法训练GTM 模型参数［22］。

1.2 NLF 短舱反设计流程

图1 定义了非轴对称短舱的4 个基准面，以短舱的顶部截面为0°基准面，沿顺时针方向依次为90°基准面、180°基准面和270°基准面。

图1 非轴对称NLF 短舱各个基准面的定义Fig.1 Reference planes of non-axisymmetric NLF nacelles

相隔180°的2 个基准面之间的层流区相互影响很小，因此可以先对其中一组对向的2 个基准面分别进行反设计，然后在此基础上设计另一组对向的2 个基准面。2 组基准面的设计顺序不影响设计结果。图2 给出了NLF 短舱的设计思路：

图2 非轴对称NLF 短舱的设计思路Fig.2 Design idea of non-axisymmetric NLF nacelles

步骤1在短舱初始外形baseline1 的基础上，分别对90°基准面和270°基准面进行压力分布反设计，获得相应的短舱外形design1.1 和design1.2。

步骤2应用FFD 技术，将design1.1 的90°基准面和design1.2 的270°基准面匹配到baseline1 外形上，获得新的短舱初始外形baseline2。

步骤3在baseline2 的基础上，分别对0°基准面和180°基准面进行压力分布反设计，获得相应的短舱外形design2.1 和design2.2。

步骤4应用FFD 方法，将design2.1 的0°基准面和design2.2 的180°基准面匹配到baseline2外形上，获得最终的短舱反设计外形design3。

图3 给出了NLF 短舱各个基准面的压力分布反设计流程，具体如下：

图3 NLF 短舱各个基准面的压力分布反设计流程Fig.3 Pressure distribution inverse design process of NLF nacelle reference planes

步骤1针对短舱外形和基准面上的设计变量进行实验设计，采用拉丁超立方方法选取样本点S，并计算对应的压力分布Cp。

步骤2根据高维数据集{S，Cp}创建并训练GTM 模型，对数据集进行降维处理，获得相应隐空间映射关系。

步骤3结合预设目标压力分布，采用全局优化算法在隐空间寻优，GTM 模型输出预测的最优设计变量与压力分布。

步骤4对最优设计变量进行CFD 校验，如果相应压力分布和层流区长度满足设计要求，则设计结束，否则需要重新设定目标压力分布或更新样本数据集，重复以上步骤。

基于GTM 模型的压力分布反设计有如下优点：

1） 目标压力分布的设定方式灵活。不必给出绕短舱基准面一周的完整压力分布，按照设计需求仅给定基准面上一段或多段压力分布即可进行设计。

2） 寻优效率高。GTM 模型将短舱构型及其压力分布组成的高维数据集进行降维处理，全局优化算法在低维隐空间寻优，极大提高了优化效率。

3） 计算成本低。GTM 模型通过样本训练，将高维数据集映射到低维隐空间，这种高精度映射关系避免了CFD 求解器的反复调用，缩短了设计时间，降低了计算成本。

2 NLF 短舱反设计算例

初始的非轴对称NLF 短舱模型尺寸参数如图4 所示，短舱外罩轴向最大长度Lmax=5.544 m，周向最大直径Dmax=3.719 m，最大直径的轴向位置L1=1.865 m，下垂角Φ1=2.49°，外罩0°和180°基准面的船尾角分别为θ1=16.11°和θ2=13.75°，短舱外罩出口直径Doutlet=2.642 m，后缘厚度thicknessoutlet=0.004 m。对短舱的0°，90°，180°和270°基准面（如图1 所示）进行压力分布反设计，目的是推迟边界层转捩的发生，扩大自然层流区面积，进而减小短舱的摩擦阻力。

图4 初始短舱模型尺寸参数Fig.4 Size of nacelle baseline model

短舱模型的计算网格如图5 所示，网格节点数为633 万，边界层内第1 层网格高度为1×10-5m，边界层内网格增长率为1.2。设计状态为高度H=11 km，自由来流马赫数Ma=0.85，雷诺数Re=3.1×107，来流初始湍流度Tu=0.1%，黏性比Rt=5.0，攻角α=5.0°，侧滑角β=0°。

图5 非轴对称层流短舱计算网格Fig.5 Computation grid of non-axisymmetric nacelles

采用自由曲面变形方法对短舱进行参数化，短舱各个基准面上的设计变量分布如图6 中实心圆点所示，设计变量数目为8，设定降维后隐空间维度L=2，实验设计选取的样本数目为30。通过求解RANS 方程和γ-转捩模型，获得短舱压力分布，实现边界层转捩的预测。采用遗传算法在隐空间寻优，种群数目取为200，最大迭代步数取为50 步，目标函数为短舱基准面的实际压力系数与目标压力系数分布之间的差异：

图6 短舱各个基准面的设计变量分布Fig.6 Distribution of design variables on nacelle reference planes

式中：M为压力系数分布的点数；Cip为第i个点的压力系数值；CipT为第i个点的目标压力系数值。

2.1 短舱90°和270°基准面反设计

顺压区是保持层流的前提，顺压区长度和顺压区内压力梯度是影响层流发展的2 个关键因素，这种影响关系还与实际的流动状态紧密相关。一般情况下，顺压区越长，层流区越长，顺压区内压力梯度越大，越有利于推迟边界层转捩的发生，延长自然层流区长度，这2 个因素相互耦合，共同作用。

在短舱初始外形baseline1 基础上，对短舱90°和270°基准面进行压力分布反设计。图7 给出了短舱90°基准面反设计结果，其中c为基准面当地弦长。图7（a）中黑色圆点表示基准面弦长10%～30%区域内给定的目标压力分布，目的是增大顺压区内的压力梯度，红色实线代表短舱设计外形design1.1 经CFD 校验后提取的90°基准面压力系数分布，其趋势向目标压力分布靠近，由于给定的目标压力分布不一定具有真实的物理意义，即不存在对应的物理解，因此也不要求设计所得的压力分布与目标压力分布完全吻合。图7（b）对比了短舱90°基准面设计前后的摩阻系数分布和构型，此处以摩阻系数作为判断边界层转捩的依据，可以看出，90°基准面的层流区长度从当地弦长的23%延长到26.8%。图8 是设计过程中真实压力与目标压力分布的差异Ep收敛示意图，蓝色三角形表示每一代所有样本的Ep，红色实心圆点线图表示每一代样本的最小Ep。可以看出，在第15 代时设计基本收敛，Ep＜0.01。

图7 短舱90°基准面反设计结果Fig.7 Design results of nacelle reference plane 90°

图8 实际压力与目标压力分布的差异收敛示意图Fig.8 Convergence of error between real and target pressure distribution

如果进行直接调用CFD 求解器的全局优化设计，其计算时间与迭代步数、样本数目成正比，即使所有样本并行计算，设计周期仍以天为单位。本文采用的反设计方法优化过程不调用CFD 求解器，因此寻优效率很高，在样本计算完成后，只需要分钟量级的时间就可获得设计外形，一轮优化在个人工作电脑（CPU 主频3.4 GHz）只需要50 s。

图9 给出了短舱270°基准面反设计结果。图9（a）中黑色圆点表示基准面弦长10%～30%区域内给定的目标压力分布，目的是增大顺压区内的压力梯度，红色实线代表短舱设计外形design1.2 经CFD 校验后提取的270°基准面压力系数分布，其与目标压力吻合较好；图9（b）对比了短舱270°基准面设计前后的摩阻系数分布和构型，此处以摩阻系数作为判断边界层转捩的依据，可以看出，270°基准面的层流区长度从当地弦长的25%延长到29.3%。

图9 短舱270°基准面反设计结果Fig.9 Design results of nacelle reference plane 270°

利用FFD 参数化方法，将短舱外形design1.1 对应的90°基准面构型和design1.2 对应的270°基准面构型匹配到初始外形baseline1 上，获得新的短舱初始外形baseline2，并利用CFD对其进行校验。图10 为计算得到的压力系数分布和摩阻系数分布，并与baseline1 外形的相应数据进行对比，由于baseline1 外形是左右对称的，因此其90°和270°基准面的压力分布基本重合，摩阻分布略有差异，但层流区长度基本一致。baseline2 外形在90°和270°基准面上均保持了与design1.1 和design1.2 一致的较长层流区，此时baseline2 左右不再对称，且270°基准面的顺压区更长，因此层流区更长。

图10 baseline2 的CFD 校 验 结 果Fig.10 Validation of baseline2 with CFD solver

2.2 短舱0°和180°基准面反设计

在新的短舱初始外形baseline2 基础上，对短舱0°和180°基准面进行压力分布反设计。图11 给出了短舱0°基准面反设计结果，与其他3 个基准面相比，由于0°基准面处在背风区，其外表面顺压区很短，边界层转捩发生的比较靠前，这里给定基准面弦长5%～30%区域内的目标压力分布，如图11（a）中黑色圆点所示，目的是减小该区域内的压力梯度，增长顺压区，红色实线代表短舱设计外形design2.1 经CFD 校验后提取的0°基准面压力系数分布，其趋势向目标压力分布靠近；图11（b）对比了短舱0°基准面设计前后的摩阻系数分布和构型，此处以摩阻系数作为判断边界层转捩的依据，可以看出，0°基准面的层流区长度从当地弦长的15.8%延长到20.3%。

图11 短舱0°基准面反设计结果Fig.11 Design results of nacelle reference plane 0°

图12给出了短舱180°基准面反设计结果，图12（a）中黑色圆点表示基准面弦长10%～30%区域内给定的目标压力分布，目的是增大顺压区内的压力梯度，红色实线代表短舱设计外形design2.2 经CFD 校验后提取的180°基准面压力系数分布，其趋势向目标压力分布靠近；图12（b）对比了短舱0°基准面设计前后的摩阻系数分布和构型，此处以摩阻系数作为判断边界层转捩的依据，可以看出，180°基准面的层流区长度从当地弦长的28.3%延长到40.5%。

图12 短舱180°基准面反设计结果Fig.12 Design results of nacelle reference plane 180°

2.3 层流短舱设计外形校验评估

利用FFD 参数化方法，将短舱外形design2.1 对应的0°基准面构型和design2.2 对应的180°基准面构型匹配到外形baseline2 上，获得最终的非轴对称短舱设计外形design3，并利用CFD 对其进行校验。图13 对比了不同视图下短舱设计前后的摩阻系数云图，可以看出，与初始外形baseline1 相比，设计外形design3 各个基准面上的层流区长度均有不同程度的延长，其中180°基准面的延长量最大，0°，90°和270°基准面的延长量相当，这与各个基准面所处的当地流动状态有关，还与其样本及目标压力分布的差异性相关，不同的样本集和目标压力分布驱动优化算法朝着不同的方向进化，因此取得的优化效果也不同。此外，由于0°基准面处在背风区域，流动特性决定其无法获得较长的顺压区来保持层流。总的来说，本文所提方法设计的非轴对称短舱能够增大外表面自然层流区域，与设计过程中的各阶段外形相比，不会损失最终成型后的短舱层流区面积。

图13 不同视图下短舱设计前后摩阻系数Cf云图对比Fig.13 Comparison of Cf contour between design3 and baseline1 with different views

3 结 论

本文提出了一种基于GTM 模型的压力分布反设计方法，利用该方法对非轴对称单独通气短舱进行了自然层流优化设计，得到以下结论：

1） GTM 模型能够对气动外形及压力分布组成的高维数据集进行有效的降维处理，仅需较少的样本点数即可实现精确的压力分布反设计，且无需迭代计算流场，所提设计方法具有高效便捷的优点。

2） 传统的非轴对称NLF 短舱设计在二维或二维假设下开展，并在三维情形下校验设计结果，由于未充分考虑复杂三维流场的影响，使得最终短舱的层流区域有所损失。提出的设计方法计算成本较低，直接在三维外形上进行设计，过程中充分考虑三维流动效应，设计所得非轴对称短舱能够保持理想的层流区域。

3） 所提设计方法具有较高的工程应用价值，在下一步研究中，拟深入探究压力分布对边界层转捩的影响规律，将该方法应用于带动力单独NLF短舱和考虑机翼机身干扰的NLF 短舱的设计中，进一步验证该优化设计平台的鲁棒性和实用性。