Rijke 管熱聲振蕩聲學(xué)特性數(shù)值模擬

張澤昊， 何 博， 聶萬勝， 董立寶

（航天工程大學(xué)宇航科學(xué)與技術(shù)系， 北京 101416）

0 引言

熱聲振蕩是指燃燒裝置中， 不穩(wěn)定熱釋放與聲學(xué)振蕩相耦合而產(chǎn)生的一種現(xiàn)象， 通常表現(xiàn)為高振幅的壓力振蕩，一直是燃燒裝置研發(fā)過程中的難題[1-3]。在研究該問題過程中，Rijke 管是最典型的系統(tǒng)， 因此國內(nèi)外學(xué)者針對(duì)該熱聲系統(tǒng)開展了大量研究。

如今， 通過數(shù)值模擬對(duì)熱聲振蕩開展研究的主要方法包括計(jì)算流體力學(xué)（CFD）、熱聲網(wǎng)絡(luò)法以及計(jì)算氣動(dòng)聲學(xué)（CAA）[4-10]。 CFD 通過求解完整的N-S 方程，可以確定Rijke 管中熱聲振蕩的發(fā)生，但計(jì)算成本高，同時(shí)難以分析熱聲振蕩產(chǎn)生的機(jī)理[11]。熱聲網(wǎng)絡(luò)法將燃燒室劃分為具有平均流動(dòng)特性的一維連續(xù)子單元，每個(gè)子單元之間通過傳遞矩陣相連接，可以有效減少計(jì)算量，但不適用于復(fù)雜燃燒室的分析[7]。 CAA 可以彌補(bǔ)以上兩種方法的不足。

目前，通過CAA 對(duì)Rijke 管開展熱聲振蕩研究，會(huì)忽略管中平均流速和熱源處的溫度梯度， 且很少對(duì)橫向模態(tài)開展研究[2，8-10]。 因此，本文使用穩(wěn)定有限元法在頻域求解帶源項(xiàng)的線性歐拉方程（LEE），在考慮平均流速及火焰處溫度梯度后，研究了火焰熱釋放、平均流馬赫數(shù)、燃燒響應(yīng)模型對(duì)Rijke 管聲學(xué)特性的影響。

1 數(shù)學(xué)模型

1.1 線性歐拉方程

從N-S 方程出發(fā)，假設(shè)氣體為理想氣體，忽略熱傳導(dǎo)及粘性，得到下面的歐拉方程及狀態(tài)方程[12]

式中，ζ—平均流場的瞬時(shí)變量， 可替換為ρ、u、p、q，‘—’—變量的平均量，‘′’—聲場波動(dòng)量。假設(shè)每個(gè)波動(dòng)量與平均量相比無限小其中，ωr—角頻率；ωi—增長率。 當(dāng)大于0 時(shí)，系統(tǒng)表現(xiàn)為不穩(wěn)定。 將式（12）～（15）帶入方程式（9）～（11）后得到頻域中的線性歐拉方程

在數(shù)值求解過程中， 采用穩(wěn)定有限元法將控制方程離散為特征值方程，穩(wěn)定技術(shù)選取伽遼金/最小二乘法，采用ARPACK 算法對(duì)離散后的特征值方程進(jìn)行求解。 離散及求解線性歐拉方程的詳細(xì)過程見文獻(xiàn)[13]。

1.2 燃燒響應(yīng)模型

假設(shè)熱釋放波動(dòng)發(fā)生在一固定薄層處， 火焰面上的熱釋放波動(dòng)主要受到未燃?xì)怏w控制， 采用Dowling 和Campa 使用的速度時(shí)滯模型[11，14-15]

2 Rijke 管熱聲振蕩特性影響因素分析

2.1 Rijke 管平均流場參數(shù)和邊界條件

圖1 給出了Rijke 管的二維幾何構(gòu)型， 是一長度為、直徑為的等截面管道，xf=0.5m 處有一厚度為δf的火焰區(qū)域，火焰左、右分別為未燃?xì)怏w、燃燒產(chǎn)物。

圖1 Rijke 管結(jié)構(gòu)示意圖

