非線性燃燒不穩(wěn)定振蕩幅值演化特征信號(hào)辨識(shí)方法研究①

甘 林,張翔宇,金秉寧*,魏少娟,劉佩進(jìn)

(1.西北工業(yè)大學(xué) 燃燒、熱結(jié)構(gòu)與內(nèi)流場(chǎng)重點(diǎn)實(shí)驗(yàn)室,西安 710072;2.中國(guó)航天科技集團(tuán)有限公司四院四十一所,西安 710025;3.空軍工程大學(xué) 裝備管理與無(wú)人機(jī)工程學(xué)院,西安 710000)

0 引言

非線性燃燒不穩(wěn)定問(wèn)題一直是固體火箭發(fā)動(dòng)機(jī)研制過(guò)程中常遇到的棘手問(wèn)題之一[1-5],其典型特征是在發(fā)動(dòng)機(jī)燃燒中產(chǎn)生非線性壓強(qiáng)振蕩現(xiàn)象[6-10]。通常,這種非線性壓強(qiáng)振蕩中包含一個(gè)或多個(gè)聲振蕩模態(tài),且各階模態(tài)受增益阻尼因素的影響,能量從低階模態(tài)向高階模態(tài)傳遞,激發(fā)出高階模態(tài),使得不穩(wěn)定問(wèn)題愈加嚴(yán)重[11-15]。然而,目前仍存在對(duì)能量傳遞演化機(jī)理認(rèn)識(shí)不清的問(wèn)題,無(wú)法對(duì)非線性燃燒不穩(wěn)定的預(yù)估和抑制提供理論指導(dǎo)。因此,需要開(kāi)展能量傳遞演化機(jī)理的研究,而進(jìn)行這項(xiàng)研究的前提則是從非線性振蕩數(shù)據(jù)中準(zhǔn)確的得到各階模態(tài)的關(guān)鍵參數(shù),如振蕩頻率和幅值等。

對(duì)于固體火箭發(fā)動(dòng)機(jī)非線性不穩(wěn)定振蕩數(shù)據(jù),科研人員多采用商業(yè)軟件平臺(tái)進(jìn)行分析,其特點(diǎn)是功能強(qiáng)大,操作簡(jiǎn)單,經(jīng)過(guò)長(zhǎng)期探索,可以形成一套數(shù)據(jù)分析流程。但在使用過(guò)程中發(fā)現(xiàn),商業(yè)軟件中很多模塊方法和參數(shù)不明確,使得數(shù)據(jù)處理結(jié)果出現(xiàn)失真;此外,受主觀經(jīng)驗(yàn)判斷的影響,分析結(jié)果也存在一定的誤差。在處理類(lèi)似信號(hào)如電力系統(tǒng)低頻振蕩信號(hào)、機(jī)械系統(tǒng)振動(dòng)衰減信號(hào)等方面,更多的是采用時(shí)頻分析方法,如經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)[16]和變分模態(tài)分解(VMD)[17-18]等,并且取得了較好的效果。綜上,有必要借鑒處理類(lèi)似信號(hào)的分析方法建立一種方便而準(zhǔn)確的分析固體火箭發(fā)動(dòng)機(jī)非線性燃燒不穩(wěn)定振蕩數(shù)據(jù)的方法,以獲得非線性不穩(wěn)定振蕩數(shù)據(jù)的關(guān)鍵參數(shù)。

