1. 清华大学 能源与动力工程系，北京 100084;
2. 山东科技大学 土木工程与建筑学院，青岛 266590
收稿日期：2023-02-20
基金项目：国家科技重大专项(Y2019-I-0022-0021)
作者简介：耿俊杰(1998—)，男，博士研究生
通讯作者：祁海鹰，教授，E-mail: hyqi@tsinghua.edu.cn
摘要：为了构建适用于燃烧室冷态流动混合过程的代理模型(surrogate model, SM)的方法，该文研究了构建过程中的关键步骤，通过使用Latin超立方抽样(Latin hypercube sampling, LHS)方法进行样本选取，在完成数值模拟后使用本征正交分解提取样本间的主要特征进行降维，再通过Kriging插值法完成工况插值。结果表明：该构建方法能够处理冷态高速、高湍流度、强旋流动和燃料/空气掺混，精度高于国际平均水平。同时，该文提出了构建方法的应用准则，为后续更复杂、包含燃烧反应过程以及结构变化的燃烧室SM构建奠定基础。
关键词：燃烧室代理模型Latin超立方抽样本征正交分解Kriging插值
Surrogate model of combustor flow mixing process
GENG Junjie1, WANG Xingjian1, LI Jialu2, FEI Teng1, QI Haiying1
1. Department of Energy and Power Engineering, Tsinghua University, Beijing 100084, China;
2. College of Civil Engineering and Architecture, Shandong University of Science and Technology, Qingdao 266590, China
Abstract: Objective Numerical simulation of a gas turbine combustor is an important step in its design process. Due to the complexity of the physical and chemical processes, the calculation cost is high. The calculation cost can be reduced by constructing surrogate model of combustor. This paper focuses on the key steps in the construction of a surrogate model suitable for cold gas flow and the mixing process of combustor. Furthermore, this paper proposes a surrogate model for the central nozzle and the subsequent combustor space of a heavy-duty gas turbine. Methods The construction of the surrogate model includes several key steps: design of experiments (DOE), numerical simulation, dimensionality reduction, and an interpolation process. Two parameters are selected as the input parameters for the surrogate model: The fuel mass flow rate Gf and the combustor inlet air pressure p2. Latin hypercube sampling is used in the DOE to determine 12 operating conditions for computational fluid dynamics (CFD) simulations, and the results are used to build the surrogate model. Proper orthogonal decomposition is used for dimensionality reduction, wherein a set of basis functions and corresponding coefficients are extracted. The basis functions reflects the main characteristics of the combustor flow field. Moreover, the data dimensionality is reduced from the number of grid nodes of the combustor to the number of basis functions, which do not exceed the number of operating conditions. The Kriging model is used to interpolate the coefficients of the basis function with the input parameters of the surrogate model. Four verification conditions are set up to determine the accuracy of the surrogate model through a comparison of the surrogate model results with the CFD simulation results. The outlet cross section of the central nozzle and the longitudinal section of the combustor are selected to compare multiple key parameter distributions, including axial velocity, radial velocity, tangential velocity, CH4 concentration, turbulent kinetic energy, and pressure. The vector operations are used to compare the distributions of various parameters, which can simultaneously reflect the differences in the numerical and spatial distributions of various parameters. Results The results showed that the error in most parameters was ~1%. The results also revealed that the construction method of the surrogate model could be applied to cold high-speed and high-turbulence strong rotational flow and fuel/air mixing. The accuracy was higher than the international average level, and the application criteria of the construction method were proposed. The influences of interpolation methods, sample numbers, and basis function numbers on the accuracy of the surrogate model were analyzed. The accuracy of SM was higher than extrapolation. Increasing the number of sample operations and basis functions could improve the accuracy of the surrogate model, but also increased the computational cost. Conclusions The SM construction method (POD & Kriging) is suitable for the cold gas flow and mixing process in the combustor. The paper lays the foundation for subsequent research on the construction method of combustor SM, which includes combustion reactions and geometric structure changes.
Key words: surrogate model of combustorLatin hypercube samplingproper orthogonal decompositionKriging interpolation
对燃气轮机燃烧室进行数值模拟是设计过程中的重要环节。然而，全尺寸燃烧室的数值模拟具有空间尺度跨度大，以及流动、传质传热、燃烧和污染物生成等多个复杂过程的耦合特性，需要使用多种数学模型，并且具有网格复杂、收敛难度大、耗时长的特点[1]。例如，对于一个约600万网格的算例，在960 CPU(单核主频2.5 GHz)的条件下，需要运行3 h才能获得收敛解。在燃烧室设计过程中，需要对设计和结构参数等进行调整和反复迭代，因此需要进行大规模数值验算，会导致大量的CFD算例。目前，国际上发展出一种代理模型(surrogate model, SM)方法，采用在有限数量的计算流体动力学(computational fluid dynamics, CFD)算例中的插值算法，可迅速、完整地得出多个中间工况的流场参数，取代原本需要的三维CFD算例[2-4]。对于同样的算例，SM在数秒内即可完成，与CFD相比，显著提高了工作效率。
SM是一种基于一定量样本数据来描述输入变量和输出数据之间关系且计算代价较小的数学模型[5-6]。SM的构建方法包含3个基本步骤：样本选取、数据降维和工况插值[2]。
样本选取即形成一个算例集，如同插值所需的一套离散数据点集。样本的选取质量对SM的精度有较大影响，因此需要在变量空间中均匀地分布样本点，以提高精度。目前有以下几种选取方法，例如正交法、均匀设计法、中心复合法，以及最常用的Latin超立方抽样方法(Latin hypercube sampling, LHS)等。其中，LHS是一种在多维变量空间中的分层抽样方法，随机性强，可以根据需要任意选择样本点的数量[7-8]。
在完成样本算例的CFD计算后，如果对每一个网格点上的流场参数进行插值，计算量会非常大，因此需要对流场数据进行降维。本征正交分解(proper orthogonal decomposition, POD)是一种常用的方法，通过POD提取出一组代表流场主要特征的基函数矩阵(基函数个数不超过样本数)，使用这组基函数以级数叠加的形式来计算流场数据，从而将计算维数降低到与基函数的个数一致，以大幅提高计算效率[3-4, 9]。
插值方法包括多项式法、径向基函数法、神经网络法和Kriging插值法等。其中Kriging插值法是最常用的最优线性无偏估计方法之一，被广泛应用于地质学、环境科学、空间遥感和计算统计等领域[10-13]。
目前，SM已广泛应用于多个领域。例如，国内外分别构建了飞行器气动外形参数与气动特性之间的SM，包括升阻比、升力系数等。SM也应用于火电厂空冷凝汽器内部流场、温度场与环境风速的关系建模，以及发动机旋流喷嘴冷态流场的建模等[2, 7, 14]。然而，文献研究中只涉及单一流体流动，并且未见有关燃烧室SM的研究，而燃烧室燃烧是一种多过程耦合的复杂系统。例如，三维高湍流度的旋流流动、燃料射流、燃料/空气混合过程、传热、燃烧化学反应和污染物生成等。因此，从设计角度还需对SM的结构进行改变。同时，多参数和多过程的耦合可能会对数据降维、插值精度等产生影响，进而影响SM的精度。
本文先构建燃烧室冷态物理混合过程的SM，检验SM方法的可行性，为后续构建包含燃烧过程、更复杂条件下的模型奠定基础。为了避免燃烧室复杂过程带来的干扰，拟先做冷态且仅关注状态参数的变化，即不涉及燃烧反应过程，不考虑结构参数的变化。
为构建适用于燃烧室冷态流动混合过程的SM的方法，本文对构建过程中的样本选取、数据降维和工况插值等关键步骤进行了研究，并制定了应用准则，为后续研究更复杂的燃烧反应过程和结构变化的模型构建奠定基础。
1 燃烧室数值模拟方法1.1 基本结构和工作模式燃烧室的基本结构如图 1所示，图中数值的单位均为mm。由导流衬套、头部喷嘴组件、火焰筒和过渡段组成。燃烧室采用径向分级的结构，分为中心喷嘴和环形区；环形区内包含5个沿周向分布的自带旋流器的喷嘴，且在环形区外壁有2排主燃孔，在火焰筒下游周向均布3个掺混孔。
 |
