航空发动机内流全场流动的大涡模拟

航空发动机是一个前后流通、各部件耦合工作的整体。在常规的发动机设计过程中，首先采用低维、时均的计算方法对发动机总体性能进行分析，获得压气机、燃烧室、涡轮等部件在交界面处的气动热力参数，这些参数通常以一维集总的“点”或二维分布的“线”的形式表征；然后根据交界面的气动热力参数要求，通过边界条件赋值将整机分解为各个独立的研究对象，分别开展压气机、燃烧室、涡轮三大核心部件的设计分析与试验测试；完成后再装配成核心机或整机开展部件匹配调试与整机性能测试。

在各部件独立开展的仿真分析中，作为边界条件的气动热力参数空间维度低且不随时间发生变化，无法准确反映贯穿于部件之间并随发动机运行瞬变的三维、非定常特征。例如，压气机喘振造成的燃烧室进口流量脉动与不稳定燃烧、燃烧室出口旋流高温热斑在涡轮叶片中造成的非均匀热负荷等。这种由于部件“边界效应”而人为引入的缺陷，使得各部件计算获得的流场与真实发动机整机流场存在偏差,极大地影响了部件独立预测分析结果与整机试验结果的一致性。

笔者及研究团队以航空发动机核心机为研究对象，采用湍流燃烧大涡模拟技术对高压压气机、燃烧室、高压涡轮耦合通流流场开展数值仿真，介绍了几何模型、数值计算模型中的一些要素，并分析了整机内跨部件传播的非定常特征以及上下游部件的相互影响效应，以期初探整机内流场数值仿真技术的可行性及其在工程应用中的技术价值。

从20世纪60—70年代开始，工业界逐渐采用计算流体力学（CFD）方法对航空器外流进行仿真，并将计算结果作为工程设计的参考[1]。进入21世纪后，科学研究进一步重视湍流的非定常效应，试图通过优化或控制流动过程以提升航空器内、外流的气动热力性能。因此，一种对湍流流动时间、空间尺度均足够精确的大涡模拟方法（LES）在业界逐步推广，成为目前分析航空发动机内部气动热力特征的先进工具之一。

然而，航空发动机各部件之间的气动热力状态，包括温度、压力、马赫数（Ma）、雷诺数(Re)等，差异极大，导致多部件耦合的气动热力仿真除了要具备宽速域、可压缩的求解方法外，还须结合实际物理特征，建立恰当的数学模型。在叶轮机械中，叶片表面边界层转捩、分离以及通道中二次流、端壁间隙流是主要流动现象，因此数值仿真中须建立恰当的湍流模型与近壁面条件；在燃烧室中，大尺度旋流、剪切与回流用于强化燃料与空气掺混与稳定火焰，因此数值仿真中须充分评估流动、混合与化学反应时间尺度的差异，建立微尺度下流动与燃烧耦合作用的燃烧模型。当前，上述主要计算方法在各部件的独立仿真中均有着长足发展、日趋成熟。例如，法国欧洲科学计算研究中心（CERFACS）在2009年开展了环形燃烧室大涡模拟[2]，在2019开展了3级压气机的大涡模拟[3]。

进入21世纪以来，为进一步提高发动机整机内流的认识，科学研究率先尝试进行了发动机整机气动热力流场的仿真。2003—2006年，斯坦福大学针对PW6000整机内流开展仿真计算[4]，在其研究中，采用可压缩的雷诺时间平均方法（URANS）模拟压气机和涡轮内流，采用不可压缩大涡模拟方法模拟燃烧室流动。这一尝试在当时是突破性的技术研究，但是由于需要在旋转部件与燃烧室之间进行仿真方法的切换，导致部件之间的湍流特征时间尺度并不一致，因此该工作所开展的多部件耦合仿真，只是几何流道耦合，而不是流场的物理过程耦合。

