吸气式高超声速飞行器推进系统建模及不确定性研究

 电话：[phone] 手机：[phone]

微信：tech99999 邮箱：[e-mail]

来源：战术导弹技术

作者：赵黄达，杨超，吴志刚

摘要：针对吸气式高超声速飞行器推进系统存在的不确定性问题，采用随机多项式展开(PCE) 和蒙特卡罗模拟(MCS)两种方法开展了推进系统不确定性定量研究，得到了推力的概率密度函数以及置信区间，以对推力的不确定性进行评估。首先对一典型吸气式高超声速飞行器的推进系统进行了建模，将模型分为外压缩段及内流道段，考虑溢流效应的影响，得到了气流参数分布及推力结果。随后考虑马赫数、迎角及燃油当量比的不确定性，采用PCE方法进行不确定分析，将得到的推力估计值、置信区间及计算效率与MCS结果进行对比。结果表明：所建立的推进系统模型可快速评估吸气式高超声速飞行器的推力特性，PCE和MCS两种方法的推力统计分布结果相吻合，但PCE方法计算效率更高，可以在初步设计阶段快速评估推进系统的不确定性。

关键词：高超声速飞行器; 不确定性；超燃冲压发动机；随机多项式展开方法；蒙特卡罗模拟

      吸气式高超声速飞行器作为一种重要的战略平台，得到了越来越多的关注，但该类飞行器存在着较强的机体/推进系统耦合问题，推进系统对全机的力学特性产生了重要影响，因此需要有效的建模仿真方法。然而，推进系统的复杂流动特性，使其建模和分析过程中存在很大不确定性，因此在分析中需要考虑推进系统的不确定性影响。

      国内外学者对推进系统的建模仿真工作进行了大量研究。Chavez等建立了一个低精度超燃冲压发动机模型。该模型包括扩散段、燃烧室以及内喷管三部分，建模中假设扩散段和喷管为等熵流动，通过总温的增加模拟燃烧过程。在一维模型基础上，MASIV( Michigan-AFRL Scramjet In Vehicle) 、SRGULL等二维分析程序考虑了隔离段预燃烧激波串、燃料混合情况、壁面摩擦等因素的影响，对超燃冲压发动机进行了更高精度的建模仿真。三维发动机建模分析方法精度最高，可以通过CFD 方法以及燃烧模拟软件CHEMKIN等得到流场情况，但通常需要耗费巨大的计算时间。近年来，使用影响系数法建立超燃冲压发动机模型的工作得到了一定的发展，该方法原理简单、计算效率高，得到了越来越多的关注。Shapiro对影响系数法进行了详细介绍，分析了摩擦、截面变化、热传导对一维流动的影响。文献[8]使用影响系数法分析了超燃冲压发动机内流场，并与三维CFD数值计算进行了对比，结果表明该方法准确性较好，可对发动机进行快速仿真计算。

 推进系统涉及复杂的流体力学、燃烧动力学问题，首先需要对推进系统的不确定性进行建模分析。传统的蒙特卡罗模拟方法计算直接，但计算量巨大。对于随机不确定参数的分析，随机多项式展开方法比较适用，且计算效率高。文献[9]使用非浸入式随机多项式展开方法进行了概率颤振分析，得到了舵面在高风险速度范围内的失稳概率。文献[10]基于随机响应面方法进行了结构可靠度分析，对随机变量相关情况下的方法有效性进行了验证。

 针对上述推进系统建模和不确定性研究工作可以发现，现有的推进系统建模仿真方法大多面向推进系统专业需求，可以实现对复杂燃烧流动的模拟，但对于气动弹性力学、飞行力学等学科，分析时往往需要多次改变推进系统的工况条件，现有的高精度的模型过于复杂，计算效率太低，不能满足上述学科的需求。此外，推进系统本身的不确定性分析也需要大量计算工作。这就要求推进系统模型在处理上述问题时，要在反应流动基本特性的同时，兼具一定的计算效率。另一方面，在现有的不确定性相关研究工作中，PCE方法得到了一定的应用，但是尚未应用于吸气式高超声速飞行器推力不确定性的研究工作中，在处理该类问题时的精度和效率需要进一步的讨论。

 针对以上问题，本文基于影响系数法对吸气式高超声速飞行器推进系统进行了建模仿真，在此基础上，考虑来流马赫数、迎角、以及燃油当量比的不确定性，通过PCE方法进行了不确定性分析，并与MCS结果对比，得到推力的概率密度函数、推力估计值以及置信区间。

1. 推进系统建模 

考虑典型的吸气式高超声速飞行器构型，将飞行器下表面作为推进系统的研究对象进行建模，如图1所示。区域Ⅰ为外压缩段，区域Ⅱ为内压缩段，区域Ⅲ包含发动机隔离段、燃烧室以及后体/尾喷管。