固體火箭發(fā)動(dòng)機(jī)非線性燃燒不穩(wěn)定振蕩數(shù)據(jù)是典型的非線性信號(hào)。在非線性信號(hào)處理方面,典型的時(shí)頻方法有經(jīng)驗(yàn)?zāi)B(tài)分解(EMD)[19]、局部均值分解(LMD)[20]、本征時(shí)間尺度分解(ITD)[21]和變分模態(tài)分解(VMD)[22]等。其中,前三種基于遞歸迭代的思想,能夠?qū)π盘?hào)模態(tài)進(jìn)行分解。然而,EMD存在著一些問(wèn)題,如缺乏嚴(yán)格的數(shù)學(xué)公式,易受模態(tài)混疊、端點(diǎn)效應(yīng)、過(guò)擬合和欠采樣的影響;LMD在抑制端點(diǎn)效應(yīng)和模態(tài)混疊方面優(yōu)于EMD,但其主要問(wèn)題是計(jì)算復(fù)雜度高,仍容易產(chǎn)生模態(tài)混疊;ITD存在末端效應(yīng)和時(shí)間分辨率低的問(wèn)題;而VMD方法基于求解一個(gè)完全非遞歸變分問(wèn)題的思想,具有理論上嚴(yán)謹(jǐn)?shù)臄?shù)學(xué)原理和對(duì)噪聲及采樣的魯棒性等優(yōu)點(diǎn)[23]。但VMD方法依賴于經(jīng)驗(yàn),需要人工預(yù)先設(shè)定模態(tài)分解個(gè)數(shù)K、二次懲罰項(xiàng)α和帶寬τ,參數(shù)設(shè)置不當(dāng)可能導(dǎo)致模態(tài)混疊,影響分解結(jié)果的準(zhǔn)確性[24]。

綜合考慮,本文引入了VMD[24]算法,并采用標(biāo)準(zhǔn)粒子群優(yōu)化算法(SPSO 2011,以下簡(jiǎn)稱SPSO)[25]對(duì)VMD算法所需的預(yù)置參數(shù)進(jìn)行優(yōu)化,建立了SPSO-VMD方法。將SPSO-VMD方法應(yīng)用于典型特征的燃燒不穩(wěn)定聲振蕩信號(hào)處理,準(zhǔn)確地獲得了各階模態(tài)的頻率和幅值,驗(yàn)證了該方法的可行性。并將其應(yīng)用在了真實(shí)發(fā)動(dòng)機(jī)振蕩數(shù)據(jù)的處理上,準(zhǔn)確地分離出了各階模態(tài)。

1 SPSO-VMD方法的建立

VMD方法將原始信號(hào)X(t)分解為K(模態(tài)分解個(gè)數(shù))個(gè)緊繞其中心頻率ω的本征模態(tài)函數(shù)。在該算法中,定義本征模態(tài)函數(shù)(IMF)為一個(gè)調(diào)幅-調(diào)頻信號(hào),VMD分解過(guò)程實(shí)質(zhì)上是變分問(wèn)題的構(gòu)造和求解過(guò)程,詳見(jiàn)文獻(xiàn)[24]。

SPSO算法是粒子群算法(Particle Swarm Optimization,PSO)[26]的一種改進(jìn)。PSO算法作為一種群體智能優(yōu)化算法,具有原理簡(jiǎn)單、易于實(shí)現(xiàn)和搜索速度快等優(yōu)勢(shì),在PSO算法中,群體中的每個(gè)個(gè)體稱為一個(gè)粒子,代表優(yōu)化問(wèn)題的一個(gè)潛在解。對(duì)于一個(gè)n維優(yōu)化問(wèn)題,每個(gè)粒子被賦予一個(gè)n維速度向量和n維位置向量。在每一次迭代也即搜索的過(guò)程中,每個(gè)粒子根據(jù)自身的搜索經(jīng)驗(yàn)和群體中其他粒子的搜索經(jīng)驗(yàn)動(dòng)態(tài)調(diào)整各自的運(yùn)動(dòng)狀態(tài),從而保證每個(gè)粒子逐步向問(wèn)題的最優(yōu)解方向靠近[27-28]。但是,基本PSO算法存在兩個(gè)缺陷:(1)不能很好地動(dòng)態(tài)調(diào)整粒子的全局和局部搜索能力;(2)粒子容易陷入搜索迭代停滯[29]。RATNAWEERA等[30]為克服PSO算法的第一個(gè)缺陷,提出了線性自適應(yīng)參數(shù)調(diào)整法則,即根據(jù)當(dāng)前迭代次數(shù)和最大迭代次數(shù),線性地調(diào)整參數(shù)的值;CLERC等[25]為克服基本PSO算法過(guò)早陷入迭代停滯的缺陷,提出了SPSO算法。