表1 給出了Rijke 管未燃?xì)怏w的平均流場參數(shù)值。 假設(shè)未燃?xì)怏w、燃燒產(chǎn)物在管道中的流場參數(shù)保持均勻，未燃?xì)怏w通過火焰瞬間轉(zhuǎn)變?yōu)槿紵a(chǎn)物，且不考慮兩者之間物質(zhì)組分變化。若此時(shí)未燃?xì)怏w溫度T1上升為燃燒產(chǎn)物溫度T2，則根據(jù)質(zhì)量、動(dòng)量守恒和理想氣體狀態(tài)方程，可以得到火焰下游燃燒產(chǎn)物速度ux2、壓力p2、密度ρ2與火焰上游未燃?xì)怏w速度ux1、 壓力p1、密度ρ1之間的關(guān)系：

表1 Rijke 管幾何模型及火焰上游平均流場參數(shù)

為了使熵波僅在火焰處產(chǎn)生，Rijke 管入口及出口聲邊界條件被設(shè)置為

該邊界條件表示：入口無速度及熵波動(dòng)，出口無壓力波動(dòng)。 其余邊界為硬壁面，網(wǎng)格采用非結(jié)構(gòu)化網(wǎng)格。

2.2 結(jié)果分析

2.2.1 火焰區(qū)域熱釋放的影響

火焰區(qū)域熱釋放的變化會(huì)使火焰上、 下游溫度產(chǎn)生變化。 假設(shè)此時(shí)上游溫度T1不變，通過改變火焰下游區(qū)域溫度T2與上游區(qū)域溫度T1比值K，研究火焰區(qū)域熱釋放對(duì)Rijke 管熱聲振蕩聲學(xué)特性的影響，K 越大， 表明火焰區(qū)熱釋放越強(qiáng)。

圖2 給出了入口馬赫數(shù)Ma1=0.02 時(shí)，溫度比K 對(duì)Rijke 管前四階縱向模態(tài)本征頻率及增長率的影響。

圖2 Ma=0.02 時(shí)前四階本征模態(tài)增長率隨溫度比的變化，箭頭表示溫度比K 由1 至3.5 的變化，ΔK=0.5，‘▲、▼、●、◆’分別為一階、二階、三階、四階模態(tài)

可以看出， 當(dāng)溫度比K=1 時(shí)，Rijke 管前四階縱向模態(tài)的增長率為0，系統(tǒng)處于邊緣穩(wěn)定狀態(tài)；隨著K 的增加，同一階模態(tài)的本征頻率逐漸增大，而其增長率逐漸減小。Rijke 管本征頻率的增大是因?yàn)闇囟缺菿 增加，導(dǎo)致火焰區(qū)域下游溫度升高，下游區(qū)域的聲速也相應(yīng)提高。而其增長率的變化與熵波的產(chǎn)生和傳播有關(guān)，火焰區(qū)域的溫度梯度使熵波產(chǎn)生，并通過平均流進(jìn)行傳播。當(dāng)溫度比K=1時(shí)，火焰區(qū)域熱釋放為0， 無熵波產(chǎn)生， Rijke 管系統(tǒng)不會(huì)產(chǎn)生阻尼；當(dāng)溫度比K＞1 時(shí)，火焰區(qū)域產(chǎn)生溫度梯度，熵波產(chǎn)生，此時(shí)平均流運(yùn)輸熵波，系統(tǒng)產(chǎn)生阻尼，穩(wěn)定性增強(qiáng)。

Rijke 管的一階橫向模態(tài)本征頻率始終穩(wěn)定于4056Hz，且增長率為正，大小幾乎不變。 說明對(duì)流熵波的中熵波幅值的增加對(duì)該模態(tài)產(chǎn)生的影響較小。

楊亞晶[16]在Rijke 管熱聲振蕩實(shí)驗(yàn)中， 得到Rijke 管的本征頻率隨熱源功率增加而逐漸提高的結(jié)論，本節(jié)的結(jié)果與該結(jié)論相一致，驗(yàn)證了該方法的準(zhǔn)確性。 但實(shí)驗(yàn)中只能獲得熱聲振蕩發(fā)生時(shí)的頻率，不能對(duì)穩(wěn)定性裕度進(jìn)行評(píng)估，使用該方法可以有效解決這一問題。