图1 推进系统示意图推进系统分析流程如图2所示。整个模型的输入参数为来流条件及燃油当量比，输出为推力。外压缩段（区域Ⅰ）的气流参数可以通过斜激波关系式给出，内流道段（区域Ⅱ、Ⅲ）通过影响系数法进行建模。下面分别对两部分的建模过程进行介绍。

图2 推进系统分析流程图（12CM）

1.1 外压缩段

吸气式高超声速飞行器前体下表面通常为预压缩面，产生一道或多道激波，从而改变气流的流动参数及方向。一般在设计状态下，前体产生的斜激波将汇聚于唇口。对于本文的研究对象，在Ma=7、零迎角飞行条件时，飞行器前体处于设计状态，此时前体产生的一道激波恰好与唇口相交。当马赫数小于设计马赫数时，激波角增大，使得激波位置前移，产生溢流，影响发动机进气道捕获流量。

图3 考虑溢流影响的外压缩段模型

      溢流影响下的外压缩段模型如图3所示，此时气流转角为迎角和几何转角之和，即。前体激波将不再与唇口相交，捕获流量为对应的气流高度。通过坐标变换及几何方程求解，得到点的纵坐标为：

当飞行器横向宽度为 时，捕获流量为：

在确定气流转角后，可由式（3）计算激波角，并带入斜激波关系式（4）~式（6）求得外压缩段的气流参数。

    其中，为比热比，分别为来流和波后气流的马赫数、温度及压力。

1.2 发动机内流道段

通过影响系数法对发动机内流段（区域Ⅱ、Ⅲ）进行建模，并使用区域Ⅰ的气流参数作为入口气流条件。

影响系数法假设流动是一维定常的，可以考虑内流道截面变化、壁面摩擦等因素，通过给定的释热规律来模拟化学燃烧反应。控制方程的微分形式如式（7）~式（11）所示，具体推导流程见文献[7]。

     连续性方程：

     状态方程：

     马赫数定义式：

     动量方程：

     能量方程：

     其中，为流体密度、速度以及摩尔质量，为摩擦系数，为水力直径，为流体比热容，为燃料流动方向参数，为热传导、辐射传热引起的能量变化，为气体做功引起的能量变化，为化学反应引起的能量变化。

     计算时将控制方程进行整理，忽略比热比、摩尔质量以及质量流量的变化，将马赫数、温度及压力的导数整理成式（12）形式，其中矩阵表示各独立变量（面积、能量等）对于非独立变量（马赫数、温度等）的影响，通常将其中的元素称为影响系数。     

     采用4阶龙格库塔法，将上述微分方程组按发动机轴向位置推进求解，得到内流道段的气流参数分布情况。在得到发动机出口气流参数后，采用式（13）计算发动机推力。

2. 不确定性分析方法 

2.1 随机多项式展开方法

随机多项式展开方法的核心思想是，将随机响应展开为随机变量的一系列正交基底函数的线性组合，随后求解基底系数，完成响应函数的构建，利用该代理模型进行不确定性分析。

推力的多项式展开为：

     其中，为随机变量向量，为第项Hermit多项式基底函数，为该基底对应的待求解系数。

为了得到基底系数，采用配点法进行求解。按一定的配点选取原则选择计算点，通常配点数目要达到待定系数的2倍以上。选取配点后，以基底系数为未知量组成方程组，通过最小二乘法对该过约束方程组进行求解。

2.2 蒙特卡洛模拟

蒙特卡罗模拟是一种抽样试验方法，通过大量的样本数据，来估算变量的统计分布结果。该方法直观简便，便于计算机编程处理。

推力与马赫数、迎角以及燃油当量的关系可以表示为：

分析时假设各随机变量满足随机概率分布，对各变量进行抽样并计算其响应值。当试验样本足够大时，认为响应数据平均值即为推力的估计值。

3. 仿真与结果讨论 

     在所建立的推进系统模型基础上，考虑来流马赫数、迎角以及燃油当量比的不确定性，使用随机多项式展开方法对推力进行不确定性定量分析，并通过蒙特卡罗模拟验证方法正确性。

3.1 标称状态仿真结果

本文分析时使用的发动机模型如图1所示，主要几何数据为：飞行器全长4.448m，内流道入口高度0.1015m，喉道高度0.0425m，发动机出口高度0.55m，发动机横向宽度0.228m。

参考典型的吸气式高超声速飞行器飞行条件，选取标称状态的计算参数，如表1所示，此时前体处于溢流状态。

