楊尚榮, 吳林龍, 于 涵, 楊寶娥
(西安航天動力研究所 液體火箭發(fā)動機(jī)技術(shù)重點(diǎn)實(shí)驗(yàn)室, 西安 710100)
燃燒不穩(wěn)定是航空航天發(fā)動機(jī)研制中需要重點(diǎn)關(guān)注的問題,通常表現(xiàn)為燃燒室聲模態(tài)頻率附近的高幅值壓力振蕩。國內(nèi)外在其產(chǎn)生機(jī)理、評估和預(yù)測方法、以及控制手段等方面開展了大量的研究[1-2]。目前預(yù)測方法在工程上的應(yīng)用還未成熟,尚不能在發(fā)動機(jī)試驗(yàn)前對其穩(wěn)定性進(jìn)行可靠的預(yù)測。因此一旦出現(xiàn)燃燒不穩(wěn)定,通常從兩方面出發(fā)去解決:一方面減少激勵能量;另一方面增加燃燒室聲學(xué)阻尼。前者需要考慮燃燒和聲場間的耦合作用過程,由于燃燒現(xiàn)象本身的復(fù)雜性,實(shí)施改進(jìn)方案難度較大。相比較而言,被動阻尼裝置的設(shè)計(jì)要容易一些。
從線性角度,增長率表征小擾動下燃燒室振蕩能量的變化率。如果其值大于零,說明系統(tǒng)不穩(wěn)定,反之,說明系統(tǒng)穩(wěn)定,此時增長率也稱為衰減率。對于出現(xiàn)燃燒不穩(wěn)定的燃燒室,新設(shè)計(jì)阻尼裝置的衰減率應(yīng)大于燃燒室的增長率并保留一定的裕度,即需要燃燒室不穩(wěn)定聲模態(tài)的增長率作為設(shè)計(jì)輸入?yún)?shù)。
學(xué)者們提出了多種方法來識別增長率。主要分為外部激勵方法(輸入-輸出識別)和無外部激勵方法(僅輸出識別)。外部激勵方法中,對于Hopf類型的分叉,Lee等[3]利用分叉前的噪聲相干共振特性[4],識別系統(tǒng)的非線性燃燒響應(yīng),外推出系統(tǒng)分叉點(diǎn)的位置和類型、增長率、極限環(huán)振蕩幅值等,進(jìn)一步可以對系統(tǒng)穩(wěn)定性邊界進(jìn)行預(yù)測[5]。另一類方法采用主動控制[6],在不穩(wěn)定工況點(diǎn)施加和關(guān)閉控制,獲得振蕩幅值從小擾動增加到極限環(huán)的過程。利用振蕩幅值初始增長階段(滿足線性假設(shè))數(shù)據(jù),擬合指數(shù)函數(shù)或采用其他數(shù)據(jù)處理方法(如DMD(dynamic mode decomposition)[7])識別增長率。實(shí)際發(fā)動機(jī)試驗(yàn)中加入特定形式的外部激勵或主動控制比較困難,因此外部激勵方法在工程應(yīng)用中不方便。
無外部激勵的方法,第一類與主動控制方法中關(guān)閉控制后的處理方法相同,區(qū)別是利用自發(fā)燃燒不穩(wěn)定振蕩幅值初始增長階段數(shù)據(jù)[8]。另一類方法基于噪聲激勵作用下燃燒室聲模態(tài)的隨機(jī)動力學(xué)理論[9],從極限環(huán)振蕩數(shù)據(jù)中提取線性增長率。根據(jù)噪聲強(qiáng)度等級,Noiray等[10]提出對應(yīng)的3種識別方法,并利用燃?xì)廨啓C(jī)中壓力振蕩數(shù)據(jù)進(jìn)行了比較。進(jìn)一步,利用數(shù)值模擬[11]和燃燒器試驗(yàn)[12]對第三種方法進(jìn)行了驗(yàn)證。Bonciolini等[13]比較了彩色噪聲和白噪聲的激勵效果,發(fā)現(xiàn)在振蕩頻率附近采用合適的帶通濾波,兩者的區(qū)別可以忽略。若帶通濾波過窄,則兩者的識別結(jié)果差別較大。Boujo等[14]通過引入伴隨優(yōu)化方法解決了上述問題。即使考慮燃燒響應(yīng)相對于壓力振蕩的時間延遲,該方法依然可以較準(zhǔn)確的識別線性增長率[15]。Hummel等[16]將方法推廣到周向不穩(wěn)定對聲模態(tài)增長率的識別。
本文開展同軸離心式噴嘴穩(wěn)定性試驗(yàn)研究,采用無外部激勵方法識別不同工況下的線性增長率,并對其參數(shù)敏感性開展分析,為液體火箭發(fā)動機(jī)燃燒穩(wěn)定性分析提供參考。
燃燒試驗(yàn)系統(tǒng)的詳細(xì)介紹參見文獻(xiàn)[17]。模擬燃燒室為敞口圓筒形結(jié)構(gòu)。試驗(yàn)用噴嘴(如圖1所示)為帶縮進(jìn)室的同軸離心噴嘴,外噴嘴為煤油離心噴嘴,內(nèi)噴嘴為直流氧化劑噴嘴。
1. 直流噴嘴; 2. 離心噴嘴; 3. 節(jié)流嘴; 4. 富氧氣通道; 5. 縮進(jìn)室。圖1 帶縮進(jìn)室的同軸離心噴嘴Fig.1 Recessed swirl coaxial injector
表1 試驗(yàn)工況Tab.1 Test operating conditions
當(dāng)燃燒室只有單階聲模態(tài)振蕩時,可以不考慮模態(tài)間的耦合作用。壓力振蕩可以用非線性振動方程描述
式中:η為脈動壓力;f為由氣動力和燃燒過程產(chǎn)生的激勵源項(xiàng);α為聲模態(tài)衰減率;ω0為聲模態(tài)角頻率;ξ為燃燒噪聲,假定為高斯白噪聲,有〈ξξτ〉=Γδ(τ),Γ為噪聲強(qiáng)度。
(2)
假定壓力振蕩接近諧振,可以將脈動壓力寫成如下幅頻形式
η(t)=A(t)cos[ωt+σ(t)]
(3)
式中,A(t)和σ(t)分別為壓力振蕩的幅值和相位。代入式(2),利用隨機(jī)平均方法,式(2)可以寫成幅值A(chǔ)和相位σ的一階隨機(jī)微分方程。由于A和σ可以解耦,僅列出幅值A(chǔ)滿足的方程
(4)
(5)
駐定情況下,FP方程存在解析解
(6)
(7)
進(jìn)一步通過非線性擬合識別出線性增長率ν、飽和系數(shù)κ和噪聲強(qiáng)度Γ。
不同混合比條件下的壓力時間序列及其對應(yīng)的相空間軌跡如圖2所示。當(dāng)混合比為15.0時,壓力振蕩幅值較小,認(rèn)為是燃燒噪聲階段,相軌跡圖為無序軌跡線的聚集。當(dāng)混合比為17.5時,振蕩幅值增加,振幅波動也較大,軌跡在不同幅值周期振蕩之間轉(zhuǎn)變,中心為空心。當(dāng)混合比為20.0和22.5時,達(dá)到極限環(huán)振蕩,軌跡圖為圓環(huán)面,環(huán)的寬帶代表噪聲對幅值的影響。時頻分析顯示振蕩燃燒頻率在(2 340±100)Hz內(nèi)。
圖2 壓力時間序列和相空間軌跡圖Fig.2 Pressure time series and phase portrait
試驗(yàn)中壓力數(shù)據(jù)有單位,如Pa,kPa等。本節(jié)分析壓力數(shù)據(jù)取不同單位時對識別結(jié)果的影響。令幅值A(chǔ)(t)=kB(t),k為單位間的換算關(guān)系,B(t)為換算后的幅值時間序列。代入式(5),化簡并結(jié)合可得
(8)
比較式(5)和式(8),可知單位變換后,增長率不變,飽和系數(shù)為原來的k2倍,噪聲強(qiáng)度為原來的1/k2倍。
表2 不同幅值下參數(shù)識別結(jié)果
本文壓力單位選用kPa進(jìn)行數(shù)據(jù)分析,其他單位下的飽和系數(shù)和噪聲強(qiáng)度可通過2.2節(jié)換算關(guān)系得到。4種工況下壓力時間序列識別結(jié)果如圖3所示。隨著混合比增加(總流量也增加),壓力振蕩均方根(root mean square,RMS)幅值增加,增長率ν由負(fù)值轉(zhuǎn)變?yōu)檎?飽和系數(shù)κ下降,噪聲強(qiáng)度Γ增加。具體地,當(dāng)混合比為15.0時,其增長率為負(fù)值,系統(tǒng)是穩(wěn)定的。當(dāng)混合比為17.5時,增長率也為負(fù)值,從線性角度分析系統(tǒng)應(yīng)該是穩(wěn)定的。但從圖2壓力曲線和相軌跡可知,其周期性振蕩幅值也較大。說明存在噪聲條件下(實(shí)際發(fā)動機(jī)中燃燒噪聲不可避免),不僅要保證系統(tǒng)增長率為負(fù),還需要留有一定裕度,才能保證系統(tǒng)不發(fā)生大幅值周期振蕩。當(dāng)混合比為20.0和22.5時,增長率均為正值,且幅值越高,增長率越大。該結(jié)果不能作為通用結(jié)論來推廣,理論上,以式(2)作為控制方程的系統(tǒng),壓力振蕩幅值由增長率ν和飽和系數(shù)κ共同決定。不同工況下,飽和系數(shù)κ非定值。此外試驗(yàn)中也發(fā)現(xiàn)較高幅值的壓力振蕩并非一定對應(yīng)較大的增長率[20]。
圖3 試驗(yàn)壓力時間序列RMS值和識別結(jié)果Fig.3 RMS values of test pressure time series and identification results
將識別結(jié)果代入FP方程解析解式(6),與試驗(yàn)壓力振蕩幅值概率密度分布進(jìn)行比較,如圖4所示。當(dāng)混合比為20時,試驗(yàn)結(jié)果和識別結(jié)果差別較大,其他工況結(jié)果相對較好,一定程度上驗(yàn)證了燃燒響應(yīng)模型和識別方法的合理性。Oriray研究中的數(shù)值試驗(yàn)發(fā)現(xiàn)增長率在-5~5 rad/s內(nèi)時,準(zhǔn)確識別需要約30 s時長的數(shù)據(jù)。本文試驗(yàn)采樣時間T=5 s,相比準(zhǔn)確識別所需的數(shù)據(jù)長度而言較短,因此需要對識別結(jié)果的不確定性進(jìn)行分析。
圖4 振蕩幅值概率密度函數(shù)Fig.4 Probability density function of oscillation amplitude
本節(jié)分析識別結(jié)果的不確定性以及數(shù)據(jù)采樣時間和采樣頻率對識別結(jié)果的影響。每個工況下試驗(yàn)數(shù)據(jù)只有一次,且系統(tǒng)真實(shí)的增長率未知。因此,為了量化識別方法的不確定性,將每個工況下的識別值作為設(shè)定值,數(shù)值求解隨機(jī)微分式(2)生成100組壓力時間序列,識別其增長率并取平均值和標(biāo)準(zhǔn)差,來量化該工況下識別結(jié)果的不確定性。
4個工況下的計(jì)算結(jié)果如圖5所示。采樣時間T和采樣率sf與試驗(yàn)中取值相同,即T=5 s,sf=25.6 kHz。結(jié)果表明,當(dāng)混合比為20.0和22.5時,增長率均值與設(shè)定值偏差分別為11%和15%,標(biāo)準(zhǔn)差分別為17%和20%。當(dāng)混合比為15.0和17.5時,增長率均值與設(shè)定值偏差分別為76%和50%,標(biāo)準(zhǔn)差分別為62%和56%。工程應(yīng)用時可以確定一個閾值,設(shè)定值(試驗(yàn)數(shù)據(jù)識別值)與識別值偏差小于該閥值時,則可以將識別值的標(biāo)準(zhǔn)差作為設(shè)定值的標(biāo)準(zhǔn)差,從而確定識別結(jié)果的不確定性。
圖5 采樣率對增長率識別的影響,采樣時間T=5 sFig.5 Effect of sampling rate on growth rate identification, sampling time T=5 s
采樣頻率對結(jié)果的影響見圖5,采樣時間均為5 s,采樣率分別為25.6 kHz,50 kHz,1 000 kHz。均值和標(biāo)準(zhǔn)差均由式(2)生成的100組壓力時間序列計(jì)算得到。對識別值的橫坐標(biāo)作了平移,以便圖形能清晰表達(dá)。增加采樣率,當(dāng)混合比為15.0和17.5的偏差(設(shè)定值與識別均值之差)減小,混合比20.0和22.5的偏差在采樣率50 kHz時最小,4個混合比下其識別值的標(biāo)準(zhǔn)差均增加。
采樣時間對結(jié)果的影響如圖6所示,采樣頻率均為1 000 kHz,采樣時間分別為1 s,3 s,5 s。每個均值和標(biāo)準(zhǔn)差均由式(2)生成的100組壓力時間序列計(jì)算得到。同樣對識別值的橫坐標(biāo)作了平移。從圖6可知,增加采樣時間,4個混合比下識別值的偏差和標(biāo)準(zhǔn)差均減小。由此可知,為了增加識別的精度,應(yīng)該增加采樣時間,同時提高采樣頻率,以便在數(shù)據(jù)分析時可選出對應(yīng)工況下較優(yōu)的采樣頻率。
圖6 采樣時間對增長率識別的影響,采樣率sf=1 000 kHzFig.6 Effect of sampling time on growth rate identification, sampling rate sf=1 000 kHz
增長率也可以從壓力振蕩幅值線性增長階段識別得到。由于壓力振蕩幅值在線性增長階段以指數(shù)形式增長[21],因此可對壓力振蕩時間序列取包絡(luò),利用指數(shù)函數(shù)去擬合,或?qū)j(luò)取對數(shù)后利用線性函數(shù)去擬合。本章分析噪聲對增長率識別結(jié)果的影響。
圖7 無噪聲壓力時間序列和增長率識別Fig.7 Noise free pressure time series and growth rate identification
圖8 含噪聲壓力時間序列和增長率識別Fig.8 Pressure time series with noise and growth rate identification
圖9 噪聲對增長率識別的影響Fig.9 Effect of noise on growth rate identification
(1)壓力振蕩幅值的單位不改變增長率的值,但會影響飽和系數(shù)和噪聲強(qiáng)度,理論上得到了相似關(guān)系并獲得了數(shù)值試驗(yàn)的驗(yàn)證。
(2)即使增長率小于零,如果噪聲強(qiáng)度較大,系統(tǒng)也有可能發(fā)生較大幅值的壓力振蕩,因此系統(tǒng)增長率應(yīng)該留有一定裕度。
(3)提高采樣時間可以提高識別精度,不同工況下存在各自較優(yōu)的采樣率。
(4)當(dāng)噪聲強(qiáng)度較大時,從振幅線性增長階段識別增長率存在較大誤差,需要考慮噪聲導(dǎo)致的不確定性。