在SPSO算法中,首先根據(jù)每個(gè)粒子的當(dāng)前位置Xi(t)、粒子的個(gè)體最優(yōu)位置Pbesti(t)以及群體的全局最優(yōu)位置Gbest,為每個(gè)粒子定義一個(gè)超球面;然后,在定義的超球面中,為每個(gè)粒子隨機(jī)產(chǎn)生一個(gè)新的位置點(diǎn);最后,將這個(gè)隨機(jī)產(chǎn)生的位置視為擾動(dòng)添加到粒子的速度更新公式中,所以在迭代初期,SPSO算法可從最大程度上保證每個(gè)粒子的飛行速度不為0?？梢员M量保證在迭代的初期每個(gè)粒子都能以非零速度在解空間中一直搜索,從而避免算法過(guò)早陷入迭代停滯。

本文按照線性自適應(yīng)參數(shù)調(diào)整法則的SPSO算法優(yōu)化VMD分解,建立了SPSO-VMD方法,其思路如下:影響VMD算法的三個(gè)需要提前輸入的參數(shù)為模態(tài)分解個(gè)數(shù)K、二次懲罰因子α和帶寬τ。在處理的過(guò)程中,模態(tài)分解個(gè)數(shù)K是通過(guò)對(duì)振蕩信號(hào)進(jìn)行傅里葉變換并濾波從而確定的;二次懲罰因子α和帶寬τ是通過(guò)SPSO算法來(lái)優(yōu)化的。

SPSO算法優(yōu)化需要一個(gè)適應(yīng)度函數(shù)f(x),其定義如下:首先,對(duì)振蕩信號(hào)進(jìn)行傅里葉變換并濾波,得到K個(gè)模態(tài)頻率及對(duì)應(yīng)的K個(gè)幅值A(chǔ)k;再使用VMD算法對(duì)振蕩信號(hào)進(jìn)行分解得到K個(gè)模態(tài)分量,對(duì)這些模態(tài)分量進(jìn)行傅里葉變換,得到K個(gè)模態(tài)分量對(duì)應(yīng)的K個(gè)幅值A(chǔ)k,VMD(k=1,2,…,K)。各模態(tài)幅值差值ΔA計(jì)算公式為

ΔAk=|Ak-Ak,VMD|

ΔA={ΔA1,ΔA2,…,ΔAk}

(1)

將Max(ΔA)定義為SPSO算法的適應(yīng)度函數(shù)f(x),也即SPSO優(yōu)化VMD算法是朝著值Max(ΔA)變小的方向去優(yōu)化的。

應(yīng)用SPSO算法優(yōu)化VMD算法流程如下:首先設(shè)定兩個(gè)參數(shù)α和τ優(yōu)化的具體范圍,粒子在這個(gè)范圍內(nèi)進(jìn)行全局尋優(yōu),α和τ的值由粒子的位置來(lái)表征,每個(gè)粒子每次循環(huán)將自身的位置代入VMD算法對(duì)振蕩信號(hào)進(jìn)行一次分解,并計(jì)算一次適應(yīng)度值,通過(guò)對(duì)比新粒子的適應(yīng)度值進(jìn)行位置的更新,使每個(gè)粒子向著最優(yōu)解方向更新位置,通過(guò)不斷的迭代更新,最終找到使整個(gè)種群全局適應(yīng)度最優(yōu)的位置,也即所需優(yōu)化參數(shù)的最佳值。

SPSO-VMD方法流程圖如圖1所示。

2 SPSO-VMD方法的驗(yàn)證