最近10年，大规模高性能并行计算技术的快速发展为发动机整机耦合仿真带来了新的契机。2020年，CERFACS采用20亿网格单元、14400核计算节点、仿真分析了DGEN380商用飞机发动机内部的流动与燃烧过程，并揭示了贯穿整机内流场的非定常现象[5]。

笔者及研究团队采用一核心机计算模型（见图1）进行压气机、燃烧室和涡轮的耦合仿真。其中压气机为10级高压压气机、燃烧室为头部中心分级的短环形燃烧室、涡轮为两级高压涡轮。在整机模型中，燃烧室周向共有18个头部，每个头部对应的周向角度为20°。为节省计算资源，数值计算域周向角度取值20°，这一周向角取值可包含一个完整的单扇区燃烧室。为保证非定常计算中叶轮机械转子叶片、静子叶片交界面尾迹传递不受周向简化的影响，采用了稠度不变的原则，将原方案中每排压气机、涡轮叶片数调整为燃烧室头部数量的整数倍。因此，本文采用的数值计算域可在发动机周向进行周期性拓展。

图1  核心机计算模型示意

在压气机与涡轮的几何建模中，保留由叶型、端壁、机匣所形成的主流流道，省略级间引气、盘腔引气、端壁泄漏、气膜冷却等次流流道与对应空气流量的变化。在燃烧室的几何建模中，保留火焰筒内、外环腔流道，火焰筒内流道以及伸入燃烧室内的燃油喷杆，省略火焰筒气膜孔、火焰筒悬挂、支撑结构等几何细节，省略了燃烧室扩压器、火焰筒内环腔的引气流路与对应空气流量的变化。但是，在燃烧室实际工作中，燃烧室内、外环腔空气经火焰筒壁面冷却孔流入火焰筒内流道，汇合头部燃气流向高压涡轮。该部分冷却空气量占燃烧室总气量的30%以上，因此不可忽略。在研究中，采用自行开发的火焰筒壁面源项法模拟燃烧室内、外环腔流入燃烧室的冷却空气[6]（见图2）。

图2  壁面源项法对火焰筒壁面气冷结构的建模示意

整个计算域（见图3）进口为高压压气机支板进口，设置均匀分布的总温、总压，并忽略进口空气的湍流脉动特征。计算域第一处出口为燃烧室火焰筒内、外环腔出口，根据计算出的燃烧室流量分配设置环腔流出流量，并以此流量作为火焰筒壁面源项法中的输入气量添加入火焰筒内。计算域第二处出口为高压涡轮第二级转子叶片出口，边界条件设置为静压出口。作为燃料的航空煤油以液态形式通过燃烧室内若干离散的喷嘴进行喷射。

图3  包含压气机、燃烧室与涡轮的多部件耦合流道示意

在压气机、涡轮的静子叶片流道与燃烧室中，除火焰筒内侧壁面外，其余固体壁面为绝热无滑移壁面。在压气机、涡轮的转子叶片流道中，设置转子叶片局部流道为旋转域并根据实际转速设置域的旋转速度。叶轮机内流场采用正交网格进行划分。燃烧室头部旋流器、火焰筒内流道采用正交网格进行划分，二者之间的过渡区域采用非结构网格进行划分，计算域网格总数约为3000万个。

为兼顾压气机、涡轮中的高马赫数流动以及燃烧室中的低马赫数流动，选择可压缩的流场求解方式。湍流模式采用Dynamic Smagorinsky-Lilly模型并结合标准壁面函数处理近壁面流动，气、液两相流模型采用项目团队自行开发的二次雾化模型（见图4），燃烧模型为增厚火焰面模型，化学反应机理为两步简化的航空煤油代用燃料（C10H23）与空气的化学动力学机理。计算域空间离散为二阶迎风格式并添加了限制器，时间步长为2×10-7s，求解器选择双精度。上述主要计算模型在项目团队的前期工作中，已通过基础试验数据完成了精度与适用性验证。

图4  完善后的燃油二次雾化模型