| 图 1 燃烧室结构示意图 |
| 图选项 |
中心喷嘴的工作模式为：空气从喷嘴外环进入预混通道，燃料则直接喷入，与空气相互掺混，形成预混可燃混合物，并通过轴向旋流器进入燃烧区，轴线处设有扩散燃烧的值班喷嘴。环形区有2种工作模式：当燃烧室相对负荷Ne≤0.8时，发生燃烧；当Ne>0.8时，仅进行预混。此前，该燃烧室已有详细的数值模拟和试验验证，SM在制作样本算例时，使用原来的网格[15]。
1.2 计算域选取本文选取从中心喷嘴入口到过渡段入口截面的三维燃烧空间(见图 1)作为计算域，主要考虑以下几方面：1) 中心喷嘴是燃烧室的典型核心部件，在气体动力学方面有双重作用，一是通过旋流器在其下游的燃烧空间形成用于稳定火焰的回流区；二是实现燃料喷射及其与空气的预混过程，可反映与高速流动、高湍流度和燃料/空气的预混均匀度等方面有关的问题。因此，SM能展示燃烧室的特殊性和典型性。2) 在过渡段入口截面之前，回流区已结束，流场较均匀，因此可作为计算域的出口。3) 暂不考虑环形区和掺混孔，简化不必要的影响因素。4) 可同时考察和评价SM再现横截面和纵剖面的CFD结果的能力。
中心喷嘴出口截面将燃烧室分为2个区域，前面为流动混合区，后面为燃烧反应区。如果需要详细模拟燃烧的过程，可通过SM提供燃烧区的入口边界条件。同时，纵剖面可以反映心回流区、角回流区等主要流场特征。因此，选取中心喷嘴出口截面和纵剖面作为SM的考察截面。
1.3 数值模拟设置1) 计算网格。
本文采用原先工作的网格文件[15]。在整体燃烧室网格基础上去掉中心喷嘴前的导流衬套、环形区和掺混孔的网格，仅保留计算域。网格类型为非结构四面体网格，网格数为182.3万，如图 2所示。
 |
