导读：航空发动机折流燃烧室几何结构复杂，其高保真数值模拟需要高效的网格与边界条件处理方法。采用曲线坐标系隐式浸没边界方法结合大涡模拟-概率密度函数输运方程湍流燃烧模型开发自研软件，并实现 WP11 中折流燃烧室的高保真模拟。流动模拟中准确解析了该燃烧室中的三股主要气流，且三股气流分别约占进口流量的75%、12.5%和12.5%。两相燃烧模拟中针对拉格朗日框架下的液滴运动和欧拉框架下的湍流燃烧采用不同网格标记，模拟得到的出口径向温度分布规律与实验一致，平均相对误差为17.95%，表明基于本方法开发的自研软件能准确模拟折流燃烧室中的两相湍流燃烧现象。 一、写在文前
离心甩油折流燃烧室广泛应用于中小型航空发动机中，利用中小型发动机转子的高转速通过离心甩油盘实现燃油雾化，具有供油压力低、雾化质量稳定、燃烧效率高和高空性能好等优点[1]。美国 TeledyneCAE 公司、Williams 公司及法国 Turbomeca 公司均生产了多种型号的采用折流燃烧室的航空发动机[2]，国内涡轴系列发动机也大量采用了折流燃烧室。航空发动机燃烧室几何结构复杂，其中的湍流燃烧现象存在强烈非线性关系，迫切需要高保真、高效率、高精度的数值模拟技术[3]。受限于真实折流燃烧室复杂的几何结构，目前对折流燃烧室的数值模拟研究大多基于采用非结构网格或在将大幅简化燃烧室几何后采用贴体结构网格。而随着近年来计算机硬件的发展，应用 IBM 和LES-TPDF 模型模拟真实折流燃烧室成为可能。AECSC-IBM（Aero Engine Combustor Simulation Code based on Immersed Boundary Method）是北京航空航天大学（北航）能源与动力工程学院航空发动机数值仿真中心（仿真中心）在 LESTPDF湍流燃烧模型算法基础[14-16]上结合曲线坐标系隐式 IBM[17]研发的适用于复杂结构燃烧室内高保真模拟的两相湍流燃烧数值模拟软件，其中包括几何模型导入、网格标记生成、边界条件模化以及两相湍流燃烧的求解等功能的实现均采用自研程序和算法。本文首先介绍该软件中网格标记的生成、IBM 结合 LES-TPDF 模型的实现以及针对折流燃烧室中甩油盘甩出的高速液滴粒子轨迹求解方法等算法，并基于该软件模拟某型折流燃烧室地面实验工况下的两相湍流燃烧。2 研究对象与计算方法2.1 研究对象基于 AECSC-IBM 软件模拟和研究某型真实结构的折流燃烧室，其结构与内部流动情况如图1 所示，经压气机压缩过的高压空气经轴向扩压器进入燃烧室，并分为三股流入火焰筒，参加燃烧反应和与燃气掺混，液相煤油由高速旋转的甩油盘进入主燃区。本文中所模拟该全环形燃烧室的1/10模型如图 2 所示。二、研究对象与计算方法
基于 AECSC-IBM 软件模拟和研究某型真实结构的折流燃烧室，其结构与内部流动情况如图1 所示，经压气机压缩过的高压空气经轴向扩压器进入燃烧室，并分为三股流入火焰筒，参加燃烧反应和与燃气掺混，液相煤油由高速旋转的甩油盘进入主燃区。本文中所模拟该全环形燃烧室的1/10模型如图 2 所示。Fig. 1 Schematic diagram of the structure and internal flow of the slinger combustor
Fig. 2 Geometric model of the slinger combustor为将 IBM 与 LES-TPDF 湍流燃烧模型结合，应用于折流燃烧室中两相湍流燃烧的数值模拟中，首先需要将初始的计算域网格（背景网格）的各个点区分和标记为流体域和固体域，从而将真实几何模型映射到计算域网格中，之后，基于网格标记对流体域边界的网格点设置边界条件，从而将真实几何模型的固体壁面转化为控制方程的边界条件。本文参考光线追踪中的常用的Möller-Trumbore 求交算法[18]，基于背景网格生成无限长扫描线，再与三角面网格格式（STereoLithography，STL）的几何模型中的所有三角面求交，得到对应于折流燃烧室模型的网格标记[17]。大涡模拟采用滤波函数过滤流场中的变量，直接求解尺度大于滤波尺度的变量并用亚网格模型模化尺度小于滤波尺度的变量，经过坐标变换和空间滤波后积分形式的守恒方程组为：
上式中，𝐺𝑘为控制体某方向上的流量，𝑆表示面积，𝑉表示体积。采用散度形式近似得到方程中的对流项。坐标系变换中非正交性产生的交叉导数项通过显式方法计算，并加入源项。通过前一时间步的质量守恒方程的解计算非线性项。以上所采用的离散方法产生了对于通用变量的准线性方程组：
其中，𝑎表示系数，𝜙̃表示密度加权滤波后的通用变量， 
表示六面体网格中当前网格点紧邻的六个网格点的对流扩散项，需要结合浸没边界网格标记求解。下标P表示当前网格点，下标𝛼可取S,N, W, E, L, R，表示与当前网格点相邻的网格点。𝑆P代表源项，包含所有不能用面通量表示的项，并且依赖于 
。本文基于曲线坐标系下的隐式浸没边界方法[18]，通过处理离散方程系数阵实现 Neumann 边界条件。根据浸没边界网格标记寻找紧贴固体网格的流体网格。若某一流体边界网格的南方向（S）网格为固体壁面，对于离散形式的动量方程，首先向其中加入彻体力源项:
中包含了未知量
，可在求解离散方程前令𝑎S = 0实现。这种处理相当于求解过程中恒为 0，即无滑移边界条件。此外，为修正壁面附近的速度分布，用壁面函数法求解并加入由壁面切应力形成的源项。对于包含流场中各物质浓度和混合物焓值在内的标量输运方程，其边界均为零梯度边界条件。可通过向离散方程中加入源项
，离散形式的标量输运方程变为：
中包含未知量和，可以通过在求解离散方程前先令 
，再令 
实现。这种处理使离散方程中原有的 
变为 
，相当于在求解过程中恒等于，即实现零梯度边界条件。为了求解各物质浓度和混合物焓的亚网格分布，需求解随机变量集 合 
的输运方程：其中，为第𝑚个欧拉随机场中的随机变量集 合，其样本空间为标量集 合 
， 
为对流扩散项，
为随机项， 
为小尺度混合项， 
为化学反应源项。随机变量不存在亚网格分布，其输运方程与标量输运方程相比多出了随机项和小尺度混合项，这两项以源项的形式加入到式（6）中。其中小尺度混合项是亚网格的，与浸没边界网格标记的分布无关，而随机项以网格间随机输运的形式形成显式源项。
其中 
表示第𝑚个 Weiner 过程 
的分量， 
为网格间随机输运速度，由 Weiner 过程定义，与流场无关。因此在定义随机输运时，若𝑖方向相邻网格被标记为固体域，则令 
。航空发动机燃烧室中，航空煤油经初始雾化形成大量煤油液滴，煤油液滴蒸发为气相煤油后与空气掺混，并在一定条件下发生燃烧反应。本文采用欧拉-拉格朗日两相模型模拟液滴粒子与气流的质量和动量相互作用。两相模拟中，液滴的破碎、蒸发以及与流体的作用力求解均为亚网格的，与燃烧室壁面的位置无直接关系。但在求解液滴粒子的运动轨迹时，则需要根据壁面位置和法向量等信息计算粒子与壁面碰撞后的反弹或黏附情况。当采用如图3所示的贴体网格时，壁面位置对应于计算域边界，而固体域一定位于计算域外。因此对于𝑡时刻某粒子的位置 
，只需要判断该粒子按当前速度运动在下一时间步将到达的位置 
是否位于计算域外，若位于计算域外，则粒子将在𝑡~(𝑡 + Δ𝑡)这段时间内碰撞到壁面。计算其反弹时可求解与 
关于壁面对称的位置 
，即为该粒子在(𝑡 + Δ𝑡)时刻的位置，并且将粒子速度向量 
设置为反弹前速度向量𝒖𝑡关于壁面的对称向量。Fig. 3 Schematic diagram of the scan line penetrating into the model via the edge
基于浸没边界法模拟燃烧室中的两相燃烧时，燃烧室几何模型被映射为网格上的流体或固体标记，用网格标记表示的固体壁面与真实壁面一般是不完全一致的，如图 4 所示。因此在用网格标记表示壁面时丢失了真实壁面的部分信息，这对于用欧拉法求解的流动与气相燃烧影响不大，但对于用拉格朗日法求解的液滴粒子运动，则影响了粒子反弹轨迹和速度的求解，并可能导致粒子与壁面碰撞判断错误。
Fig. 4 Schematic diagram of the relationship between the grid marker and the real wallwhen using IBM首先，假设粒子不能在单一时间步Δ𝑡内跨越被标记为固体域，若粒子在𝑡~(𝑡 + Δ𝑡)这段时间内碰撞到壁面，则粒子下一时间步的位置
应当位于固体域内，如图 5 所示。因此，只需判断所在的网格标记是否为固体标记，若是，则需求解粒子的反弹轨迹和反弹后的速度。
Fig. 5 Schematic diagram of particle bounce after collision with the wall when using IBM
粒子反弹的求解需要求解点 关于壁面的对称位置 
和速度 
关于壁面的对称矢量 
，因此需要真实壁面的法向量等信息。对于紧邻固体域的流体网格，重新调用前文中的求交算法，从而获取真实壁面的信息。将某一紧邻固体域的流体网格位置计为 
，生成三条正交的扫描线 
，其中𝑫1、𝑫2、𝑫3相互垂直。对于 STL 模型三角面集 合中第𝑚个三角面，几何文件中包含三角面的三个顶点坐标 
以及该三角面的法向量 
，因此将集 合中第𝑚个三角面记为 
。基于 Möller-Trumbore求交算法得到 
的值。因此点到三角面 
的距离 
可以表示为 
或 
或 
。遍历所有三角面并求出点到三角面的最小距离 
，则将 
作为点处的壁面法向量，用于求解粒子在位于处的网格附近的反弹轨迹和速度。折流燃烧室利用发动机转子带动离心甩油盘高速转动实现燃油雾化，在高转速工况下，甩油盘初始雾化得到的液滴速度往往很高。而折流燃烧室火焰筒为薄壁结构，采用 IBM 方法模拟其两相燃烧时，可能会出现粒子在单一时间步内穿越薄壁的情况，如图 6 所示。由于火焰筒壁对应的固体域较薄且粒子速度较高，当前位于流体域 
位置的某粒子在下一时间步的位置 仍位于流体域，但在𝑡~(𝑡 + Δ𝑡)时间段内的轨迹穿过了固体域。Fig. 6 Schematic of a particle crossing a wall in a single time step主燃烧室中，燃油液滴一般只在火焰筒中运动，因此可以对欧拉法求解湍流燃烧和拉格朗日法求解粒子轨迹采用不同的流固标记。对于液滴粒子的运动，可以将燃烧室模型作较大幅度的简化，只保留简化后的火焰筒内流体域作为粒子轨迹的计算域，从而减少求解壁面法向量等信息的计算量，并避免高速粒子在单一时间步内穿越薄壁结构，如图 7 所示，图中绿色阴影为液滴粒子轨迹求解时的固体壁面。Fig. 7 Schematic diagram of particle bounce after collision with thin wall
三、计算结果及分析
基于以上算法形成的 AECSC-IBM 软件，采用 349 万网格模拟 Sandia 实验室射流燃烧实验[19]中雷诺数分别为 33600 和 44800 的 Flame-E 与 Flame-F 两种射流火焰，以检验算法与软件应用于湍流燃烧模拟的求解精度。模拟中底部进口中心的半径3.6 mm的内喷口喷射293 K，1: 3的甲烷-空气混合气，半径9.1 mm外喷口喷射1880 K（Flame-E）或1860 K（Flame-F）的高温空气，模拟得到两种算例的瞬态温度云图如图 8 所示。
Fig. 8Transient temperature diagram of jet flame统计射流火焰的时均温度分布的模拟结果并与实验测点的温度数据对比，得到图 9 及表1，其中𝑟表示截面上的径向位置，𝑧表示截面在轴向上的位置，𝑑表示射流出口直径。Flame-E的时均温度分布的平均相对误差在 11.51%~22.38%，均值为 14.69%，Flame-F 时均温度分布的平均相对误差在 7.16%~21.23%，均值为 14.18%。射流火焰的模拟结果表明，上文所述算法与软件能够准确地模拟湍流燃烧现象，适合进一步应用于折流燃烧室中的湍流燃烧模拟。Fig. 9Time-average temperature distribution of jet flame along the radial direction
IBM 网格标记生成折流燃烧室的 STL（stereolithography）模型如图 10 所示，实际为多个三角面组成的面网格。划分计算域背景网格时，为了降低固体域在背景网格中的占比，采用外形与折流燃烧室计算域轮廓相近的曲线坐标系中的分块结构网格，如图 11 所示。
Fig. 10 Triangular surface grid of the refractory combustion chamber
Fig. 11 Z-axis central section background grid with parallel chunking schematic
按以上方法生成的网格标记在 Z 轴中央截面上的分布如图 12 所示。其中流体计算域被标记为红色，而固体域被标记为蓝色。对应于流体和固体网格标记边界的燃烧室几何壁面的三维示意图如图 13 所示，可知由网格标记表示的三维燃烧室模型与图 2 中的燃烧室模型一致，因此本文中的扫描算法能够准确地将折流燃烧室几何模型映射为网格标记。
Fig. 12 Z-axis central section background grid and grid marker schematic
Fig. 13 3D grid marker schematic
采用该折流燃烧室整机性能试验数据[2]中巡航状态设计点（相对转速95%）工况中的燃烧室进口数据作为本文折流燃烧室算例模拟的进口条件。总温为493.6K，总压为485822.5Pa的空气以1.233kg/s的质量流量经模型左上方的轴向扩压器流入折流燃烧室的1/10模型中。流场达到稳定后，三维瞬态速度矢量图如图 14 所示。Fig. 14 3D velocity vector diagram of cold flow field图 15 中展示了三维流线图向 X-Y 平面的二维投影，从图中可知空气经轴向扩压器流入燃烧室后，分为三股主要气流：第 1 股气流经过火焰筒外壳的进气斗进入火焰筒，在图中用蓝色箭头标记，这股气流的流量约占进口流量的75%，用于与在主燃区燃烧后的高温燃气掺混；第 2 股气流首先沿火焰筒外壳向离心叶轮轴心方向流动，并经前进气锥进入火焰筒主燃区，图中用黄色箭头标记了其主要流动方向，这股气流约占进口气流量的12.5%；第 3 股气流首先在燃烧室外套和火焰筒外壳之间大致沿轴向流动，然后沿径向通过空心的涡轮导向叶片流经火焰筒内壳与鼓筒轴之间的流道，最终进入主燃区，这股气流也约占进口气流量的12.5%。
Fig. 15 Schematic diagram of 3D flow lines and three streams of cold flow field在采用折流燃烧室的发动机工作时，航空煤油经高速旋转的甩油盘雾化为煤油液滴进入主燃区，煤油液滴蒸发后在主燃区与图 15 中的第 2 和第 3 股气流掺混并剧烈燃烧。图中可以发现第2 和第 3 股气流在主燃区形成了稳定的涡结构，有利于增强气相煤油与空气的掺混和保持燃烧反应稳定进行。在主燃区反应后的高温燃气与进气斗中流入的第 1 股气流掺混，进行补燃并降低了燃气温度，掺混后的燃气经涡轮导向器流出燃烧室。在图 2 中折流燃烧室模型的基础上大幅简化，仅保留火焰筒内部流道的简单结构，作为液滴粒子轨迹求解的计算域，如图 16 所示。Fig. 16 Computational domain model for solving droplet particle trajectories扫描图 16 中的简化火焰筒流道几何，获取边界处的法向量并生成用于求解粒子轨迹的网格标记。用于气相湍流燃烧模拟和粒子轨迹求解的网格标记如图 17 所示，其中灰色部分表示用欧拉法求解气相湍流燃烧时的燃烧室壁面，蓝色部分表示用拉格朗日法求解煤油液滴运动轨迹时的计算域。
Fig. 17 3D grid markers for solving gas phase combustion and droplet particle trajectories首先根据两种网格标记求解得到无燃烧反应时的液滴粒子运动轨迹，如图 18 所示。煤油液滴从甩油盘边缘沿图中黄色箭头方向进入计算域，之后受到红色箭头 1 表示的来自前进气锥气流的作用下向燃烧室后方运动至靠近火焰筒内壳的位置，其后跟随用红色箭头 2 表示的从火焰筒内壳流入的气流运动并与红色箭头 3 所示的前进气锥气流相遇，沿图 15 中所示的主燃区涡结构的边缘向火焰筒外壳运动。部分粒子运动到靠近火焰筒外壳的液滴运动的计算域边界，并在反弹后随红色箭头 4 所示的进气斗进入的气流继续运动。Fig. 18 Trajectory of kerosene droplet particles in cold working condition对于本文中模拟的折流燃烧室1/10模型，其巡航状态设计点（相对转速95%）工况的燃油流量为0.0195225kg/s。采用C12H23作为航空煤油的替代组分，并基于烷烃通用四步反应机理[20]模拟其燃烧化学反应。点火并计算一段时间后，得到充分发展的湍流燃烧场。图 19 中展示了稳定燃烧状态下瞬态的中央截面气相煤油浓度分布和速度矢量分布，气相煤油输运方程的计算域与煤油液滴运动的计算域不同，对应于未简化的折流燃烧室模型。液相煤油在主燃区内全部蒸发为气相煤油，且气相煤油在与进气斗流入的冷气流掺混前基本完全参与反应而分解。在主燃区中心位置，高速煤油液滴沿径向带动气流运动。因此速度矢量图中，在煤油浓度较高的位置，存在沿煤油液滴运动方向流动的小股气流。而在燃烧器出口位置，由于燃烧放热导致气体密度降低，喷雾燃烧状态下的气流速度明显高于图 14 中冷态流动的速度。Fig. 19 Central section transient kerosene mass fraction and velocity vector diagram
图 20 为折流燃烧室中央截面的瞬态温度云图。在温度云图中，主燃区除了由前进气锥和火焰筒内壁流入的两股冷气流产生了低温区外，燃油液滴运动轨迹上也存在蒸发吸热以及气相煤油分解引起的低温区，其分布基本与图 19 中的气相煤油浓度区域一致。Fig. 20 Central section transient temperature cloud map
真实折流燃烧室的非定常燃烧流场中，湍流掺混作用强烈。图 21 中为混合分数与温度对应关系的散点图，散点的颜色表示释热率（HRR）。由于瞬态燃烧场的部分区域中的燃烧化学反应未达到平衡，在湍流掺混充分的区域仍存在温度偏低而释热率高的点，这些点对应控制体内正在发生强烈的氧化放热反应。
Fig. 21 Scatter plot of mixing fraction, temperature and heat release rate折流燃烧室两相湍流燃烧模拟中，采用式(9)-(12)所示的燃烧反应机理求解煤油燃烧化学反应，其中𝑛 = 12，𝑚 = 23。
上式中的四步燃烧机理应用于航空煤油燃烧模拟时，参与反应的组分包括 
。其中含碳元素的包 
，三种物质反应率对应的散点图如图 22所示，碳元素的守恒性使散点形成了一个平面。图中散点的颜色表示含碳组分浓度对应的释热率，可知释热率随CO和CO2反应率增大而增大，随煤油分解的反应率增大而降低。
(b) Projection of the carbon conservation planeFig. 22 Scatter plot of heat release rate distribution on the carbon conservation plane同样地，机理中含氢元素的组分包括C12H23、H2和H2O，其反应率对应的散点图如图 23 所示，图中的散点形成了代表氢元素守恒的平面。含氢组分反应率散点图中的释热率分布与含碳组分反应率散点图不同。因为CO和CO2生成放热，生成率增大均使释热率增大，而H2生成吸热，H2O生成放热，二者生成率增大对释热率的影响相反，因此含氢与含碳组分反应率散点图上的释热率分布不同。对于含氢组分，释热率随H2O反应率增大而增大，随H2生成和煤油分解的反应率增大而降低。Fig. 23 Scatter plot of heat release rate distribution on the hydrogen conservation plane整机实验中的燃烧室出口温度测点如图 24 及图 25 所示，等距分布的镍铬-镍硅热电偶被镶嵌于空心涡轮导向器叶片前缘作为温度测点。统计模拟得到的温度场，得到折流燃烧室中的时均温度分布，燃烧室中央截面时均温度如图 26 所示。
Fig. 24 Total temperature probes on a turbine guide vane[2]
Fig. 25 Distribution of temperature measuring points[2]
Fig. 26 Time-averaged temperature cloud of central section将涡轮导向器叶片前截面上的流体计算域的时均温度作周向平均，得到时均径向温度分布曲线，并与实验中的涡轮导向器叶片上测点测量得到数据对比，得到燃烧室出口径向温度分布曲线，如图 27 所示。在涡轮导向器叶片叶尖和叶根附近温度较低，而在叶片中部靠近叶尖的位置较高，呈倒“C”形状。模拟得到的温度分布规律与实验测量的结果一致，这种温度分布规律有利于延长热端部件的工作寿命，符合燃烧室设计的要求。
Fig. 27 Time-averaged radial temperature distribution of the combustion chamber exit measurementpoint表 2 中在模拟得到的燃烧室出口时均径向温度分布曲线上对应于实验测点的位置取值，并将其与实验测量得到的燃烧室出口温度对比。模拟得到的出口温度总体偏高，在涡轮导向器叶根附近与实验值相差不大，而在叶尖附近误差偏高，最大相对误差为 69.51%，五个测点位置的平均相对误差为17.95%。模拟得到的温度分布相对实验数据整体偏高，并且相对温度较高的高温段与实验相差较大，这可能是因为实验中的测点位于空心涡轮叶片上，测得的温度受到涡轮叶片温度的影响，而涡轮叶片中央通过的由上至下的冷气流对叶片高温段冷却效果明显，如图 28 所示。相对高度较低的测点附近的涡轮叶片壁面受叶片内部冷气流影响较小，而模拟结果在这些位置的误差较低表明模拟结果及误差分析是合理的。一方面，考虑到实验测点被镶嵌于空心的涡轮导向器叶片上，测点测得的温度受到叶片本身温度的影响，因此测得的温度可能略低于流场中的温度。另一方面，模拟中采用绝热壁面边界条件，且不考虑辐射传热，不同于实际情况。并且，本文中采用C12H23作为航空煤油的替代组分，其物性和热值与实验中采用的航空煤油可能存在差异。另外，本文模拟折流燃烧室时将甩油盘初始雾化出的煤油液雾简化为从一个固定的环形喷口喷出的液雾，也可能对温度分布造成影响。四、结论
本文基于隐式浸没边界方法，结合 LES-TPDF 模型，采用北航仿真中心自研软件 AECSCIBM，针对真实几何结构的某型号折流燃烧室 1/10 模型开展两相湍流燃烧模拟研究，得到以下结论：（1）本文介绍的曲线坐标系下的隐式浸没边界方法及网格标记生成方法可将折流燃烧室几何模型精确映射为 IBM 网格标记，并基于网格标记为湍流燃烧模拟中的控制方程组提供边界条件。（2）本文中，对基于欧拉框架的湍流燃烧模拟和基于拉格朗日框架的喷雾粒子轨迹求解分别采用不同的 IBM 网格标记，能够实现复杂壁面附近的粒子反弹求解，从而实现真实几何燃烧室内的两相燃烧模拟。（3）由流场模拟结果可知，高压空气流入燃烧室后分为三股主要气流，分别约占进口流量的75%、12.5%和12.5%，用于参与高温燃气掺混和气相燃料的燃烧反应。（4）燃烧场中，反应机理中含氢（或碳）元素组分集 合的反应率在散点图中形成表征该元素守恒性质的平面，在这种平面上可直观表现出组分集 合中的某一组分对释热率等参数的影响。（5）基于自研软件 AECSC-IBM 的模型和算法能够高精度模拟折流燃烧室中的两相湍流燃烧现象，对涡轮导向器叶片前缘温度径向分布模拟的平均误差为17.95%。因此 AECSC-IBM 软件适合应用于结构复杂的真实折流燃烧室中的两相湍流燃烧高保真模拟，能够为航空发动机燃烧室精细化设计提供数据参考。[1] 宋双文, 蒋荣伟, 陈剑, 等. 离心甩油折流环形燃烧室的性能试验[J]. 航空动力学报, 2007, 22(9): 1401-1404.
[2] 刘邓欢, 金捷, 王方, 等. 涡轮喷气式发动机整机环境下折流燃烧室性能实验[J]. 航空动力学报, 2015, 30(10):2410-2415.
[3] 王方, 王煜栋, 姜胜利, 等. AECSC-JASMIN 湍流燃烧仿真软件研发和检验[J]. 航空学报, 2021, 42(12): 128-140.
[4] 周玉珍,邵文清,王菊金,等. 涡扇发动机燃烧室三维冷态流场的计算[J].推进技术, 2001, 06: 480-482.（ZHOUYu-zhen, SHAO Wen-qing, WANG Ju-jin, et al. Flow-field Computation of Turbofan Engine Combustion Chamber[J].Journal of Propulsion Technology, 2001, 06: 128-140.）
[5] 余宗明, 黄勇,王方. 方腔折流燃烧室冷态流场研究[J].推进技术, 2010, 05: 533-538.（YU Zong-ming, HUANGYong, WANG Fang. A Study on the Flowfield of a Z-flowpath Combustor[J]. Journal of Propulsion Technology,2010,05: 533-538.）
[6] 颜应文, 宋双文,胡好生, 等. 折流燃烧室两相喷雾流场的数值模拟[J].航空动力学报, 2011, 26(5): 1003-1010.
[7] 解晓东,李本威,伍恒,等.折流燃烧室燃烧特性研究[J]. 兵器装备工程学报, 2018, (2): 74-79.
[8] Sivaramakrishna G, Kumar S K. CFD analysis of gas turbine slinger combustor[C]. Bangalore: ACFD Symposium,2012.
[9] Rana R, PrajwalA, Sivaramakrishna G, et al. Conjugate Heat Transfer Analysis of a Small Annular Combustor withCentrifugal Fuel Injection System[R]. ASME GTINDIA2019-2356.
[10] 阎超, 于剑, 徐晶磊, 等. CFD 模拟方法的发展成就与展望[J]. 力学进展, 2011, 41(5): 562-589.[11] Peskin C S. Flow Patterns Around Heart Values: A Numerical Method[J]. Journal of Computational Physics, 1972,10(2): 252-271.
[12] 张漫, 王铮钧, 王晶, 等. 航空发动机内流全场流动的大涡模拟[J]. 航空动力, 2021, 19(2): 57-60.
[13] 金捷, 刘邓欢. 航空发动机燃烧室湍流两相燃烧模型发展现状[J]. 南京航空航天大学学报, 2016, 03: 303-309．
[14] Jones W P, Prasad V N. Large Eddy Simulation of the Sandia Flame Series (D-F) Using the Eulerian Stochastic FieldMethod[J]. Combustion and Flame, 2010, 157(9): 1621-1636.
[15] Jones W P, Marquis A J, Wang F. Large Eddy Simulation of a Premixed Propane Turbulent Bluff Body Flame Usingthe Eulerian Stochastic Field Method[J].Fuel, 2015, 140: 514-525.
[16] Wang F, Liu R, Dou L, et al. A Dual Timescale Model for Micromixing and Its Application in LES/TPDF Simulationsof Turbulent Nonpremixed Flames[J]. Chinese Journal of Aeronautics, 2019, 32(4): 52-64.
[17]王煜栋,王方,周佳伟,等.AECSC-IBM 航空发动机燃烧室数值模拟软件研发与检验[J/OL].航空动力学报:1-11.DOI:10.13224/j.cnki.jasp.20220216.
[18] Lagae A, Dutre P. An Efficient Ray-Quadrilateral Intersection Test[J]. Journal of Graphics Tools, 2005, 10(4): 23-32.
[19]Meier U, Heinze J, Freitag S, et al. Spray and Flame Structure of a Generic Injector at Aeroengine Conditions[J].Journal of Engineering for Gas Turbines and Power, 2012, 134(3): 031503.1- 031503.9.
[20] KunduKP, Deur J. A Simplified Reaction Mechanism for Calculation of Emissions in Hydrocarbon (Jet-A)Combustion[C]. Monterey:29th Joint Propulsion Conference and Exhibit, 1993.
声明：本文是《推进技术》网络首发 论文《基于曲线坐标系浸没边界方法的折流燃烧室模拟》，作者是王煜栋 1，王方 1,2,3，甘甜 1，金捷 1,2,3（1.北京航空航天大学 能源与动力工程学院，北京 100191；2.北京航空航天大学 江西研究院，江西 南昌 330096；3.北京航空航天大学 成都航空动力创新研究院，四川 成都 611930）