计算中，通过监测压气机、燃烧室内预置测点的速度脉动与涡轮出口流量判断流场的计算结果是否达到统计学稳定。计算稳定后，继续计算7个非定常脉动周期并在期间保存有效计算数据进行统计平均，随后终止计算并进行结果分析。数值计算工作在中国航发商发高性能计算平台HPCCII上开展。

图5为计算获得的压气机与涡轮内流中径叶高处熵分布与燃烧室中2000K的温度等值面分布。由图5可以看出，从压气机第一排转子叶片开始，每一级叶片均产生了明显的非定常尾迹，这一尾迹流入下一排叶片中，与该排叶片自身的尾迹相互叠加，逐级向下游发展。为定量表征非定常脉动特征在发动机中的传播，在计算域中取3个典型位置监测速度随时间的变化，分别位于压气机末级导向叶片内部、燃烧室扩压气出口、高压涡轮第一级导向叶片进口（如图5中P1、P2、P3位置）。

图5  整机内流场瞬时熵与燃烧室中的火焰面分布

图6为P1、P2、P3测点速度随时间的脉动特征。从图中可以看出，压气机末级导向叶片通道中，有着明显的非定常速度脉动，这是由于上游多级压气机转子所产生的。在前置扩压器出口，由于气流在自由通道中突扩，上游速度扰动被显著衰减，到燃烧室出口时，速度脉动幅值进一步衰减，但是在火焰筒出口仍能观测到上游非定常流动对下游的影响。

图6  发动机内3个典型测点气流速度随时间的变化特征

为分析多级压气机尾迹对燃烧室出口热斑的影响，对比计算了考虑上游压气机与省略上游压气机时的火焰筒温度分布特征（见图7）。可以看出，当计算中耦合上游压气机时，燃烧室出口高温区在周向、径向所占区域比忽略上游压气机所呈现出的结果更宽一些，燃烧室中时均的高温火焰面更紧凑。分析认为产生该现象的原因在于压气机转子尾迹扫掠对燃烧火焰带来了非定常扰动，强化了燃料与空气的混合程度，使得火焰面更为紧凑，也使得旋流火焰的扩张角度有所提高。

图7  燃烧室时均火焰面形状与出口温度分布云图

通过仿真分析可知，发动机上游产生的非定常流动会向下游传播，并影响下游部件的流动特征。通过初步研究可以看出，发动机内流场多部件耦合仿真，可消除部件之间的边界效应，完整地分析空气从发动机流入至流出的全部过程，并揭示之前无法分析的多部件内流耦合作用关系，有助于进一步提升各部件性能设计水平。同时，整机内流仿真也可在数字孪生、虚拟样机等技术中发挥作用，配合其他学科的仿真方法，进一步助推数字化发动机的实现。

（张漫，中国航发商发，研究员，主要从事航空发动机燃烧室设计与仿真）

参考文献

[1] Jameson A. The present status, challenges, and future developments in Computational Fluid Dynamics[J]. Agard fluid dynamics panel symposium, 1995.

[2] Laborderie J,Duchaine F,Gicquel L,et al. Wall-modeled large-eddy simulations of a multistage high-pressure compressor[J]. Flow turbulence and combustion, 2019.

[3] Staffelbach G,Gicquel L, Boudier G,et al. Large eddy simulation of self excited azimuthal modes in annular combustors[J]. Proceedings of the combustion institute, 2009, 32(2):2909-2916.

[4] Medic G , You D , Kalitzin G . On coupling of RANS and LES for integrated computations of jet engines[J]. ASME Turbo Expo 2007: power for land, sea, and air - Montreal, 2007(6):1181-1187.

[5] Perez C A,Dombard J,Duchaine F, et al. Large-eddy simulation of an integrated high-pressure compressor and combustion chamber of a typical turbine engine architecture - ASME Turbo Expo GT2020-16288[C]. ASME Turbo Expo,2020.

[6] Wang J, Mao R, Tao W, et al. Numerical simulation of combustor effusion cooling flow based on source term method[C].ASME Turbo Expo 2019: turbomachinery technical conference and exposition, 2019.