2.2.2 平均流馬赫數(shù)的影響

本節(jié)選取溫度比K=3 反映典型發(fā)動(dòng)機(jī)的工況， 通過改變?nèi)肟隈R赫數(shù)Ma1研究平均流馬赫數(shù)對(duì)Rijke 管聲學(xué)特性的影響。

圖3、 圖4 分別給出了Rijke 管前四階縱向模態(tài)、一階橫向模態(tài)本征頻率和增長率隨入口馬赫數(shù)Ma1增長的變化規(guī)律。其中縱向模態(tài)的變化規(guī)律與Nicoud 等[2，11]人相一致，但他們沒有研究橫向模態(tài)所受影響。

由圖3 可以看出，當(dāng)入口馬赫數(shù)Ma1為0 時(shí)，增長率為0；隨著入口馬赫數(shù)Ma1的增大，系統(tǒng)前四階模態(tài)的本征頻率逐漸減小，增長率也逐漸減小， 系統(tǒng)的穩(wěn)定性增強(qiáng)；另外，馬赫數(shù)對(duì)增長率的影響更為顯著。 這是因?yàn)楫?dāng)入口馬赫數(shù)Ma1=0 時(shí)，Rijke 管中無平均流動(dòng)， 火焰區(qū)域雖有熵波存在，但是不會(huì)傳播到其他區(qū)域，無阻尼產(chǎn)生；當(dāng)入口馬赫數(shù)Ma1＞0 時(shí)，平均流會(huì)運(yùn)輸火焰區(qū)域產(chǎn)生的熵波，同時(shí)所設(shè)置的邊界條件允許聲能輻射到計(jì)算域外， 系統(tǒng)中聲能衰減，處于穩(wěn)定狀態(tài)；馬赫數(shù)越大，穩(wěn)定性越強(qiáng)。

圖3 K=3 時(shí)前四階縱向模態(tài)增長率隨入口馬赫數(shù)的變化，箭頭表示入口馬赫數(shù)Ma1 由0 至0.16的變化，ΔMa1=0.02，‘▲、▼、●、◆’分別為一階、二階、三階、四階模態(tài)

由圖4 可以看出，隨著入口馬赫數(shù)Ma1的增加，一階橫向模態(tài)的本征頻率不斷減小，而增長率呈先增加，后減小的趨勢(shì)，且?guī)缀醵即笥?。這說明平均流速的提高可能使Rijke管發(fā)生橫向振蕩。

圖4 一階橫向模態(tài)受入口馬赫數(shù)的影響

2.2.3 火焰區(qū)域熱釋放波動(dòng)的影響

本節(jié)設(shè)置溫度比K=3，入口馬赫數(shù)Ma1=0.02，其他流場條件不變。 在火焰區(qū)域添加燃燒響應(yīng)模型，通過改變模型中的參數(shù)n和τ，來改變聲與火焰之間的耦合情況，預(yù)測(cè)熱聲振蕩現(xiàn)象。

（1）相互作用指數(shù)影響。 圖5 和圖6 分別給出了時(shí)滯τ=0 時(shí)，Rijke 管前四階縱向模態(tài)本征頻率及增長率隨相互作用指數(shù)n 的變化。 由圖5 可以看出， 隨著n 的增加，Rijke 管一階、四階縱向模態(tài)的本征頻率逐漸增大，二階、三階縱向模態(tài)的本征頻率逐漸減小。 由圖6 可以看出，隨著n 的增加，一階及三階縱向模態(tài)的增長率逐漸增大，系統(tǒng)穩(wěn)定性降低，但二階縱向模態(tài)的增長率逐漸減小，四階模態(tài)的增長率幾乎不變。

圖5 Ma=0.02， 時(shí)滯τ=0 時(shí)前四階模態(tài)本征頻率隨相互作用指數(shù)的變化