表1 标称状态计算参数标称状态下，前体及内流道的气流参数分布如图4所示。飞行器前体气流参数由斜激波关系式给出，因此不随位置变化。内压缩段及隔离段气流继续减速，马赫数下降，温度、静压上升。燃烧室由于燃料燃烧释热，气体总能量增加，温度压力均大幅上升。后体/喷管段气体膨胀加速喷出，产生推力，标称状态下计算得到推力为1678.67N。文献[8, 13]对某吸气式高超声速飞行器模型进行了气动分析，本文流动参数变化趋势与之相同。

图4 标称状态下气流参数分布

3.2 不确定性分析结果

3.2.1 随机多项式展开方法结果

考虑来流马赫数、迎角以及燃油当量比的不确定性，将三个随机参数表示为：

     其中，为常数影响因子，本文取0.02，为随机变量，满足标准正态分布。

按3.1节方法构建响应多项式函数。对于上述三个随机参数，采用7阶Hermit多项式时，响应过程的多项式展开共有120项。选取300个计算点作为配点，进行多项式系数的求解，从而构建出响应函数。

在使用构建的响应函数进行不确定分析前，首先讨论响应函数与随机参数的关系，三个随机参数的影响中，燃油当量比为主要影响因素，与推力成正相关关系，影响机理比较清晰，不单独列出分析。此处重点分析燃油当量比不变时，推力与迎角及马赫数的关系，如图5所示。  

图5 推力响应曲面（燃油当量比0.4）

从图5可以看出，燃油当量比不变时，推力随马赫数增大而减小，随迎角增大而增大。进一步对数值的比较可以发现：马赫数增大时，发动机出口处温度、压力均增大，根据推力表达式中，占主要部分，该部分对推力的贡献减小，的贡献增大，两者综合使得推力减小。迎角增大时，出口参数并无显著变化，推力增加主要源于捕获流量的增加。

在分析了响应函数与随机参数的关系之后，利用构建好的响应函数产生新的统计样本，进而分析推力的分布情况。将计算结果利用核密度估计方法进行统计，得到推力的概率密度函数。样本数量为5000时，得到概率密度函数如图6所示。从结果可以看出，当随机参数满足正态分布时，推力响应也近似满足正态分布。推力均值为1678.57N，95%置信区间为[1604.85N，1753.45N]。

图6 标称状态下推力分布及概率密度函数

3.2.2 蒙特卡洛模拟结果为了验证随机多项式展开分析方法的正确性，对推进系统模型进行蒙特卡罗模拟，分别计算了1000、2000、5000个样本点数据下，推力与随机参数关系，计算结果如图7所示。

图7 蒙特卡罗模拟推力分布

图7中颜色越深的地方推力数值越大，从结果可以看出，推力大小随燃油当量比的改变变化比较明显。同样采用核密度估计方法，得到推力均值为1678.23N，95%置信区间为[1601.56N，1751.28N]。概率密度函数结果与图6类似，不重复列出。
3.2.3 结果对比讨论将两种方法产生的推力数据进行整理对比，结果如表2所示。

表2 PCE与MCS方法计算结果对比

从表2结果可以看出，两种方法一致性较好，验证了PCE方法的合理性。通过对样本均值的分析可以发现，MCS的样本均值精度均随样本数量增加而提高，而PCE方法得到的样本均值呈现出一定的随机性，这是由于PCE方法构建响应函数的信息源自配点处的推力值，构建完成后，统计样本数量的增加并没有对构建的响应函数产生改变，因此结果不具有通用性规律。

通过对置信区间的对比可以发现，在样本数量足够时，样本标准差以及95%置信区间主要与单次分析时随机参数样本有关，与样本数量无关。置信区间的变化量约为均值的4%-5%。

通过对计算时间的对比可以发现，相比于MCS方法，PCE方法耗时大幅减少，计算效率高。PCE方法的耗时主要集中在对配点处样本的计算，且构建出响应函数后，产生新的统计样本时，仅需要对响应函数进行求值计算，因此计算时间大幅度减小。

3.2.4 多工况计算从4.2.1~4.2.3节结果中可以看出，PCE方法在计算效率上具有较大的优势，因此可以快速得到不同工况条件下推力计算估计值及置信区间。为了说明该特点，本节选取飞行轨迹中的14个状态点进行计算。对于每一个状态，采用式（17）构建配点及样本点，其中，根据飞行工况确定，，均取1和0.4，样本数为5000，全部工况总的计算耗时约为4200s。

飞行工况的不同反映在马赫数和高度的变化，其中马赫数直接对推力计算产生影响，作为不确定性来源；高度的变化体现在飞行密度的变化上，在计算时同样对推力结果产生影响，但不作为不确定来源。分析时假定初始时刻（t=0s）飞行高度为21.7km，飞行马赫数为5.7，结束时刻(t=130s)飞行高度23km，飞行马赫数为7.0。中间计算点的高度和马赫数按线性插值给出，其中第四个状态点，为本文前述标况计算点。

