同軸離心式噴嘴熱聲不穩(wěn)定性遞歸網(wǎng)絡分析

楊尚榮, 周晨初 ,王牧昕

(1.西安航天動力研究所 液體火箭發(fā)動機技術重點實驗室, 陜西 西安 710100； 2.中國航天科工集團公司六院210所，陜西 西安 710065)

0 引言

熱聲系統(tǒng)中火焰、流場和聲場間存在復雜的耦合關系，因此，熱聲系統(tǒng)可以作為復雜系統(tǒng)處理[1]。復雜網(wǎng)絡方法常用來研究復雜系統(tǒng)[2]，研究思路為首先將研究對象轉換為網(wǎng)絡結構形式(僅包含節(jié)點和邊)，然后利用網(wǎng)絡測度(包括節(jié)點、邊、全局測度)進行量化分析。近來被應用到熱聲系統(tǒng)的動力學研究中，可以獲得更多非線性動力學特性方面的信息，此外，基于復雜網(wǎng)絡方法，可以發(fā)展出熱聲系統(tǒng)中熱聲振蕩和貧燃熄滅的預測方法[3]。

Murugesan利用可視化網(wǎng)絡(visibility networks)方法將壓力時間序列轉變?yōu)閺碗s網(wǎng)絡，發(fā)現(xiàn)燃燒噪聲具有無標度(scale—free)網(wǎng)絡特征[4]，而熱聲不穩(wěn)定具有規(guī)則的分布特征。王金華利用可視化網(wǎng)絡構建了湍流火焰面的拓撲結構，標記出對湍流火焰面有較大影響的關鍵褶皺結構，分析了湍流與火焰的相互作用規(guī)律[5]。Okuno采用周期網(wǎng)絡(cycle networks)和k—鄰近網(wǎng)絡(k—nearest neighbors)研究了燃氣輪機模型燃燒器中的熱聲不穩(wěn)定現(xiàn)象，發(fā)現(xiàn)熱聲不穩(wěn)定具有準周期和高維度特征，節(jié)點度顯示出冪律型分布，同時較小的平均路徑長度和較大的聚類系數(shù)說明系統(tǒng)存在小世界(small—world)特征[6]。Gotoda利用改進的可視化網(wǎng)絡和遞歸網(wǎng)絡研究了湍流受限火焰的熄滅過程，發(fā)現(xiàn)接近火焰熄滅時的燃燒狀態(tài)具有小世界特征，網(wǎng)絡的平均度可以作為火焰熄滅的在線預測工具[7]。Murugesan在鈍體穩(wěn)焰湍流燃燒器的動力學研究中，發(fā)現(xiàn)復雜網(wǎng)絡的一些測度，如特征路徑長度、聚類系數(shù)、網(wǎng)絡直徑、全局有效性可作為熱聲振蕩和貧燃熄滅的預測指標[8-9]。Godavarthi利用遞歸網(wǎng)絡方法[10]，發(fā)現(xiàn)特征路徑長度、介數(shù)中心性可對湍流燃燒器中的動力學轉變過程進行預測。Murayama利用旋渦間相互作用構造了湍流網(wǎng)絡，研究了預混火焰燃燒噪聲階段速度場的動力學特征，發(fā)現(xiàn)其具有無標度網(wǎng)絡特征[11]。

同軸離心式噴嘴是液體火箭發(fā)動機常用的雙組元噴嘴之一，對其燃燒穩(wěn)定性的研究一直受到學者們關注。Miller開展了同軸離心式噴嘴燃燒試驗，發(fā)現(xiàn)對于特定噴嘴，燃燒不穩(wěn)定只發(fā)生在一定的頻率范圍內(nèi)[12]。Wierman研究了氣中心同軸離心噴嘴的燃燒響應，發(fā)現(xiàn)氧噴嘴長度對振蕩幅值的影響較大[13]。王楓研究了氣液同軸噴注器的結構尺寸和工作參數(shù)對燃燒穩(wěn)定性的影響，發(fā)現(xiàn)噴嘴縮進長度對燃燒穩(wěn)定性裕量影響很大且存在相對最佳值[14]。王延濤試驗發(fā)現(xiàn)氣氣同軸離心噴嘴的自發(fā)燃燒不穩(wěn)定過程出現(xiàn)了“滯后”現(xiàn)象[15]。王迪開展了氣液同軸離心式單噴嘴模型發(fā)動機燃燒試驗，發(fā)現(xiàn)增加縮進長度會對縱向燃燒不穩(wěn)定產(chǎn)生阻尼作用，但在某些特定的燃燒室長度下，上述阻尼作用可以忽略[16]。