圖6 Ma=0.02，時(shí)滯τ=0 時(shí)前四階模態(tài)增長率隨相互作用指數(shù)的變化

（2）時(shí)滯的影響。時(shí)滯τ 反映了熱釋放波動(dòng)q′ 與速度波動(dòng)u′之間的相位差，該相位差可表示為

當(dāng)t=1 時(shí)，時(shí)滯變?yōu)橐越穷l率ω 振蕩的一個(gè)周期，此時(shí)熱釋放波動(dòng)與速度波動(dòng)之間的相位差為2π。

圖7（a）～（d）分別給出了n=0.5 時(shí)，Rijke 管前四階縱向模態(tài)隨歸一化時(shí)滯的變化。 可以發(fā)現(xiàn)，不同模態(tài)的本征頻率及增長率隨著的變化呈周期性變化。 當(dāng)0＜t＜0.5 時(shí)，一階及四階模態(tài)本征頻率降低，增長率大于0，系統(tǒng)處于不穩(wěn)定狀態(tài)，二階及三階模態(tài)則相反；當(dāng)0.5＜t＜1 時(shí)，一階、四階模態(tài)進(jìn)入穩(wěn)定狀態(tài)。 同時(shí)， 這些變化存在一個(gè)共同點(diǎn)，當(dāng)本征頻率增大時(shí)，系統(tǒng)處于穩(wěn)定狀態(tài)，當(dāng)本征頻率減小時(shí)，系統(tǒng)處于不穩(wěn)定狀態(tài)。

圖7 Ma=0.02，n=0.5 時(shí)前四階模態(tài)隨歸一化時(shí)滯的變化

系統(tǒng)穩(wěn)定狀態(tài)與不穩(wěn)定狀態(tài)之間的切換， 在各類時(shí)滯系統(tǒng)中廣泛存在。由Rayleigh 準(zhǔn)則可知，當(dāng)熱釋放波動(dòng)q′ 與壓力波動(dòng)p′ 處于同相位， 即兩者相位差為0°～90°[17]時(shí)，熱聲振蕩發(fā)生。在燃燒響應(yīng)模型中，改變，熱釋放波動(dòng)與速度波動(dòng)之間的相位差發(fā)生改變， 而速度波動(dòng)與壓力波動(dòng)的相位相差90°[14]。 由此可知，如果速度波動(dòng)滯后壓力波動(dòng)90°， 那么當(dāng)熱釋放波動(dòng)相位超前速度波動(dòng)0°～180°時(shí)，熱聲振蕩發(fā)生，對(duì)應(yīng)圖7 中二階、三階縱向模態(tài)t＞0.5 時(shí)；如果壓力波動(dòng)滯后速度波動(dòng)90°，那么當(dāng)熱釋放波動(dòng)相位滯后速度波動(dòng)0°～180°時(shí)，熱聲振蕩發(fā)生，對(duì)應(yīng)圖7 中一階、四階縱向模態(tài)t＜0.5 時(shí)。

Rijke 管一階橫向模態(tài)對(duì)火焰區(qū)域的熱釋放波動(dòng)不敏感，其本征頻率始終穩(wěn)定在4056Hz，增長率變化也很小。

3 結(jié)論

本文以Rijke 管熱聲系統(tǒng)為研究對(duì)象，在不同情況下對(duì)Rijke 管熱聲系統(tǒng)進(jìn)行線性穩(wěn)定性分析，得出以下結(jié)論：

采用穩(wěn)定有限元法求解頻域線性歐拉方程，可快速獲得Rijke 管的聲學(xué)特性，適合燃燒裝置的熱聲振蕩數(shù)值模擬。

火焰熱釋放對(duì)Rijke 管縱向模態(tài)的本征頻率影響顯著，本征頻率隨火焰熱釋放的提高而增大；平均流馬赫數(shù)更易對(duì)增長率產(chǎn)生影響， 其隨平均流馬赫數(shù)的增大而下降； 縱向模態(tài)的本征頻率及增長率隨燃燒響應(yīng)模型中時(shí)滯的變化呈周期性變化。

Rijke 管一階橫向模態(tài)受平均流馬赫數(shù)影響明顯。