王治宇,陈朋,林伟,仝毅恒,郭康康,黄卫东,聂万胜
(航天工程大学宇航科学与技术系,北京 101416)
高频不稳定燃烧问题在液体火箭发动机、固体火箭发动机、燃气轮机等热力装置中广泛存在[1-2]。其诱发原因是非定常热释放与燃烧室内声波的相互耦合。高频不稳定燃烧会引发高幅值周期性压力振荡和热负荷急剧增加,对燃烧室造成巨大破坏。不稳定燃烧问题的发生涉及湍流燃烧和传热、化学反应、火焰和声波的相互作用等[3-4],而这些过程又有着不同的时间尺度。因此,长久以来并没有普适的数学物理模型将其彻底描述清楚。在液体火箭发动机中,工程上采用隔板、声腔等被动控制装置对不稳定燃烧进行抑制。但这通常以牺牲燃烧效率、增加燃烧室重量和结构复杂性为代价,甚至需要进行多次全尺寸试车才能确定被动控制装置的有效性。主动控制技术则是通过传感器和控制机构等,对推进剂流量等变量进行控制,避免不稳定燃烧的发生。不稳定燃烧的“早期辨识”和“提前预测”可为未来在液体火箭发动机中发展主动控制提供技术支撑。近年来,基于数据驱动的不稳定燃烧早期辨识与预测越来越具备可行性[5-7]。
分析燃烧过程的混沌动力学特性变化是研究不稳定燃烧早期辨识和预测的可行途径。高频不稳定燃烧的发展过程通常表现为从燃烧噪声到热声耦合极限环振荡的演变[8],该过程伴随着非线性混沌动力学特性的变化。文献[9]表明,燃烧噪声阶段往往体现了较强的混沌特性,而不稳定燃烧的发展阶段则经历了混沌特性的减弱。燃烧噪声阶段的压力扰动是非周期性的,在经过Hopf 分岔后,将演变为等振幅的极限环震荡[10]。Murayama等[11]和Gotoda等[12]指出,燃烧噪声阶段的混沌动力学特性能预测燃烧系统后续的重要变化。为此,Kasthuri等[13]基于普渡大学的矩形发动机试验数据,利用递归定量分析和多重分形理论等,研究了多种非线性特征在热声耦合早期辨识中的可行性。Lyu等[14]采用统计学指标和人工神经网络对环形燃烧器的热声耦合预测问题进行了研究。然而,上述研究主要采用单次试验数据进行验证,较难考察这些方法的泛化能力。采用的指标较为单一,这在一定程度上限制了模型非线性特征的表达能力和通用性。单一的指标难以刻画燃烧室内动力学特征的非线性,导致无法细致区分不同燃烧状态。这些指标往往与燃烧室几何构型、推进剂组合等因素密切相关,开展早期辨识或预测时采用的阈值也随之变化,严重限制了模型的通用性。比如,当改变燃烧室几何构型后,原先适用的阈值和指标可能不再适用。
高频不稳定燃烧预测的两个挑战是早期辨识的有效性和预测模型的鲁棒性。早期辨识的有效性是指不稳定燃烧尚处于发展阶段时的预测。液体火箭发动机可在极短的时间内(100 ms)[1]从稳定运行状态发展出高频不稳定燃烧。在如此短的时间内捕捉热声振荡的相关特征并进行预测是极大的挑战。预测模型的鲁棒性指深度学习模型的泛化能力。比如,燃烧室几何构型调整后预测模型的可应用性。
为应对这两个挑战,引入机器学习模型建立预测模型。这既能增强模型对燃烧信号中非线性特征的刻画能力,又便于建立更为通用的不稳定燃烧预测框架。Sarkar等[15]针对高速摄影的图片序列,基于深度学习开展了燃烧不稳定性的实时诊断。Kobayashi等[16]研究了利用压力信号和释热速率等信息来开展燃烧不稳定性预测的可行性。Cellier等[5]基于旋流燃烧器的试验数据,开展了基于深度学习的燃烧不稳定性预测,取得了较好的效果。Lyu等[6]利用长短时神经网络对环形燃烧室的声压进行了预测,实现了不稳定燃烧的早期辨识,证明了该方法提前预测热声耦合的有效性。
上述研究表明,燃烧室的动态压力信号、释热信号等,能揭示燃烧噪声的非线性动力学特征,具备开展燃烧不稳定性预测的潜力。注意到,相比于原始信号,递归矩阵[13]中包含了丰富的混沌动力学特征,且便于与深度学习模型相结合进行特征的提取与计算。因此,本文提出采用递归矩阵作为深度学习模型的前处理步骤,利用深度学习模型对这些矩阵特征进行训练,基于大量的热试车试验数据开展有监督学习,从而对高频不稳定燃烧进行提前预测。为验证早期辨识的有效性和预测框架的鲁棒性,采用包含不同运行工况和燃烧室几何构型的热试车试验数据进行训练和预测[17],既能考察早期辨识的有效性,又能对预测模型的泛化能力进行验证。
单喷嘴火箭模型发动机采用的推进剂为氧气/甲烷,实物图如图1 所示[17]。发动机结构的细节尺寸如图2 所示,其中氧喷嘴总长度70 mm,靠近氧腔的上游部分内径为9.6 mm,靠近燃烧室的下游部分内径为8 mm,甲烷环缝宽度为1 mm,喷嘴缩进长度为2 mm。该发动机采用同轴剪切式喷嘴,氧气经由氧腔沿轴向通过氧喷嘴进入燃烧室,甲烷通过8 个径向孔进入喷嘴环缝,氧气和甲烷在喷嘴缩进室发生最初的相互作用,之后进入燃烧室。燃烧室为圆柱形燃烧室,直径35 mm,可变长度在180 mm到452 mm 之间。喷管内型面为带有收缩和扩张段的拉瓦尔喷管,喷管喉部直径21 mm,扩张比为2.7。发动机采用火花塞点火方案,火花塞安装在燃烧室侧壁,距离喷注面60 mm,单次点火能量2 J,点火频率50 Hz。发动机热试车试验过程中的动态压力测量主要通过安装在燃烧室侧壁面的压电式动态压力传感器实现。该传感器型号为Kistler 6043A,量程为±10 MPa,数据采样频率为100 kHz。两支动态压力传感器的安装位置分别为距离喷注面37 mm处和距离燃烧室尾部端面37 mm 处。热试车试验中,获得的高频燃烧不稳定性类型包括一阶和二阶纵向不稳定。采用高频压力传感器能完整捕捉到燃烧室内动力学系统的变化。因此,高频动态压力信号可作为预测框架的输入信号。
图1 气氧甲烷单喷嘴火箭模型发动机[17]Fig.1 O2/CH4 single-injector model rocket combustor[17]
图2 发动机结构示意图[17]Fig.2 Structure schematic of the rocket combustor[17]
图3 为高频不稳定燃烧的试验数据结果,其中图3(a)为靠近燃烧室头部的高频压力传感器记录的压力-时间轨迹,在等幅振荡阶段,压力振荡峰-峰值约为0.25 MPa,占比平均室压25%。图3(b)为图3(a)中压力信号的时频谱图。时频谱图中,短时傅里叶变换的每个计算窗口大小为4 096 个数据样本,重叠数据样本为2 048个。由于采用的是火花塞点火,图中所出现的压力峰为点火压力峰。大约经过50 ms后,燃烧室发展为典型的热声不稳定燃烧。
图3 前三组热试车试验的动态压力及对应时频谱图Fig.3 Dynamic pressure signals and corresponding spectrogram of the first three groups of hot tests
标准运行工况条件下,氧气和甲烷处于室温(约为295 K),流量分别为0.17 kg/s 和0.08 kg/s。通过热力计算程序CEA(Chemical equilibrium with applications)计算可得燃烧室压力为1.225 MPa,燃烧室温度为2 692 K。氧气和甲烷流量通过推进剂供应管路中的音速喷嘴控制。试验时控制系统根据设定时序,控制电磁阀和火花塞的开关来控制发动机点火和熄火,并采用富氧点火和富燃关机方案,燃烧室采用热沉冷却,因此每次试验点火时间较短。25 组试验中均发生了明显的热声耦合不稳定。如表1所示,不同热试车试验之间,几何构型和运行工况略有变化,因此高频不稳定燃烧的主频和振型也会相应发生变化。试验过程中一共出现了两种高频燃烧不稳定性,分别对应一阶纵向和二阶纵向模态。二阶纵向模态都出现在燃烧室长度为316.5 mm的试验中,其余均为一阶纵向模态。
表1 各组试验操作设置情况Table 1 Overview of used test runs
模型的训练和预测都基于滑动窗口的形式,如图4 所示。若当前时间窗动态压力的振荡峰-峰值低于平均室压的10%,则认为该时间窗样本是稳定燃烧样本,否则为不稳定燃烧样本。样本的标签为:若下个滑动时间窗样本为不稳定燃烧,则当前时间窗标记为阳性(positive),否则为阴性(negative)。据此对剩下25 组试验数据进行滑动窗口划分和标签的标定。采用交叉验证的方式来评估模型泛化能力。以第一组试验数据作为测试组为例,使用剩下的24 组试验数据进行训练。输入信号为各个组动态压力时序的时间窗样本,输出为阴性或阳性的二进制标签,即二分类的输出。为评估模型的准确率,采用若干评价指标,包括准确率、受试者工作特性(Receiver operating characteristic,ROC)曲线等[7]。
图4 通过时间切片的滑动窗口对不同试验样本进行操作Fig.4 Manipulating different experimental samples through sliding windows of time slicing
总体框架如图5 所示,包括两大部分。第一部分是基于混沌动力学分析的前处理,通过混沌分析将压力时序转化为递归矩阵(Recurrence matrix),从而提取压力信号的非线性特征。第二部分是基于深度学习的预测计算模块[18]。骨架网络为ResNet-18模型和CBAM(Convolutional block attention module)注意力模块[19]。其中,Resnet-18 模型的主干网络由18层网络组成,包括多种卷积层(Convolution)、池化层(Maxpooling)、归一化层等。ResNet-18 模型是深度学习领域代表性的模型,可用于提取递归矩阵的特征。CBAM 注意力模块为嵌入式模块,用于提升网络模型的准确率。预测框架的输入信号为传感器记录的高频动态压力。输出变量为:0 表示稳定燃烧(阴性),而1 表示不稳定燃烧(阳性)。为便于分析,将该框架记为RRC。
图5 RRC模型的燃烧不稳定预测总体框架Fig.5 The technical framework of the RRC model for combustion instability prediction
相空间重构计算递归矩阵的第一步,通过将压力振荡信号组成的时间序列拓展到高维空间,在一定程度上恢复吸引子结构,从而获得系统的非线性特征。根据Takens 定理[20],对于N个离散数据对应的一维压力时间序列{p'(i)|i=1,2,…,N},引入延迟时间τ和适当的嵌入维度m可构建相空间矢量:
式中:t=1,2,…,N-(m-1)τ。
根据Takens 的理论,时间延迟和嵌入维数是互不相关的,因此可以对其独立求解。本文采用自相关函数法对时间延迟进行计算[21],采用Cao 方法计算嵌入维数[22]。确定了最佳时间延迟和嵌入维数后,可计算相应的递归矩阵[20]。
首先进行相空间重构,得到维度为m,时间延迟为τ的时间序列。按照式(2)计算重构时间序列中第i向量、第j向量的距离:
进一步根据式(3)计算递归矩阵:
式中:ε是阈值常数,H1为Heaviside函数。
基于滑动窗口划分,得到所有试验数据集对应的递归矩阵即可完成深度学习模型的前处理。随机抽取80%的数据集作为训练集,开展深度学习模型训练。其余为验证集,以监测训练过程。训练环境为搭载Nvidia 3080 显卡的Cuda 11.2,并在该环境下基于Pytorch搭建深度学习框架。
标签集由0 和1 组成,分别表示阴性(稳定燃烧)和阳性(不稳定燃烧)。将测试集输入到训练好的模型后,输出的是0 到1 之间的概率,即表示在下一个时间窗口内发生不稳定燃烧的概率。输出与标签集之间的误差则用于评估模型预测的准确性。当模型输出超过0.5,表明下个时间窗内会发生燃烧不稳定,反之则不发生燃烧不稳定。准确率指标a定义为:
式中:NS表示需要测试的数据总量;ni表示其中第i条数据的准确率,具体计算方式如下:
式中:H2为Heaviside函数为RRC的决策值,y为真实的标签值。
以第一组热试车数据为例时,在交叉验证过程中,使用剩下的24组数据对预测模型进行训练。以训练后的模型对第一组试验数据进行测试,与真实标签值对比后得到准确率。进行验证测试的试验数据对预测模型而言是未知的,可以较好地验证模型的泛化能力。特别是,25 组试验中包含了不同燃烧室几何构型和运行工况的数据,这更加考验模型的泛化能力。完成交叉验证后,得到25组测试的准确率,如图6所示。可以看到,所有测试评估组的准确率都稳定在95%以上,且其中大部分组的准确率超过了98%。因此,交叉验证的准确率结果较好地验证了预测模型的有效性和泛化能力。
图6 交叉测试组的准确率结果Fig.6 Accuracy results of the cross testing group
注意到试验数据集阴性和阳性的样本分布是不平衡的。在不稳定燃烧的预测中,最重要的是尽可能辨识到那些可能发生燃烧不稳定的样本,从而采取相关控制措施提前避免。因此,为进一步刻画预测模型在不同标签数据下的表现,引入混淆矩阵对阴性和阳性样本进行量化分析。该矩阵可描述预测类别和实际类别之间的差异。对于二分类的模型,结果可分为以下4 种。真阳性样本(True positive,sTP):模型预测为发生不稳定燃烧的样本,实际也为不稳定燃烧;假阳性样本(False positive,sFP):模型预测为发生不稳定燃烧,但实际情况是稳定燃烧的;真阴性样本(True negative,sTN):模型预测为稳定燃烧,实际情况也是如此;假阴性样本(False negative,sFN):模型预测为稳定燃烧,但实际情况中却发生了不稳定燃烧。
当模型预测样本为阳性,即发生不稳定燃烧时,只有相对较少的一部分属于误判。这意味着预测模型应用在实际发动机中时,做出阳性决策的可信度比较高。可以用精确率φ来表示预测结果为阳性的样本中正确预测的比例:φ=sTP(sTP+sFP)。
精确率越高,模型的阳性预测的可信度也越高。前4 组数据测试后的精确率分别为:92.00%,100%,92.51%,99.65%,这表明了预测模型辨识不稳定燃烧的可信度是非常高的。另一个指标是假阳性率φfpr,可以形象地理解为误诊率,即错误地将原本属于稳定燃烧的样本误判为发生不稳定燃烧,计算方式为:φfpr=sFP(sFP+sTN)。
前4 组数据试验测试后的假阳性率分别为14%,0.00%,16%,0.78%。可见第二组和第四组的假阳性率非常低,而第一组和第三组相对更高。这意味着在第一组和第三组中,有一部分稳定的样本被模型误判为会发生不稳定燃烧。从不稳定燃烧的实际意义来分析,更为关键的是不能丢失对不稳定燃烧样本的预测,所以假阳性率偏高依然是可以接受的。此外,引入召回率η来表征模型对所有阳性样本的识别成功率。计算方式为:η=sTP(sTP+sFN)。
召回率的现实意义在于预测模型有多大的概率预测出全部的不稳定燃烧样本。当没有成功预测出后续发生的不稳定燃烧时,发动机将面临不稳定燃烧的危害。因此预测模型应追求更高的召回率,从而避免在预测中出现不稳定燃烧样本的丢失。经计算,前4 组试验数据测试后的召回率分别为100.00%,96.66%,100.00%,96.56%。由此可见,预测模型对不稳定燃烧(阳性样本)的预测成功率非常高,其中只有很少一部分阳性样本被误判为了稳定燃烧样本。在不稳定燃烧的预测中,为捕捉到更多的不稳定燃烧样本,即追求更高的召回率,可以牺牲一定的假阳性率为代价,这可以通过降低一定的预测阈值来实现。此时模型将更加倾向于将样本预测为不稳定燃烧,从而采取主动控制措施以避免不稳定燃烧。在实际发动机中,假阳性率偏高的代价是模型会错误地将稳定样本也预测为不稳定燃烧。召回率偏低则表明许多实际发生了不稳定燃烧的样本并没有被成功预测,导致系统并未对不稳定燃烧采取措施。因此相对而言,召回率偏低的代价会造成更显著的危害。由此可见,以假阳性率为代价换更高的召回率,通常是可行的。
召回率又被称之为真阳性率φtpr,与假阳性率φfpr组合后可进行ROC曲线的计算[5],能很好地刻画非平衡样本中模型的分类性能。假阳性率为ROC曲线的x轴,真阳性率为y轴。通过调整模型的阈值可获得不同的φtpr/φfpr数据坐标点,也就是模型对应的ROC曲线[5]。
ROC 曲线可以动态评估预测模型将以多大的误判代价来实现不稳定燃烧样本的成功预测。当假阳性率很低时,性能优秀的模型仍能保持较高的真阳性率。因此性能越好的模型,对应的ROC 曲线应更逼近ROC图中的左上方。以前4组试验数据为例,计算ROC曲线,如图7所示。
李庚先生自幼喜爱中国古典文学艺术,对绘画更加痴迷,秉承其父李可染训教,钻研水墨,取得了丰硕成果,在国际上声誉日隆。清华大学教授、李可染画院副院长王鲁湘赞其“血液中流淌着父亲李可染先生对神秘崇高的宇宙精神悠然神往的可贵品质。他对自然神性有一种发自灵府的会意……山水风云的无言大美,水墨语言的玄化鸿蒙,在他看来都是道通天地的律动,他所做的就是‘从而和之’,让气韵发自笔端。李庚的水墨画把父亲早年的墨戏和晚年的抽象推到了一个风神飘举的境界,其绘画作品玄韵淡泊,情思淹济,风气韵度似父亲而更恣纵迈达”。
图7 前4组热试车试验数据预测结果的ROC曲线Fig.7 ROC curves of the first four hot-fire tests
进一步分析图7,注意到所有的ROC 曲线都经过坐标点(0,0)和坐标点(1,1)。当阈值设定为最高时,所有预测结果都为阴性,伪阳性率φfpr和真阳性率φfpr都为0,此时对应坐标点(0,0)。当阈值设定为最低时,所有预测结果都为阳性,此时伪阳性率和真阳性率都为1,对应了ROC 曲线的坐标点(1,1)。图中的虚线表示随机模型的ROC曲线。
图7中,4 组试验数据的ROC 曲线体现了非常理想的预测结果。为进一步考察预测模型在不稳定燃烧发展过程中早期辨识的能力,以这4 组试验数据为例,将不稳定燃烧发展过程中模型的决策函数和实时压力数据进行对比,如图8所示。
图8 前4组的压力信号与对应的RRC模型预测结果对比Fig.8 Dynamic pressure signals of the first four verified groups and corresponding prediction results obtained by RRC model
热试车试验中,高频不稳定燃烧发展速度极快,只需要100 ms 左右的时间。对早期辨识的准确性和有效性提出了较高的要求。图8中的预测点表示RRC 模型的实时预测值,而黑色实线为实时动态压力信号。为方便与压力信号进行对比,将输出值阈值从0.5 调整为0。低于阈值的输出表示下一个时间窗口内不会出现不稳定燃烧,反之则表示会发生不稳定燃烧。
在0.5 s出现的点火压力峰之前,模型的预测稳定维持在阈值下方。发动机未点火启动之前,传感器所测的系统噪声与后续是否发生不稳定燃烧是无关的。但预测模型也需要能够对这一过程的物理信号进行辨识和预测,以保证预测运行的连续性。燃烧室内点火后形成了极高的压力峰,压力的剧烈振荡改变了燃烧室的动力学特性。预测模型的输出值也随之产生变化,在阈值下方产生了逐渐增加的趋势。但由于此时模型仍然无法判断系统后续是否会发生不稳定燃烧,一段时间内仍处于阈值以下。随着时间窗口继续滑动,预测系统采集到更多的动力学特征,此时的预测值超过了阈值,表明下个时间窗内发生不稳定燃烧的概率超过了0.5。
在阈值附近区域振荡一段时间后,最终稳定在阈值上方的水平。此时表明了预测系统判断下一个时间窗内发生不稳定燃烧的概率非常大。而在之后的不稳定燃烧阶段,模型输出值一直稳定保持在阈值上方的水平。这是因为发生了不稳定燃烧后,动力学特征容易区分和预测。由此可见,燃烧状态监测的实现是简单且稳定的,较难预测的是由稳定燃烧噪声阶段到热声耦合阶段发展的过程。
在不稳定发展阶段(从点火后到热声振荡阶段),当预测的输出值稳定在最大值附近时,认为模型已完成热声耦合动力学特征的早期辨识,实现了不稳定燃烧的提前预测。将这一时间节点称为稳定预测点,如图8 所示。稳定预测点到最终试验中发展出不稳定燃烧的时间,是模型预测该组试验不稳定燃烧的提前时间。4 组的预测提前时间分别为:31.52 ms、34.58 ms、47.16 ms、20.15 ms。这4 组热试车的平均提前预测时间为33.35 ms。经计算25组预测结果后,平均提前时间为35.28 ms。
综上,基于前4 组试验数据对预测框架进行了有效性和泛化能力的评估,考察了预测准确率、精确率和ROC 等指标,并分析了各组的提前预测时间,可为主动控制系统提供一定时间裕度。
本文提出了一种较为通用的高频不稳定燃烧早期辨识和预测分析框架。该框架融合了递归矩阵和深度学习模型。递归矩阵作为ResNet-18 深度学习模型的前处理,可提供关于燃烧噪声的非线性动力学特征。ResNet-18 模型可基于递归矩阵进行大量试验数据的有监督学习,从而对未知的热试车试验进行燃烧不稳定性的早期辨识和预测。基于25 组热试车试验的交叉验证,评估了预测框架的有效性和鲁棒性。结果表明,本文提出的高频不稳定燃烧预测框架可实现热声耦合的早期辨识和预测,平均提前时间为35.28 ms。