本文開展同軸離心式噴嘴穩(wěn)定性試驗研究，采用遞歸網(wǎng)絡方法分析同軸離心式噴嘴火焰動力學特性，為液體火箭發(fā)動機穩(wěn)定性評估和預測提供參考。

1 實驗系統(tǒng)和數(shù)據(jù)處理方法

1.1 實驗裝置

燃燒試驗系統(tǒng)原理如圖1所示[15]，系統(tǒng)用來進行單噴注器燃燒模擬實驗研究[17]。

1—噴嘴；2—支架；3—模擬噴注器面板；4—模擬燃燒室；5、6—換熱器；7—通風系統(tǒng)；8—脈動壓力傳感器。圖1 燃燒試驗系統(tǒng)原理圖Fig.1 Schematic of test bench

模擬燃燒室為敞口圓筒形結構，直徑156 mm，高277 mm，垂直安放在模擬噴注器面板上。試驗用噴嘴(見圖2)[15]為帶縮進室的同軸離心噴嘴，外噴嘴為離心煤油噴嘴，出口直徑16.4 mm，內(nèi)噴嘴為直流氧化劑噴嘴，出口直徑13.4 mm，縮進室長度為10.5 mm。噴嘴安裝于模擬燃燒室中間位置，出口與模擬噴注器面板平齊。氧化劑為400 ℃的富氧氣(氧氣88%+空氣12%)，燃料為410 ℃的煤油蒸氣，均由換熱器加熱至相應溫度。

1—直流噴嘴；2—離心噴嘴；3—節(jié)流嘴；4—富氧氣通道；5—縮進室。圖2 帶縮進室的同軸離心噴嘴Fig.2 Recessed coaxial injector

1.2 試驗方法和工況

試驗中，燃料流量不變，逐步增加氧化劑流量至各工況點，如表1所示[18]。脈動壓力傳感器安裝在燃燒室壁面距噴注器面板20 mm的位置(見圖1)，采樣頻率12.8 kHz。

表1 試驗工況點

1.3 遞歸網(wǎng)絡

遞歸網(wǎng)絡是將時間序列轉變?yōu)閺碗s網(wǎng)絡的其中一種方法。對于壓力時間序列{x(i):i=1,2，…,m}，由時間延遲重構技術可以得到延遲向量[19]

式中：d為嵌入維數(shù)，通過虛假最鄰近法[20]確定；τ為延遲時間，取采樣間隔的整數(shù)倍，通過互信息法[20]確定。

式中：Θ為Heaviside函數(shù)；ε稱為遞歸閾值；N=m-(d-1)τ為重構相空間中向量的數(shù)量；2條軌跡間的距離通過L2范數(shù)進行計算。

遞歸網(wǎng)絡的鄰近矩陣元素Aij利用遞歸矩陣元素Rij由下式計算得到[20]

Aij=Rij-δij

式中δij為Kronecker delta。當i=j時，δij=1；當i≠j時，δij=0。

鄰近Aij提供了節(jié)點和節(jié)點間邊的信息，如果Aij=1，則節(jié)點i和j有邊相連；如果Aij=0，則節(jié)點i和j不相連。從定義可以看到，遞歸網(wǎng)絡是無向、對稱、無權重網(wǎng)絡。

2 實驗結果及分析

2.1 壓力時間序列分析

試驗中燃料流量不變，隨著氧化劑流量的增加(混合比增加)，系統(tǒng)從燃燒噪聲轉變?yōu)闊崧暡环€(wěn)定狀態(tài)。圖3中列出了典型工況下(混合比0.81、2.97、4.32)的壓力時間序列?？梢钥吹剑旌媳?.81時系統(tǒng)振蕩幅值很低，為燃燒噪聲狀態(tài)?；旌媳?.97時，系統(tǒng)出現(xiàn)大幅值壓力振蕩和小幅值壓力振蕩間隔出現(xiàn)的現(xiàn)象，稱為陣發(fā)現(xiàn)象?；旌媳?.32時，系統(tǒng)出現(xiàn)了幅值較一致的高幅值壓力振蕩。

不同混合比下的壓力均方根值如圖4所示。在混合比為2.35時，模擬燃燒室出現(xiàn)了自激燃燒不穩(wěn)定，混合比為2.97時，振蕩幅值達到最大，繼續(xù)增加混合比，振蕩幅值有所下降。壓力均方根值常被用來在線評估燃燒系統(tǒng)的穩(wěn)定性。但由于燃燒不穩(wěn)定的幅值在出現(xiàn)不穩(wěn)定之前是未知的，因此將其作為不穩(wěn)定預測的參數(shù)不合適。此外，即使壓力均方根值較高，也有可能是噪聲引起的[10]。因此需要其他可以評估和預測燃燒不穩(wěn)定的參數(shù)。

