首页/文章/ 详情

云讲堂|基于LES-TPDF方法的航空发动机折流燃烧室中的两相湍流燃烧模拟

4月前浏览4300
   

   

 导读:航空发动机折流燃烧室几何结构复杂,其高保真数值模拟需要高效的网格与边界条件处理方法。采用曲线坐标系隐式浸没边界方法结合大涡模拟-概率密度函数输运方程湍流燃烧模型开发自研软件,并实现 WP11 中折流燃烧室的高保真模拟。流动模拟中准确解析了该燃烧室中的三股主要气流,且三股气流分别约占进口流量的75%、12.5%和12.5%。两相燃烧模拟中针对拉格朗日框架下的液滴运动和欧拉框架下的湍流燃烧采用不同网格标记,模拟得到的出口径向温度分布规律与实验一致,平均相对误差为17.95%,表明基于本方法开发的自研软件能准确模拟折流燃烧室中的两相湍流燃烧现象。 
SIMPOP
一、写在文前  
离心甩油折流燃烧室广泛应用于中小型航空发动机中,利用中小型发动机转子的高转速通过离心甩油盘实现燃油雾化,具有供油压力低、雾化质量稳定、燃烧效率高和高空性能好等优点[1]。美国 TeledyneCAE 公司、Williams 公司及法国 Turbomeca 公司均生产了多种型号的采用折流燃烧室的航空发动机[2],国内涡轴系列发动机也大量采用了折流燃烧室。
航空发动机燃烧室几何结构复杂,其中的湍流燃烧现象存在强烈非线性关系,迫切需要高保真、高效率、高精度的数值模拟技术[3]。受限于真实折流燃烧室复杂的几何结构,目前对折流燃烧室的数值模拟研究大多基于采用非结构网格或在将大幅简化燃烧室几何后采用贴体结构网格。
而随着近年来计算机硬件的发展,应用 IBM 和LES-TPDF 模型模拟真实折流燃烧室成为可能。AECSC-IBM(Aero Engine Combustor Simulation Code based on Immersed Boundary Method)是北京航空航天大学(北航)能源与动力工程学院航空发动机数值仿真中心(仿真中心)在 LES-TPDF湍流燃烧模型算法基础[14-16]上结合曲线坐标系隐式 IBM[17]研发的适用于复杂结构燃烧室内高保真模拟的两相湍流燃烧数值模拟软件,其中包括几何模型导入、网格标记生成、边界条件模化以及两相湍流燃烧的求解等功能的实现均采用自研程序和算法。
本文首先介绍该软件中网格标记的生成、IBM 结合 LES-TPDF 模型的实现以及针对折流燃烧室中甩油盘甩出的高速液滴粒子轨迹求解方法等算法,并基于该软件模拟某型折流燃烧室地面实验工况下的两相湍流燃烧。2 研究对象与计算方法2.1 研究对象基于 AECSC-IBM 软件模拟和研究某型真实结构的折流燃烧室,其结构与内部流动情况如图1 所示,经压气机压缩过的高压空气经轴向扩压器进入燃烧室,并分为三股流入火焰筒,参加燃烧反应和与燃气掺混,液相煤油由高速旋转的甩油盘进入主燃区。本文中所模拟该全环形燃烧室的1/10模型如图 2 所示。
二、研究对象与计算方法  

1、研究对象

基于 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

2、IBM 模型与 LES-TPDF 的结合

为将 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、欧拉-拉格朗日两相模型与 IBM 的结合

航空发动机燃烧室中,航空煤油经初始雾化形成大量煤油液滴,煤油液滴蒸发为气相煤油后与空气掺混,并在一定条件下发生燃烧反应。本文采用欧拉-拉格朗日两相模型模拟液滴粒子与气流的质量和动量相互作用。两相模拟中,液滴的破碎、蒸发以及与流体的作用力求解均为亚网格的,与燃烧室壁面的位置无直接关系。但在求解液滴粒子的运动轨迹时,则需要根据壁面位置和法向量等信息计算粒子与壁面碰撞后的反弹或黏附情况。
当采用如图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

三、计算结果及分析  

1、数值模型验证

基于以上算法形成的 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

2、折流燃烧室 

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

3、折流燃烧室模拟与分析