2.1 典型特征信號(hào)構(gòu)建

固體火箭發(fā)動(dòng)機(jī)不穩(wěn)定振蕩數(shù)據(jù)中包含振蕩增長(zhǎng)段、極限環(huán)段以及振蕩衰減段。并且在整個(gè)振蕩數(shù)據(jù)中包含多階模態(tài),各模態(tài)之間存在能量傳遞,各模態(tài)幅值變化率也隨時(shí)間變化。因此,構(gòu)造符合上述特征的典型的振蕩信號(hào)用于校驗(yàn)SPSO-VMD模態(tài)分解的準(zhǔn)確性。此時(shí),各階模態(tài)壓強(qiáng)振蕩幅值Rm滿足的非線性微分方程[13]為

(2)

式中αm為m階模態(tài)壓強(qiáng)振蕩的總增益系數(shù);A為非線性系數(shù);ωm為m階聲振模態(tài)的角頻率;Emnl為表示各階模態(tài)之間的非線性耦合作用。

按上式構(gòu)造含20階模態(tài)的仿真信號(hào),前5階模態(tài)的幅值以及幅值變化率如表1所示,模態(tài)數(shù)m=20。6～20階初始幅值均為0,其幅值變化率滿足αm=(m-3)×α4。在表1中,正值αm表示能量通過(guò)線性作用從燃燒場(chǎng)流入振蕩系統(tǒng)中的第m階模態(tài),負(fù)值αm則表示能量從振蕩系統(tǒng)中流出,高階模態(tài)的αm的絕對(duì)值越大,表明高階模態(tài)越難被激發(fā)以及越容易被耗散。按上式構(gòu)造含20階模態(tài)的仿真信號(hào)如圖2所示。從圖2可看出,上升段振蕩幅值出現(xiàn)了明顯線性向非線性轉(zhuǎn)變的過(guò)程,即各階模態(tài)的幅值并不是以線性的方式持續(xù)增長(zhǎng),當(dāng)各階模態(tài)能量傳遞達(dá)到動(dòng)態(tài)平衡過(guò)程時(shí),出現(xiàn)了極限環(huán)狀態(tài);下降段各階模態(tài)的幅值變化率也是隨時(shí)變化各不相同;圖2(b)為整個(gè)振蕩段頻譜分析的前五階結(jié)果,前三階振蕩模態(tài)的幅值占主導(dǎo),所以這里僅分解前三階模態(tài)來(lái)驗(yàn)證分解結(jié)果的準(zhǔn)確性。

表1 前五階模態(tài)的幅值以及幅值變化率

(a)Amplitude of the simulated signal

2.2 典型特征信號(hào)的SPSO-VMD方法分解

對(duì)上述信號(hào)在給定α和τ范圍為(0,3000]和[0,0.2]內(nèi)進(jìn)行SPSO-VMD方法的分解,優(yōu)化所得α和τ分別為2983.78和0.000 4。各階模態(tài)幅值分解結(jié)果如圖3(b)所示。從圖3(b)可以看分解后的各階模態(tài)幅值與圖2(c)幅值輪廓線高度一致;表2展示了SPSO-VMD方法分解該特征信號(hào)的頻率和幅值結(jié)果。

表2 SPSO-VMD方法分解特征信號(hào)頻率、幅值結(jié)果

(a)Iterative process

這里三個(gè)模態(tài)的分解所得極限環(huán)幅值是通過(guò)統(tǒng)計(jì)的方式計(jì)算所得,其與標(biāo)準(zhǔn)極限環(huán)幅值誤差均低于0.5%。從表中結(jié)果可以看出,對(duì)于非線性燃燒不穩(wěn)定典型特征信號(hào),SPSO-VMD方法能夠準(zhǔn)確的分離各個(gè)模態(tài)、辨識(shí)各階頻率和幅值。比較了以下兩種方法:

(1)采用VMD算法對(duì)特征信號(hào)進(jìn)行分解,分別取值α和τ為16和0.01;

(2)采用PSO算法以同樣思路在同樣參數(shù)范圍內(nèi)對(duì)特征信號(hào)進(jìn)行優(yōu)化VMD參數(shù)的分解,對(duì)比結(jié)果如表3和圖4所示。

圖4 Mode 1極限環(huán)段振蕩某點(diǎn)誤差示意圖

表3 三種算法對(duì)比

為求出全局最小適應(yīng)度函數(shù)值f(x),采用VMD法對(duì)特征信號(hào)進(jìn)行分解。其中,α和τ的取值范圍分別為(0,3000]和[0,0.2],α的步長(zhǎng)為1,τ的步長(zhǎng)為0.01,全局求解共耗時(shí)約10 h,結(jié)果如圖5所示。

圖5 f(x)的全局求解

如表3和圖4所示,PSO-VMD方法和SPSO-VMD方法分解的誤差遠(yuǎn)小于VMD方法。從圖4中可知,SPSO-VMD方法的誤差略小于PSO-VMD方法,但SPSO-VMD方法可以全局優(yōu)化到與最小適應(yīng)度函數(shù)值相對(duì)應(yīng)的參數(shù),且SPSO-VMD方法所需時(shí)間更短,效率更高,如圖3(a)所示,在本次計(jì)算中,SPSO-VMD方法僅用3步就取得了比PSO方法更好的效果。通過(guò)兩者誤差的對(duì)比分析,可以認(rèn)為在適應(yīng)度函數(shù)為10-5數(shù)量級(jí)時(shí),分解便已經(jīng)達(dá)到了較好的效果。

3 SPSO-VMD方法的應(yīng)用

某發(fā)動(dòng)機(jī)在工作過(guò)程中發(fā)生了非線性燃燒不穩(wěn)定,其振蕩幅值數(shù)據(jù)如圖6(a)所示,對(duì)其進(jìn)行FFT處理。頻譜結(jié)果表明,該發(fā)動(dòng)機(jī)在工作過(guò)程中出現(xiàn)了典型的非線性燃燒不穩(wěn)定現(xiàn)象,存在多階模態(tài),且前8階模態(tài)幅值明顯,通過(guò)計(jì)算分析,分解前8階模態(tài)重構(gòu)該信號(hào)的精度能達(dá)到95%,于是這里使用SPSO-VMD方法對(duì)其振蕩數(shù)據(jù)進(jìn)行前8階模態(tài)分解,所用計(jì)算資源CPU型號(hào)為intel Xeon(R)Platinum 8259CL、48核,迭代在第9步時(shí),適應(yīng)度函數(shù)便達(dá)到了7.763 5×10-5,通過(guò)上一節(jié)的誤差分析,可以認(rèn)為該適應(yīng)度值表明分解已經(jīng)達(dá)到了較好的效果,對(duì)應(yīng)用時(shí)約為1880 s。

(a)Oscillation signal

8階模態(tài)的分解結(jié)果如圖7(a)～(h)所示,各分解模態(tài)的FFT結(jié)果表明,各階模態(tài)被準(zhǔn)確分離。

(a)Mode 1 and its FFT result (b)Mode 2 and its FFT result

針對(duì)于類(lèi)似第6階模態(tài),頻譜圖上出現(xiàn)兩個(gè)臨近頻率峰值的情況,這里作一定的分析:對(duì)第6階模態(tài)進(jìn)行STFT處理,所得結(jié)果如圖8所示,結(jié)果表明第6階模態(tài)在發(fā)動(dòng)機(jī)發(fā)生燃燒不穩(wěn)定的時(shí)間內(nèi)出現(xiàn)了頻率的變化,頻譜圖中兩個(gè)峰值恰好對(duì)應(yīng)STFT結(jié)果中振蕩幅值最大的0.29 s和0.46 s兩個(gè)時(shí)刻下的頻率,也即這兩個(gè)頻率是第6階模態(tài)不同時(shí)刻的振蕩頻率,它們屬于一個(gè)模態(tài)。這種頻率變化的現(xiàn)象主要是由燃燒室裝藥燃燒,燃面退移,導(dǎo)致燃燒室內(nèi)聲腔構(gòu)型發(fā)生變化造成的。