圖3 不同混合比條件下的壓力時間序列Fig.3 Pressure time series for different mixture ratios

圖4 不同混合比下最大功率譜幅值Fig.4 Maximum amplitudes of PSD for different mixture ratios

2.2 網(wǎng)絡拓撲結構

2.2.1 最優(yōu)延遲時間

采用交互信息法求最優(yōu)延遲時間τ時，取互信息法的第一個局部最小值為最優(yōu)延遲時間。求得混合比0.81、2.97、4.32時最優(yōu)延遲時間τ分別為0.11、0.07、0.07 ms。圖5顯示了混合比2.97時最優(yōu)延遲時間計算結果。

圖5 混合比k=2.97時最優(yōu)延遲時間Fig.5 Optimal delay time for mixture ratio k=2.97

2.2.2 嵌入維數(shù)

獲得最優(yōu)延遲時間后，利用虛假最鄰近法求解嵌入維數(shù)。求得混合比0.81、2.97、4.32時嵌入維數(shù)d分別為15、9、8。圖6顯示了混合比2.97時的最小嵌入維數(shù)計算結果。

圖6 混合比k=2.97時最小嵌入維數(shù)Fig.6 Minimum embedding dimension for mixture ratio k=2.97

2.2.3 網(wǎng)絡拓撲結構

圖7中列出了典型工況下(混合比0.81、2.97、4.32)的遞歸網(wǎng)絡拓撲圖?；旌媳?.81時為燃燒噪聲狀態(tài)，體現(xiàn)出隨機網(wǎng)絡特征?；旌媳?.97時為陣發(fā)狀態(tài)，網(wǎng)絡圖顯示出中央聚集和周期邊界特點，中央聚集表示系統(tǒng)為低幅值振蕩，由遞歸圖分析可知[18]，系統(tǒng)為準周期振蕩，周期邊界代表系統(tǒng)為大幅值極限環(huán)振蕩狀態(tài)?；旌媳?.32時為準周期振蕩狀態(tài)，顯示出寬邊界圓環(huán)特征。可以看出，遞歸網(wǎng)絡拓撲圖體現(xiàn)出了不同燃燒動力學狀態(tài)下吸引子的幾何特征。

圖7 3種混合比下的網(wǎng)絡拓撲結構Fig.7 The topologies of recurrence networks for the mixture ratios

2.3 遞歸網(wǎng)絡分析

采用遞歸網(wǎng)絡中幾種全局測度對系統(tǒng)狀態(tài)轉變過程進行量化估計。分別利用1 000()、1 500()、2 000()個數(shù)據(jù)點計算各網(wǎng)絡測度值，見圖8?？梢钥闯?，除介數(shù)中心性(數(shù)值隨數(shù)據(jù)點個數(shù)增加而向下平移減少，但表現(xiàn)出相同的變化規(guī)律)外，其他網(wǎng)絡測度值均體現(xiàn)出了良好的收斂性。

平均度：單個節(jié)點的度定義為

代表與節(jié)點v相連的節(jié)點的數(shù)量。平均度定義為

平均度正比于遞歸量化分析中的量化指標——遞歸率[21]，遞歸率在文獻[18]中進行了分析，因此本文不再對平均度進行分析。

聚類系數(shù)：聚類系數(shù)是節(jié)點聚集性的測量，量化任意節(jié)點領域內(nèi)的節(jié)點是否互為相鄰節(jié)點。單個節(jié)點的聚類系數(shù)定義為

式中：Nv為與節(jié)點v相連的節(jié)點之間邊的數(shù)量；kv(kv-1)/2表示與節(jié)點v相連的節(jié)點之間最大可能的邊的數(shù)量。整個網(wǎng)絡的平均聚類系數(shù)定義為

圖8 網(wǎng)絡測度隨混合比的變化Fig.8 Variation of recurrence networks measures with mixture ratios

從圖8(a)看到，聚類系數(shù)在陣發(fā)振蕩時值最大，燃燒噪聲時值最小。說明燃燒噪聲階段，節(jié)點間相互獨立性較強，體現(xiàn)出隨機特性。陣發(fā)振蕩存在低幅值準周期振蕩(中央聚集)，準周期振蕩狀態(tài)存在寬邊界圓環(huán)聚集，因此聚類系數(shù)較大。