采用该折流燃烧室整机性能试验数据[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 plane
Fig. 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 软件适合应用于结构复杂的真实折流燃烧室中的两相湍流燃烧高保真模拟,能够为航空发动机燃烧室精细化设计提供数据参考。
 
7月11日20时,2024仿真产学沿用线上报告会将邀请北京航空航天大学副教授王方老师做基于LES-TPDF方法的航空发动机燃烧室数值模拟》讲座,感兴趣的读者可提前报名,收藏和分享感兴趣的朋友。
参考文献  

[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.


声明:本文是《推进技术》网络首发 论文《基于曲线坐标系浸没边界方法的折流燃烧室模拟》,作者:王煜栋,王方,甘甜,金捷



来源:多相流在线
ACTSystem碰撞非线性燃烧化学湍流通用航空航天兵器UGUM
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-07-13
最近编辑:4月前
积鼎科技
联系我们13162025768
获赞 108粉丝 109文章 299课程 0
点赞
收藏
作者推荐

积鼎CFD界面追踪方法Level Set与VOF在气泡流动模拟的效果比较

“如果说单相流CFD的终极问题是湍流,那么多相流的终极问题简直就没有完”。对于两相流模拟,模型主要分为两大类:高相分数模型和界面捕捉类模型。当我们关注水中的含气量(气泡界面及气泡形状可忽略),则采用高相分数模型,此模型适用于气泡特别多的流动问题。对于有明确边界的流体-流体问题,基本上就是看如何捕捉边界面。常用的界面捕捉模型包括LS(Level Set)方法和VOF(Volume of Fluid)方法。多相流功能的CFD软件,首先就是针对这种有边界面的问题。目前主流商业CFD软件大多采用VOF方法,多相流分析软件Virtualflow采用Level Set 方法进行界面流仿真。1Level Set 方法Level Set方法是基于空间曲面的隐函数表达。关于Level Set 方法,往期文章已做过介绍:👉🏼 基于Level Set方法的液体燃料雾化DNS模拟在LS方法中,每一个时间步都要重新初始化LS方程,在时刻 tn 求得的LS函数与控制方程一起求解得到下一时刻的LS函数,这些初始化的过程中总伴随着界面位置的移动,会造成质量损失,导致质量不守恒。而改善初始化步骤来矫正质量守恒又会增加计算时间,提升计算成本。同时,因为LS方法采用的是光滑的距离函数来捕捉相界面,各个物理量可以在界面上光滑连续地过渡,且相界面的捕捉效果好。2VOF方法在VOF(Volume of Fluid)方法中,用来划分两相界面的函数是体积分数α,表示的是单个网格内的液体体积与这个网格总体积的比值。若求出整个网格相分数,如图1(a)来构造界面,会发现体积分数α在空间上是一个阶梯函数,在空间上是不连续的,从而重构出来的相界面,如图1(b)是间断的,两个相邻网格的界面是不连续的,且物理量在通过界面时也是不连续的,这个现象称为寄生流动(parasitic current),目前VOF方法的主要工作就是缓解数值方法造成的寄生流动现象。 图1 相分数空间分布(a)及其界面重构(b)与LS方法类似,根据相分数可以得到界面上的单位法向量和曲率以及计算域中的密度和粘度。为了解决物理量在界面两端不连续的问题,引入连续表面力模型(Continuum Surface Force),通过以下公式将界面内的压力表示为压力的连续函数。式中c是界面处的位置函数(图2),其表示为: 图2 CSF模型下c位置函数(左)及连续压力函数(右)由以上公式可以推出曲面微元上的表面张力。此外,在OpenFOAM中,为了求解动量守恒方程中的压力项和体积力项,定义prgh,如下:式中为网格位置矢量,对该公式求其梯度得到:该公式可以直接带入动量守恒方程中进行计算。在OpenFOAM中,使用VOF方法后在控制方程中添加了一个求解α的相方程:为了界面的尖锐,OpenFOAM采用Waller提出的方法,在相方程中添加人工对流项,从而保证界面的清晰。其中:c为压缩因子,值为0时表示不存在人工压缩,给c赋值后有利于提高界面清晰度,但同时也会提高计算成本和产生收敛问题。界面捕捉方法优点缺点VOF质量守恒阶梯函数不连续,依赖网格的细化程度,需要进行界面重构,对复杂尖锐界面的模拟效果不理想LS用光滑的LS函数捕捉界面,界面法向的计算精度高,复杂界面处理能力强需重新初始化LS函数,导致质量不守恒,且计算成本增加 3界面捕捉效果:LS vs VOF本文主要讨论多相流分析软件Virtualflow中用到的LS方法和开源软件OpenFOAM中用到的VOF方法。 为了验证LS方法和VOF方法对界面捕捉的效果,下面展示Albadawi文献中采用这两种方法计算模拟的气泡变化过程,并与实验进行了对比分析。 在计算域底部中心一个小孔以恒定体积流率喷射气泡,由于压力、浮力和表面张力的共同作用,气泡会经历产生,变形和分离的过程。计算域及其物性参数如下: 图3 计算域及其边界条件表2 两相流物性参数取tDet 为气泡分离时间,t/tDet =0 为初始时间,此时气泡已经是轴对称的球形,实验中各个时刻气泡形状如图4所示。 图4 不同时刻气泡形状 图5 LS方法(Virtualflow软件)模拟结果与实验数据对比 图6 VOF方法(OpenFOAM)模拟结果与实验数据对比由图5与图6可以看出:与实验结果对比,LS方法对界面捕捉的效果更好,VOF方法只能模拟出气泡变形的大致趋势,各个时刻气泡的高度都要比实验值的低。此外,还可以根据气泡的重心位置和纵横比来比较LS方法和VOF方法的模拟效果。 图7 气泡重心位置随时间的变化 图8 气泡纵横比随时间的变化由图7可以看出,在气泡产生到发生一段变形时,LS方法和VOF方法都可以很好地预测气泡的重心位置,但随着气泡的继续演化,VOF模拟得到的气泡重心会不断的偏离实验值,而LS方法模拟得到的重心轨迹与实验吻合较好。气泡的纵横比是指气泡的最大高度与最大宽度的比值,由图7可以看出,在VOF方法中模拟得到的纵横比会沿着实验值震荡,表现为模拟得到的气泡会在轴向上产生周期性的膨胀和收缩,而这一现象是实验观察中不存在的,而LS方法可以很有效地捕捉气泡演化时的形状以及气泡分离的时间。 更多Virtualflow模拟案例以下为更多Virtualflow软件Level Set 界面跟追方法模拟案例: 图9 液体碰撞壁面反弹过程(LS方法) 图10 膜态沸腾(3D)来源:多相流在线

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