在上述飞行轨迹下，飞行马赫数、高度随着时间的增加而增大。图8给出了全部工况的推力估计值及其95%置信区间，其中圆点结果为推力的估计值，灰色区域的上下界为推力的置信区间范围。结果表明，马赫数增大的影响占主要因素，使得推力的估计值随着时间不断减小，与4.2.1得到的结论一致；推力置信区间宽度大小主要受不确定参数取值范围影响，参数取值范围大致相同时，各个工况下推力置信区间宽度相差不大。

图8  推力估计值及其95%置信区间

4. 结论  

     本文针对吸气式高超声速飞行器推进系统进行建模，在此基础上，通过数值模拟分析了推力随马赫数、迎角及燃油当量比的变化规律。使用随机多项式展开方法，得到了推力的估计值及置信区间，并与蒙特卡罗模拟结果进行了对比，得出了以下主要结论：

(1)吸气式高超声速飞行器前体的溢流会影响冲压发动机的捕获流量，在建模过程中需要考虑该因素。本文提出一种基于斜激波关系式和影响系数法的推进系统建模仿真方法，该方法能够快速有效地得到推进系统推力特性，适用于该类飞行器气动弹性力学、飞行力学问题的推力特性评估。

(2)对于本文建立的推进系统模型，飞行条件在标称状态下产生扰动时，推力随马赫数增大而减小，随迎角的增大而增大。其中，燃油当量比占主要影响因素，同推力成正相关关系。

(3)蒙特卡罗模拟方法精度随样本数量增加而增加。随机多项式展开方法精度与之相当，但与样本数量无关。两者相比，随机多项式展开方法计算效率较高。可以处理多工况情况。

一网打尽系列文章，请回复以下关键词查看：

创新发展：习近平 | 创新中国 | 创新创业 | 科技体制改革 | 科技创新政策 | 协同创新 | 科研管理 | 成果转化 | 新科技革命 | 基础研究 | 产学研 | 供给侧

热点专题：军民融合 | 民参军 | 工业4.0 | 商业航天 | 智库 | 国家重点研发计划 | 基金 | 装备采办 | 博士 | 摩尔定律 | 诺贝尔奖 | 国家实验室 | 国防工业 | 十三五 | 创新教育 | 军工百强 | 试验鉴定 | 影响因子 | 双一流 | 净评估

预见未来：预见2016 |预见2020  | 预见2025 | 预见2030  | 预见2035 | 预见2045  | 预见2050

前沿科技：颠覆性技术 | 生物 | 仿生 | 脑科学 | 精准医学 | 基因 |  基因编辑 | 虚拟现实 | 增强现实 | 纳米 | 人工智能 | 机器人 | 3D打印 | 4D打印 | 太赫兹 | 云计算 | 物联网 | 互联网+ | 大数据 | 石墨烯 | 能源 | 电池 | 量子 | 超材料 | 超级计算机 | 卫星 | 北斗 | 智能制造 | 不依赖GPS导航 | 通信 | 5G | MIT技术评论 | 航空发动机 | 可穿戴 | 氮化镓 | 隐身 | 半导体 | 脑机接口 | 传感器

先进武器：中国武器 | 无人机 | 轰炸机 | 预警机 | 运输机 | 直升机 | 战斗机 | 六代机 | 网络武器 | 激光武器 | 电磁炮 | 高超声速武器 | 反无人机 | 防空反导 | 潜航器

未来战争：未来战争 | 抵消战略 | 水下战 | 网络空间战 | 分布式杀伤 | 无人机蜂群 | 太空战 | 反卫星

领先国家：美国 | 俄罗斯 | 英国 | 德国 | 法国 | 日本 | 以色列 | 印度

前沿机构：战略能力办公室 | DARPA | 快响小组 | Gartner | 硅谷 | 谷歌 | 华为 | 阿里 | 俄先期研究基金会 | 军工百强

前沿人物：钱学森 | 马斯克 | 凯文凯利 | 任正非 | 马云 | 奥巴马 | 特朗普

专家专栏：黄志澄 | 许得君 | 施一公 | 王喜文 | 贺飞 | 李萍 | 刘锋 | 王煜全 | 易本胜 | 李德毅 | 游光荣 | 刘亚威 | 赵文银 | 廖孟豪 | 谭铁牛 | 于川信 | 邬贺铨

全文收录：2017文章全收录 | 2016文章全收录 | 2015文章全收录 | 2014文章全收录

其他主题系列陆续整理中，敬请期待……