特征路徑長度：節(jié)點i和節(jié)點j之間的最短路徑長度Li,j定義為連接兩節(jié)點所需要的最少邊數(shù)，對于不相連的節(jié)點，Li,j=0。特征路徑長度L定義為所有點之間最小路徑長度之和與所有節(jié)點之間可能的邊的數(shù)量的比值

從圖8(b)看到，L在陣發(fā)振蕩時值最大，準周期振蕩和燃燒噪聲時較小。陣發(fā)振蕩時，高幅值振蕩節(jié)點和低幅值振蕩節(jié)點間距離遠，路徑長，因此特征路徑長度大。噪聲狀態(tài)網(wǎng)絡總邊數(shù)較少(可從網(wǎng)絡圖的密度看出，與遞歸閾值ε的取值方法有關)，正規(guī)化后特征路徑長度小。

網(wǎng)絡直徑：網(wǎng)絡直徑定義為特征路徑長度的最大值

D=max(Li,j)

從圖8(c)看到，網(wǎng)絡直徑趨勢與特征路徑長度稍有不同，燃燒噪聲時較小，陣發(fā)振蕩和準周期振蕩時值較大，說明燃燒噪聲時節(jié)點互聯(lián)度高。

全局有效性：全局有效性E定義為

表征節(jié)點間互達的有效性。對于不相連的節(jié)點，Li,j=。全局有效性是最短路徑長度的倒數(shù)和，因此，如果特征路徑長度小，則全局有效性高，意味著網(wǎng)絡是高度互聯(lián)的。從圖8(d)看到，全局有效性與特征路徑長度趨勢相反，燃燒噪聲時值最大，說明節(jié)點高度互聯(lián)。

介數(shù)中心性：單個節(jié)點的介數(shù)中心性定義為經(jīng)過該節(jié)點的最短路徑與所有最短路徑的比值

式中：σi,j為節(jié)點i和j之間最短路徑的數(shù)量；σi,j(v)為節(jié)點i和j之間最短路徑且經(jīng)過節(jié)點v的數(shù)量。介數(shù)中心性確定相空間中連接兩個高密度區(qū)域的低密度區(qū)域。從圖8(e)看到，介數(shù)中心性在陣發(fā)振蕩時值最大，準周期振蕩時較小，燃燒噪聲時值最小。較高的介數(shù)中心性說明吸引子具有較高的局部破碎[20](fragmentation)。陣發(fā)狀態(tài)相空間中高密度區(qū)域的分布不一致，因此介數(shù)中心性在陣發(fā)狀態(tài)時數(shù)值較大。

同配性：定義為所有邊兩端節(jié)點的Pearson相關系數(shù)

式中B=∑i

從圖8(f)看到，同配性在3種狀態(tài)下均為正值，即均為同配網(wǎng)絡。在陣發(fā)振蕩時值最大，準周期振蕩時聚類系數(shù)較小，燃燒噪聲時值最小。

通過上述計算分析，發(fā)現(xiàn)復雜網(wǎng)絡的一些測度，如特征路徑長度、聚類系數(shù)、介數(shù)中心性、網(wǎng)絡直徑、全局有效性等在3種燃燒狀態(tài)下，其值存在明顯的差異，因此均可作為同軸離心式噴嘴燃燒動力學轉變過程的表征指標。進一步，與壓力均方根值一起作為系統(tǒng)穩(wěn)定性的判斷標準，可以提高預警系統(tǒng)的可靠性。

3 結論

采用遞歸網(wǎng)絡方法對同軸離心式噴嘴燃燒穩(wěn)定性進行了分析。獲得了以下結果：

1)試驗中隨著混合比的增加，系統(tǒng)依次經(jīng)歷了燃燒噪聲狀態(tài)、陣發(fā)狀態(tài)和準周期振蕩狀態(tài)。

2)網(wǎng)絡拓撲遞歸網(wǎng)絡體現(xiàn)出了不同狀態(tài)下吸引子的幾何特征。燃燒噪聲狀態(tài)顯示出隨機網(wǎng)絡特征，陣發(fā)狀態(tài)時顯示出中央聚集和周期邊界特點，準周期振蕩狀態(tài)顯示出寬邊界圓環(huán)特征。

3)網(wǎng)絡測度，如聚類系數(shù)、特征路徑長度、介數(shù)中心性、網(wǎng)絡直徑、全局有效性、同配性等均可作為同軸離心式噴嘴燃燒動力學轉變過程的表征指標，可與壓力均方根值一起作為系統(tǒng)穩(wěn)定性的判斷標準。