| 图 2 网格示意图 |
| 图选项 |
2) 数学模型和边界条件。
除了通用的质量、动量守恒方程和组分扩散方程，湍流模型使用Realizable k-ε模型，边界条件为：空气入口使用速度进口边界条件；截取该燃烧室中心喷嘴入口截面的速度分布，并设置进口压力；燃料入口使用质量流量入口边界条件；燃烧室出口使用自由出流边界条件；壁面为无滑移边界。
1.4 SM输入参数的选择选取在冷态条件下参数变化较大的燃料量Gf和进口压力p2这2个参数。燃料量代表燃烧室的负荷变化，直接体现中心喷嘴中的预混状况和出口处的燃料浓度分布。参数Gf的具体数值范围为0.040 0~0.120 0 kg/s，p2为1 330~1 510 kPa，对应的Ne范围为0.5~1.0。
1.5 流场和混合过程的数值模拟结果图 3展示了经CFD计算得到的燃烧室纵剖面流线图。图中的计算工况为：空气流量中Ga= 1.766 kg/s，温度T2=665.5 K，压力p2=1 349 kPa；燃料流量中Gf=0.074 0 kg/s，温度Tf=288.0 K。在横截面突扩处形成角回流区；由于旋流作用在喷嘴下游形成中心回流区，且回流区前端进入中心喷嘴内，在喷嘴的出口截面处中心区域的轴向速度为负值，外围区域的轴向速度为正值，如图 4所示。由于该喷嘴上侧某旋流叶片有探视孔，气流通过探视孔以直流的形式进入燃烧室，造成旋流强度下降，导致上侧的回流区长度较短，因此中心回流区上下并不对称。
 |
| 图 3 CFD燃烧室纵剖面流线图 |
| 图选项 |
 |
| 图 4 CFD燃烧室轴向速度云图 |
| 图选项 |
中心喷嘴的燃料浓度分布，如图 5a所示。XCH4为摩尔分数；CCH4为甲烷浓度，kmol/m3。燃料通过中心通道进入燃烧室，其中一部分作为值班燃料直接进入燃烧室；另一部分通过4个燃料喷嘴喷出，并在预混通道内与空气混合。由于预混通道的方式是直流，且预混距离较短，喷嘴出口截面处的燃料/空气预混均匀度并不理想。图 5b和5c展示了CH4的浓度分布，图中的外围区分别形成了4个间隔分布的高浓度区，对应燃料喷管的位置。
 |