圖8 第6階模態(tài)STFT處理結(jié)果 圖9 前8階模態(tài)重構(gòu)信號(hào)

這里不考慮9階及其高階模態(tài),將得到的前8階模態(tài)重新疊加,與原始的振蕩幅值進(jìn)行對(duì)比,如圖9所示,計(jì)算前8階模態(tài)重構(gòu)振蕩幅值的精度為95.2%。如圖10所示,將前8階模態(tài)繪制在同一張圖中,在0.16～0.4 s時(shí)間段范圍內(nèi),壓強(qiáng)信號(hào)雖然保持極限環(huán)振蕩現(xiàn)象,但是各階模態(tài)的幅值卻在不斷變化,可看到奇數(shù)階模態(tài)幅值衰減到了一個(gè)極小值,而偶數(shù)階模態(tài)出現(xiàn)了一定的增長(zhǎng),仿佛出現(xiàn)了偶數(shù)階增長(zhǎng)抑制奇數(shù)階的現(xiàn)象,據(jù)文獻(xiàn)[15]所提:“非線性能量可從一階聲模態(tài)傳遞至所有高階聲模態(tài),而從二階聲模態(tài)傳遞只能向偶數(shù)階傳遞,奇數(shù)階模態(tài)均不發(fā)生能量傳遞過(guò)程,且高階聲模態(tài)不能向低階聲模態(tài)傳遞能量。”在該發(fā)動(dòng)機(jī)振蕩信號(hào)中,由于某種原因,一階模態(tài)幅值衰減到了接近0的值,傳遞向高階聲模態(tài)的能量極少,而二階聲模態(tài)幅值較大,不斷向偶數(shù)階模態(tài)傳遞著非線性能量,導(dǎo)致偶數(shù)階維持一定幅值的振蕩,而奇數(shù)階高階模態(tài)因?yàn)檩斎氲姆蔷€性能量較少,不足以維持自身耗散的作用,導(dǎo)致幅值一直處在一個(gè)較低的值。這種現(xiàn)象恰好也提示研究人員可以從能量傳遞的角度對(duì)發(fā)動(dòng)機(jī)非線性燃燒不穩(wěn)定進(jìn)行抑制,比如切斷能量傳輸?shù)穆窂交蛘吒淖兡芰總鬏數(shù)牧康鹊?這些都有待研究。

圖10 前8階模態(tài)奇偶數(shù)階幅值演化對(duì)比

4 結(jié)論

(1)SPSO-VMD方法能夠準(zhǔn)確地分解發(fā)動(dòng)機(jī)非線性燃燒不穩(wěn)定典型特征信號(hào)的各階模態(tài),獲得各模態(tài)的頻率和幅值,可應(yīng)用于真實(shí)發(fā)動(dòng)機(jī)燃燒不穩(wěn)定振蕩信號(hào)的處理,可為研究非線性不穩(wěn)定能量傳遞演化機(jī)制提供數(shù)據(jù)支撐。

(2)SPSO-VMD方法改進(jìn)了PSO-VMD方法陷入局部最優(yōu)解的缺陷,能夠準(zhǔn)確地找到全局最優(yōu)解,分解模態(tài)誤差較小,分解效率較高。

(3)SPSO-VMD方法自適應(yīng)強(qiáng),不依賴主觀經(jīng)驗(yàn),集成性好,使用簡(jiǎn)單方便,可適用于工業(yè)部門(mén)對(duì)發(fā)動(dòng)機(jī)燃燒不穩(wěn)定問(wèn)題的定位和分析。