本发明属于航空航天推进系统领域，具体涉及一种燃烧室燃烧不稳定的在线预测方法。
背景技术：
在液体火箭发动机研制过程中，在燃烧室经常会观察到随机燃烧不稳定现象，即燃烧室内随机地出现自发的大幅值压力振荡。随机不稳定的一个特点是难以复现，即使在同一台燃烧室中采用相同的工况进行测试，这给燃烧室的稳定工作带来了潜在的安全风险，因此对其进行在线预测具有重要的理论和工程意义。
宽频低幅值燃烧噪声会触发燃烧不稳定，目前认为是发动机出现随机不稳定的机理。从非线性动力系统角度看，发动机中的燃烧不稳定是通过霍普夫分叉产生的，包括亚临界和超临界分叉，亚临界存在双稳定区。燃烧室处于双稳定区时，有限幅值的扰动和宽频低幅值噪声均可以触发燃烧不稳定。
对噪声触发燃烧不稳定的在线预测方法目前未见文献报道。与之相关的内容为：文献scheffer,m.,bascompte,j.,brock,w.a.,brovkin,v.,carpenter,s.r.,dakos,v.,held,h.,vannes,e.h.,rietkerk,m.,andsugihara,g.,2009,“early-warningsignalsforcriticaltransitions,”nature,461(7260),pp.53–59.和文献dakos,v.,carpenter,s.r.,brock,w.a.,ellison,a.m.,guttal,v.,ives,a.r.,kefi,s.,livina,v.,seekell,d.a.,vannes,e.h.,andscheffer,m.,2012,“methodsfordetectingearlywarningsofcriticaltransitionsintimeseriesillustratedusingsimulatedecologicaldata,”plosone,7(7),p.e41010。上述两个文献提出了临界慢化理论，即对于多数非线性系统，当系统接近分叉点时，系统对扰动的恢复率降低，同时发展了多种参数来预测临界转变(分叉)，包括恢复率降低、自相关系数、方差和条件异方差增加等。
文献gopalakrishnan,e.a.,sharma,y.,john,t.,dutta,p.s.,andsujith,r.i.,2016,“earlywarningsignalsforcriticaltransitionsinathermoacousticsystem,”sci.rep.,6(1),p.35310中分析发现方差和条件异方差可以作为热声不稳定预警信号(缓慢改变控制参数引起的不稳定，非噪声触发不稳定)，同时，延迟为1的自相关系数在噪声环境下鲁棒性较差。
文献an.,q.,steinberg,a.m.,jella,s.,bourque,g.&furi,m.2019earlywarningsignsofimminentthermoacousticoscillationsthroughcriticalslowingdown.j.engr.gasturb.power141:054501中，对燃气轮机中的自发燃烧不稳定进行分析(不能确定是噪声触发不稳定或是线性不稳定)，发现方差和延迟为1的自相关系数可以作为不稳定的预警信号。研究发现，上述参数对噪声触发不稳定同样具有预警作用。但是，在将其发展为在线预测方法时，发现直接采用参数值进行判断不可行，因为不同工况下出现燃烧不稳定时的参数值相差较大。该文献中，提出对压力数据的自相关系数曲线利用滚动窗口方法计算kendall秩相关系数，达到1时说明系统即将出现燃烧不稳定。但是，在实际应用时，由于自相关系数曲线为非单调上升，导致kendall秩相关系数呈现振荡状态，并且达不到1，判断不稳定存在困难。
综上所述，现有方法对燃烧室的随机燃烧不稳定难以在线预测，给燃烧室带来了潜在的安全风险问题。
技术实现要素：
本发明的目的是解决现有方法对随机燃烧不稳定难以在线预测的问题，提出一种燃烧室燃烧不稳定的在线预测方法，该方法实现对噪声触发不稳定的在线预测，且计算速度快，可靠性高。
为实现以上发明目的，本发明技术方案是：
一种燃烧室燃烧不稳定的在线预测方法，包括以下步骤：
步骤一、获取燃烧室的压力脉动时间序列
通过动态压力传感器获取燃烧室或燃烧室喷前的压力脉动时间序列
步骤二、在滚动计算窗口1内，去除压力脉动时间序列
步骤三、在滚动计算窗口1内，计算压力脉动时间序列的延迟为1自相关系数ρ1，随着计算窗口的滚动，得到各个窗口延迟为1自相关系数ρ1组成的延迟为1自相关系数序列ρ；
延迟为1自相关系数ρ1按下式计算，
其中，e为数学期望计算，pt 1为pt延迟为1的时间序列，μ为pt的平均值或数学期望，
步骤四、设定滚动计算窗口2的长度为m，当延迟为1自相关系数序列ρ数量大于滚动计算窗口2的长度时，在滚动计算窗口2内，计算延迟为1自相关系数序列ρ的kendall秩相关系数；
kendall秩相关系数τρ按下式计算，
其中，ρi和ρj均为滚动计算窗口2内延迟为1自相关系数序列ρ中的点，i和j为自相关系数序列ρ中的点的时间序号；
步骤五、根据阈值判断状态；
若kendall秩相关系数τρ大于阈值τ*，说明燃烧室内即将出现燃烧不稳定。
进一步地，步骤二中，在滚动计算窗口1内，利用移动平均方法去除压力脉动时间序列
2.1)计算趋势mt，按下式计算；
其中，k为移动平均窗口长度，
2.2)去除原始压力脉动时间序列中的趋势，按下式计算，
其中，pt为去除趋势后的压力脉动时间序列。
进一步地，步骤二中，也可以通过高斯光滑法或loess光滑法去除压力脉动时间序列中的趋势。
进一步地，步骤三中，方差
其中，n为滚动计算窗口1的长度。
进一步地，步骤五中，阈值为0.5。
本发明技术方案与现有技术相比，具有如下优点：
本发明提供一种燃烧室燃烧不稳定的在线预测方法，采用本发明提出的方法能够提前预测燃烧室内即将发生的燃烧不稳定现象，进而通过适当的方法进行控制，同时，该方法计算速度较快，可靠性高。
附图说明
图1为本发明方法中压力时间序列的示意图；
图2为本发明方法中计算自相关系数时滚动计算窗口1和计算kendall秩相关系数时滚动计算窗口2的示意图。
具体实施方式
以下结合附图和具体实施例对本发明的内容作进一步详细描述。
本发明提供一种燃烧室燃烧不稳定的在线预测方法，采用本发明提出的方法能够提前预测燃烧室内即将发生的燃烧不稳定现象，进而通过适当的方法进行控制。该方法包括：实时监测燃烧室压力脉动时间序列；利用移动平均方法去除时间序列中的趋势；采用滚动窗口方法计算压力脉动时间序列的延迟为1自相关系数；计算延迟为1自相关系数序列的kendall秩相关系数；基于kendall秩相关系数对燃烧室内脉动压力脉动时间序列的发展进行在线预测。
本发明方法原理是：对于非线性系统，当系统接近分叉点时，系统对扰动的恢复率降低。滞后为1自相关系数可评估系统对扰动的恢复快慢程度，其对噪声触发不稳定具有预警作用。kendall秩相关系数可对滞后为1自相关系数上升趋势进行量化度量，其值的大小可以判断随机燃烧不稳定出现的可能性。
本发明提供的燃烧室燃烧不稳定的在线预测方法具体包括以下步骤：
步骤一、获取燃烧室的压力脉动时间序列
通过动态压力传感器监测燃烧室或燃烧室喷前的压力脉动时间序列
步骤二、在滚动计算窗口1内，去除压力脉动时间序列
该步骤中，采用移动平均方法、高斯光滑法或loess光滑法去除压力脉动时间序列中的趋势，下面以移动平均方法为例，获取去除趋势后的压力脉动时间序列pt，移动平均方法相对于其他方法，运算速度较快；
2.1)在滚动计算窗口1内，利用移动平均方法去除压力脉动时间序列中的趋势，首先计算趋势mt，按下式计算
其中，k为移动平均窗口长度；
2.2)然后去除原始压力脉动时间序列
其中，pt为去除趋势后的压力脉动时间序列；
步骤三、在滚动计算窗口1内，计算压力脉动时间序列的延迟为1自相关系数ρ1，随着计算窗口的滚动，得到各个窗口延迟为1自相关系数ρ1组成的延迟为1自相关系数序列ρ；
延迟为1自相关系数ρ1按下式计算，
其中，e为数学期望计算，pt 1为pt延迟为1的时间序列；μ为pt的平均值或数学期望；
其中，n为滚动计算窗口长度；
步骤四、设定滚动计算窗口2长度为m，当延迟为1自相关系数序列ρ数量大于滚动计算窗口2的长度时，在滚动计算窗口2内，计算延迟为1自相关系数序列ρ的kendall秩相关系数τρ；
kendall秩相关系数τρ按下式计算
其中，ρi和ρj均为滚动计算窗口2内延迟为1自相关系数序列ρ中的点，i和j为自相关系数序列ρ中的点的时间序号；
步骤五、根据阈值判断状态；
如果kendall秩相关系数τρ大于阈值τ*，则说明燃烧室内即将出现燃烧不稳定。
下面以对本发明方法做举例性的进一步详细描述。
如图1所示，通过动态压力传感器监测燃烧室或燃烧室喷前的压力脉动时间序列，确定滚动计算窗口1的长度n，本实施例中的滚动计算窗口1长度为1000t，t为燃烧室最低阶固有振荡周期。当压力脉动时间序列数量大于滚动计算窗口1的长度时，开始后面的计算；
在滚动计算窗口1内，利用移动平均方法去除压力脉动时间序列中的趋势，本实施例中的移动平均窗口长度取为2t，t为燃烧室最低阶固有振荡周期；
在滚动计算窗口1内，计算去除趋势后压力脉动时间序列的延迟为1自相关系数，并保存计算结果；
如图2所示，确定滚动计算窗口2的长度，本实施例中的滚动计算窗口2长度取为500t，t为燃烧室最低阶固有振荡周期。当滞后为1自相关系数序列长度大于滚动计算窗口2的长度时，计算滞后为1自相关系数序列的kendall秩相关系数。随着发动机工作时间增加，压力脉动时间序列数据也逐渐增加。为了减少计算量，可以间隔数个点计算kendall秩相关系数，本实施例中的间隔点数为20。(步骤四计算kendall秩相关系数，默认是逐点计算，这儿为减少计算量，可以间隔数个点计算kendall秩相关系数)
确定判断阈值，本实施例中的阈值取值为0.5。如果某一时刻kendall秩相关系数τρ大于阈值，说明燃烧室内即将出现燃烧不稳定性，发出警告信号，同时采取控制预案。
1.一种燃烧室燃烧不稳定的在线预测方法，其特征在于，包括以下步骤：
步骤一、获取燃烧室的压力脉动时间序列
通过动态压力传感器获取燃烧室或燃烧室喷前的压力脉动时间序列
步骤二、在滚动计算窗口1内，去除压力脉动时间序列
步骤三、在滚动计算窗口1内，计算压力脉动时间序列的延迟为1自相关系数ρ1，随着计算窗口的滚动，得到各个窗口延迟为1自相关系数ρ1组成的延迟为1自相关系数序列ρ；
延迟为1自相关系数ρ1按下式计算，
其中，e为数学期望计算，pt 1为pt延迟为1的时间序列，μ为pt的平均值或数学期望，
步骤四、设定滚动计算窗口2的长度为m，当延迟为1自相关系数序列ρ数量大于滚动计算窗口2的长度时，在滚动计算窗口2内，计算延迟为1自相关系数序列ρ的kendall秩相关系数；
kendall秩相关系数τρ按下式计算，
其中，ρi和ρj均为滚动计算窗口2内延迟为1自相关系数序列ρ中的点，i和j为自相关系数序列ρ中的点的时间序号；
步骤五、根据阈值判断状态；
若kendall秩相关系数τρ大于阈值τ*，说明燃烧室内即将出现燃烧不稳定。
2.根据权利要求1所述的燃烧室燃烧不稳定的在线预测方法，其特征在于：步骤二中，在滚动计算窗口1内，利用移动平均方法去除压力脉动时间序列
2.1)计算趋势mt，按下式计算；
其中，k为移动平均窗口长度，
2.2)去除原始压力脉动时间序列中的趋势，按下式计算，
其中，pt为去除趋势后的压力脉动时间序列。
3.根据权利要求1所述的燃烧室燃烧不稳定的在线预测方法，其特征在于：步骤二中，通过高斯光滑法或loess光滑法去除压力脉动时间序列中的趋势。
4.根据权利要求1或2或3所述的燃烧室燃烧不稳定的在线预测方法，其特征在于：步骤三中，方差
其中，n为滚动计算窗口1的长度。
5.根据权利要求4所述的燃烧室燃烧不稳定的在线预测方法，其特征在于：步骤五中，阈值为0.5。
技术总结