| 图 5 CFD中心喷嘴内和出口截面上的燃料分布 |
| 图选项 |
2 结果与分析2.1 SM构建流程构建流程如图 6所示，使用LHS方法确定样本算例；在完成CFD计算后，通过POD进行数据降维和Kriging模型进行工况插值，完成初步构建；最后以验证算例的CFD计算结果为基准，对精度进行检验；如果误差较大，则需重新构建。
 |
| 图 6 构建流程 |
| 图选项 |
2.2 样本选取首先，确定样本个数Nsamp，主要取决于输入变量的个数Nvar。为保证SM的精度，应满足的准则可表示如下：
| ${N_{{\rm{samp }}}} \ge 5{N_{{\rm{var }}}}.$ | (1) |
样本选取的具体方法为：把每个变量的变化范围等分为Nsamp个区间，在每个区间内随机生成一个数值，这样每个变量都得到Nsamp个抽样值，然后将各变量的抽样值随机组合，即可得到Nsamp个样本。本文生成12个样本，同时也构成了12个CFD算例，将Gf按大小进行排序，样本的参数条件如表 1所示。
表 1 样本的参数条件
| 样本编号 | 1 | 2 | 3 | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 11 | 12 |
| Gf/(kg·s-1) | 0.043 0 | 0.050 0 | 0.057 0 | 0.063 0 | 0.070 0 | 0.077 0 | 0.083 0 | 0.090 0 | 0.097 0 | 0.103 0 | 0.110 0 | 0.117 0 |
| p2/kPa | 1 442.5 | 1 412.5 | 1 427.5 | 1 472.5 | 1 337.5 | 1 382.5 | 1 352.5 | 1 502.5 | 1 457.5 | 1 487.5 | 1 367.5 | 1 397.5 |
表选项
在表 1的样本中选择4个CFD验证算例。算例取自该燃烧室负荷Ne范围为0.5~1.0的真实运行参数，用于检验SM的误差，如表 2所示。表 1和2中的工况只有参数Gf、p2发生了改变，其他参数保持不变，与1.5节中的工况相同。
表 2 验证算例对应的负荷和参数条件
| 算例 | YZ1 | YZ2 | YZ3 | YZ4 |
| Ne | 0.5 | 0.8(切换前) | 0.8(切换后) | 1.0 |
| Gf /(kg·s-1) | 0.074 0 | 0.108 0 | 0.058 0 | 0.054 0 |
| p2 /kPa | 1 349.0 | 1 440.0 | 1 440.0 | 1 486.0 |
表选项
2.3 数据降维处理使用POD进行降维，首先将各个样本的流场数据，包括：速度、压力、燃料浓度、湍动能等，按照计算网格节点的顺序进行排列，构成m×1的矩阵Y t，表示如下：
| ${\mathit{\boldsymbol{Y}}_t} = \left[ {\begin{array}{*{20}{c}}{{Y_{t, 1}}}\\{{Y_{t, 2}}}\\ \vdots \\{{Y_{t, m}}}\end{array}} \right].$ | (2) |
将各样本的流场数据进行排列构成矩阵Y，表示如下：
| $\mathit{\boldsymbol{Y}} = \left[ {{\mathit{\boldsymbol{Y}}_1}, {\mathit{\boldsymbol{Y}}_2}, \cdots , {\mathit{\boldsymbol{Y}}_{{N_{{\rm{samp }}}}}}} \right].$ | (3) |
按照Galerkin加权残值法完成基函数的求解，具体算法为：计算样本矩阵Y的內积；计算矩阵Q的特征值λi及对应的特征向量ai(i=1, 2, …, Nsamp)，并按照特征值大小进行排列；计算各特征值λi所对应的函数φi，并进行归一化；计算基函数?i的系数βt, i。计算过程分别表示如下：
| $\begin{array}{c}\begin{array}{*{20}{c}}{\mathit{\boldsymbol{Q}} = {\mathit{\boldsymbol{Y}}^{\rm{T}}}\mathit{\boldsymbol{Y}}, }\\{{\lambda _1} > {\lambda _2} > \cdots > {\lambda _{{N_{{\rm{samp}}}}}}, }\\{{\mathit{\boldsymbol{\varphi }}_i} = \mathit{\boldsymbol{Y}}{\mathit{\boldsymbol{a}}_i}, }\\{{\mathit{\boldsymbol{\phi}} _i} = \frac{{{\mathit{\boldsymbol{\varphi }}_i}}}{{\left| {{\mathit{\boldsymbol{\varphi }}_i}} \right|}}, }\end{array}\\\left[ {\begin{array}{*{20}{l}}{{\beta _{t, 1}}}\\{{\beta _{t, 2}}}\\ \vdots \\{{\beta _{t, {N_{{\rm{samp}}}}}}}\end{array}} \right] = \left[ {\begin{array}{*{20}{l}}{\mathit{\boldsymbol{\phi}} _1^{\rm{T}}}\\{\mathit{\boldsymbol{\phi}} _2^{\rm{T}}}\\ \vdots \\{{\mathit{\boldsymbol{\phi}} _{N_{{\rm{samp}}}}^{\rm{T}}}}\end{array}} \right]{\mathit{\boldsymbol{Y}}_t}.\end{array}$ | (4) |
| ${\mathit{\boldsymbol{Y}}_t} = \sum\limits_{i = 1}^{{N_{{\rm{samp }}}}} {{\beta _{t, i}}} {\mathit{\boldsymbol{\phi}} _i}.$ | (5) |
1) 基函数是一组标准正交基，表示如下：
| $\left( {{\mathit{\boldsymbol{\phi}} _i}, {\mathit{\boldsymbol{\phi}} _j}} \right) = \left\{ {\begin{array}{*{20}{l}}{1, }&{i = j;}\\{0, }&{i \ne j.}\end{array}} \right.$ | (6) |
| $\max \sum\limits_{t = 1}^{{N_{{\rm{samp }}}}} {\sum\limits_{i = 1}^{{N_{{\rm{samp }}}}} {\left( {{\mathit{\boldsymbol{Y}}_t}, {\mathit{\boldsymbol{\phi}} _i}} \right)} } .$ | (7) |
2.4 Kriging插值在对流场数据进行POD后，需对系数βt, i和SM的输入参数(即Gf和p2)进行插值，从而得到非样本点上的系数值，并计算得到各参数的分布，表示如下：
| ${\mathit{\boldsymbol{Y}}_{{\rm{SM}}}} = \sum\limits_{i = 1}^{{N_{{\rm{samp }}}}} {{\beta _{{\rm{SM}}, i}}} {\mathit{\boldsymbol{\phi}} _i}.$ | (8) |
Kriging插值法是指，假设对一个目标函数y进行插值，共有n个数据点，定义x1、x2、…、xn为自变量，y1、y2、…、yn分别为对应的函数值。
该方法将y(x)表述成函数平均值μ与一个随机过程ε(x)的组合，表示如下：
| $y(x) = \mu + \varepsilon (x).$ | (9) |
| $E[\varepsilon (x)] = 0, $ | (10) |
| ${\mathop{\rm corr}\nolimits} \left[ {\varepsilon \left( {{x_l}} \right), \varepsilon \left( {{x_m}} \right)} \right] = \exp \left[ { - \sum\limits_{k = 1}^d {{\theta _k}} {{\left( {{x_{kl}} - {x_{km}}} \right)}^2}} \right].$ | (11) |
假设使用Kriging插值法求得在预测点xnew处的函数值为ynew，最优解是指该点处预测值的误差方差最小，即var(ynew-y0)最小，其中y0为真实值，使用最大似然估计方法进行推导，表示如下：
| ${y_{{\rm{new}}}} = \hat \mu + \mathit{\boldsymbol{r}}_{{\rm{new}}}^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}\left( {\mathit{\boldsymbol{y}} - {\mathit{\boldsymbol{A}}_n}\hat \mu } \right), $ | (12) |
| $\hat \mu = \frac{{\mathit{\boldsymbol{A}}_n^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}\mathit{\boldsymbol{y}}}}{{\mathit{\boldsymbol{A}}_n^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{A}}_n}}}, $ | (13) |
| $\mathit{\boldsymbol{y}} = {\left[ {{y_1}, {y_2}, \cdots , {y_n}} \right]^{\rm{T}}}{\rm{. }}$ | (14) |
| $\mathit{\boldsymbol{r}}(l, 1) = {\mathop{\rm corr}\nolimits} \left[ {\varepsilon \left( {{x_{{\rm{new}}}}} \right), \varepsilon \left( {{x_l}} \right)} \right].$ | (15) |
| ${s^2} = {\hat \sigma ^2}\left[ {1 - \mathit{\boldsymbol{r}}_{{\rm{new}}}^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{r}}_{{\rm{new}}}} + \frac{{{{\left( {1 - \mathit{\boldsymbol{A}}_n^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{A}}_n}} \right)}^2}}}{{\mathit{\boldsymbol{A}}_n^{\rm{T}}{\mathit{\boldsymbol{R}}^{ - 1}}{\mathit{\boldsymbol{A}}_n}}}} \right], $ | (16) |
| ${\hat \sigma ^2} = \frac{{{{\left( {\mathit{\boldsymbol{y}} - {\mathit{\boldsymbol{A}}_n}\hat \mu } \right)}^{\rm{T}}}{\mathit{\boldsymbol{R}}^{ - 1}}\left( {\mathit{\boldsymbol{y}} - {\mathit{\boldsymbol{A}}_n}\hat \mu } \right)}}{n}.$ | (17) |
2.5 SM计算结果采用上述方法构建SM，并计算出各验证算例(见表 2)在喷嘴出口横截面和纵剖面上各参数的分布。验证算例1(YZ1)的计算结果如图 7—11所示，va为轴向速度，vr为切向速度。由图可知，SM完全再现了CFD计算得到的横截面和纵剖面的轴向速度和切向速度分布。此外，其他所有参数，如k、p、CH4浓度场等的分布均与CFD计算结果相似。
 |
| 图 7 SM得到的轴向速度云图 |
| 图选项 |
 |
| 图 8 SM得到的径向速度 |
| 图选项 |
 |
| 图 9 SM得到的湍动能 |
| 图选项 |
 |
| 图 10 SM得到的压力 |
| 图选项 |
 |
| 图 11 SM得到的CH4浓度场 |
| 图选项 |
2.6 误差分析为了综合评价经SM计算的变量大小和空间分布与CFD结果之间的误差δ，需要采用向量运算进行表征[2]，表示如下：
| $\delta = \frac{{\left| {{\mathit{\boldsymbol{Y}}_{{\rm{SM}}}} - {\mathit{\boldsymbol{Y}}_{{\rm{CFD}}}}} \right|}}{{\left| {{\mathit{\boldsymbol{Y}}_{{\rm{CFD}}}}} \right|}} \times 100\% .$ | (18) |
SM计算结果的误差，如表 3所示，表中vt为径向速度，可知SM在横截面和纵剖面内的计算结果的误差均≤4%，远低于国际上10%的平均水平[2-4]，表明尽管燃烧室内的气流速度和湍流度较高，存在旋流和回流区，且存在不同气体组分的混合，在这种情况下，本文构建的SM仍有较强的处理能力和较高的准确性。
表 3 SM计算结果的误差
| 算例 | 喷嘴出口截面误差/% | 纵剖面误差/% | |||||||
| YZ1 | YZ2 | YZ3 | YZ4 | YZ1 | YZ2 | YZ3 | YZ4 | ||
| va/(m·s-1) | 0.18 | 0.08 | 0.21 | 0.15 | 0.47 | 0.61 | 0.85 | 0.81 | |
| vr/(m·s-1) | 0.10 | 0.14 | 0.18 | 0.20 | 1.30 | 1.70 | 1.26 | 1.67 | |
| vt/(m·s-1) | 0.16 | 0.13 | 0.18 | 0.14 | 1.39 | 1.23 | 1.07 | 1.00 | |
| CCH4/(kmol·m-3) | 0.50 | 0.15 | 0.96 | 1.23 | 0.41 | 0.40 | 0.57 | 0.47 | |
| k/(m2·s-2) | 1.22 | 1.01 | 3.15 | 3.96 | 2.97 | 3.04 | 3.94 | 3.18 | |
| p/MPa | 3.01 | 0.11 | 3.45 | 2.89 | 0.19 | 0.29 | 0.58 | 0.27 | |
表选项
2.7 SM构建方法的应用准则1) 采用内插值的方式。
在SM中，外插和内插是指预测点工况与样本工况区间的关系。如果预测点工况在样本工况区间之外，称为外插；如果预测点工况在样本工况区间之内，称为内插。
研究表明SM是一种数据驱动(data driven，DD)方法，由于样本集不包含参数范围外的信息，采用外插会显著降低SM的精度。本文的算例也证明了这一点，外插使SM算出的CH4浓度场的误差成倍增加，如图 12所示。图中以表 1中算例11的CH4浓度场作为验证算例，从表 1和2的算例中选取8个不同的组合构建SM，并且这些组合的下边界值为Gf, down=0.043 3 kg/s，上边界值Gf, up由0.117 0 kg/s变化至0.074 0 kg/s，可表示如下：
| $\Delta {G_{\rm{f}}} = {G_{{\rm{f}}, 11}} - {G_{{\rm{f}}, {\rm{up}}}}.$ | (19) |
 |
| 图 12 外插方式导致的误差 |
| 图选项 |
其中ΔGf为样本算例11的边界值与样本组合上边界值的差值。ΔGf为负值时，表明插值方式为内插，此时以表 1中的算例12为上边界，其余所有组合均使算例11构成了外插；ΔGf越大，表明验证算例外插得越远。
2) 样本数目足够多。
从数学原理上讲，用于插值的样本数量越多，SM计算结果的精度就越高。使用验证算例1的喷嘴出口截面的CH4场进行验证，当样本数量Nsamp由12降至7时，误差约增大2倍，如图 13所示。
 |
| 图 13 样本个数对SM误差的影响 |
| 图选项 |
3) 基函数个数选取合理。
在POD分解中，得到一组基函数?i，分别对应样本內积矩阵的特征值λi(按从大到小排列)。通常前几个特征值较大，其对应的基函数也代表了样本分布的主要特征。因此，在实际应用中，通常不会使用所有的基函数，而只选用前几阶基函数就可以达到较高的精度，以提高计算效率。定义前q个基函数的能量系数ξ[2]，表示如下：
| $\xi = \frac{{\sum\limits_{j = 1}^{{N_{{\rm{samp }}}}} {\sum\limits_{i = 1}^q {{{\left( {{\beta _{j, i}}{\mathit{\boldsymbol{\phi}} _i}} \right)}^2}} } }}{{\sum\limits_{j = 1}^{{N_{{\rm{samp }}}}} {\sum\limits_{i = 1}^{{N_{{\rm{samp }}}}} {{{\left( {{\beta _{j, i}}{\mathit{\boldsymbol{\phi}} _i}} \right)}^2}} } }} \times 100\% .$ | (20) |
表 4 SM误差、能量系数与基函数个数
| 基函数个数 | 1 | 2 | 3 | 4 | 12 | |
| 能量系数/% | 99.98 | 99.98 | 99.99 | 99.99 | 100.00 | |
| δ/% | YZ1 | 0.863 1 | 0.542 3 | 0.498 8 | 0.498 2 | 0.498 0 |
| YZ2 | 0.633 9 | 0.157 5 | 0.156 0 | 0.154 0 | 0.153 0 | |
| YZ3 | 2.035 3 | 1.029 4 | 1.012 7 | 1.010 8 | 0.963 1 | |
| YZ4 | 2.456 5 | 1.262 0 | 1.239 1 | 1.235 6 | 1.235 6 | |
表选项
3 结论本文提出了一种构建燃烧室冷态流动混合过程的SM的方法，使用了LHS方法、本征正交分解、Kriging插值等算法，并证明了其合理可行。在只改变状态参数的情况下，SM具有处理冷态高速、高湍流、强旋流动和燃料/空气掺混的能力，其精度高于国际平均水平。本文还制定了SM构建方法的应用准则，包括：1) 使用内插值的方式；2) 遵循样本个数和基函数个数选取的约束准则。为后续更复杂、包含燃烧反应过程以及结构变化的燃烧室SM构建方法研究奠定基础。
参考文献
| [1] | 焦树建. 燃气轮机燃烧室[M]. 北京: 机械工业出版社, 1981. JIAO S J. Gas turbine combustor[M]. Beijing: China Machine Press, 1981. (in Chinese) |
| [2] | MAK S, SUNG C L, WANG X J, et al. An efficient surrogate model for emulation and physics extraction of large eddy simulations[J]. Journal of the American Statistical Association, 2018, 113(524): 1443-1456. DOI:10.1080/01621459.2017.1409123 |
| [3] | WANG X J, YEH S T, CHANG Y H, et al. A high-fidelity design methodology using LES-based simulation and POD-based emulation: A case study of swirl injectors[J]. Chinese Journal of Aeronautics, 2018, 31(9): 1855-1869. DOI:10.1016/j.cja.2018.07.004 |
| [4] | YEH S T, WANG X J, SUNG C L, et al. Common proper orthogonal decomposition-based spatiotemporal emulator for design exploration[J]. AIAA Journal, 2018, 56(6): 2429-2442. DOI:10.2514/1.J056640 |
| [5] | 朱新奇. 飞行器外形优化设计方法及多输出代理模型研究[D]. 西安: 西北工业大学, 2019. ZHU X Q. Research on flight vehicle shape optimization design algorithms and multi-output surrogate model[D]. Xi'an: Northwestern Polytechnical University, 2019. (in Chinese) |
| [6] | 吕利叶. 先进代理模型方法与应用研究[D]. 大连: 大连理工大学, 2020. LV L Y. On advanced surrogate models and application[D]. Dalian: Dalian University of Technology, 2020. (in Chinese) |
| [7] | 夏陈超. 基于CFD的飞行器高保真度气动外形优化设计方法[D]. 杭州: 浙江大学, 2016. XIA C C. High-fidelity aerodynamic shape optimization method of aircraft based on computational fluid dynamics[D]. Hangzhou: Zhejiang University, 2016. (in Chinese) |
| [8] | 李俊芳, 张步涵. 基于进化算法改进拉丁超立方抽样的概率潮流计算[J]. 中国电机工程学报, 2011, 31(25): 90-96. LI J F, ZHANG B H. Probabilistic load flow based on improved Latin hypercube sampling with evolutionary algorithm[J]. Proceedings of the CSEE, 2011, 31(25): 90-96. (in Chinese) |
| [9] | WANG X J, CHANG Y H, LI Y X, et al. Surrogate-based modeling for emulation of supercritical injector flow and combustion[J]. Proceedings of the Combustion Institute, 2021, 38(4): 6393-6401. DOI:10.1016/j.proci.2020.06.303 |
| [10] | 冯俞楷. POD降维算法在传热与流动数值模拟中的应用[D]. 北京: 华北电力大学, 2017. FENG Y K. Application of POD reduced-order algorithm on heat transfer and flow simulation[D]. Beijing: North China Electric Power University, 2017. (in Chinese) |
| [11] | 黄煌. 风资源评估后处理系统的实现与并行优化[D]. 北京: 清华大学, 2015. HUANG H. The implementantion and optimization of post-processing system for wind resource assessment[D]. Beijing: Tsinghua University, 2015. (in Chinese) |
| [12] | 本杰明. 基于遥感和机器学习从地表特征对地下水深度的预测方法研究[D]. 北京: 清华大学, 2019. BEN J M. Predicting groundwater depth from surface features using remote sensing and machine learning[D]. Beijing: Tsinghua University, 2019. (in Chinese) |
| [13] | 邹林君. 基于Kriging模型的全局优化方法研究[D]. 武汉: 华中科技大学, 2011. ZOU L J. Research on global optimization algorithm based on kriging model[D]. Wuhan: Huazhong University of Science and Technology, 2011. (in Chinese) |
| [14] | 胡和敏. 火电空冷系统跨尺度热质传递的数值模拟研究[D]. 北京: 华北电力大学, 2014. HU H M. Cross scale simulation on heat and mass transfer of air-cooled condenser in power plant[D]. Beijing: North China Electric Power University, 2014. (in Chinese) |
| [15] | 祁海鹰, 田园, 耿俊杰. X重型燃机联合研制低排放燃烧室模拟计算项目报告[R]. 北京: 清华大学, 2021. QI H Y, TIAN Y, GENG J J. Report on simulation of low emission combustor for a heavy gas turbine[R]. Beijing: Tsinghua University, 2021. (in Chinese) |
| [16] | SACKS J, WELCH W J, MITCHELL T J, et al. Design and analysis of computer experiments[J]. Statistical Science, 1989, 4(4): 409-423. |
