數值解法在燃燒室出口燃氣溫度計算中的應用

韓 冰，王明瑞，李亞娟，馬 征

（中國航發沈陽發動機研究所，沈陽110015）

0 引言

隨著航空發動機推重比的不斷提升，燃燒室出口溫度越來越高，燃氣分析法作為1種間接測溫方法，在高溫升燃燒室出口溫度場的測試中應用越來越廣泛[1]。燃氣分析測溫法是集熱力學、數值解技術和程序設計于一體的工程測試技術。高溫燃燒產物是多成分的復雜混合物，根據熱力學知識利用燃氣中各成分建立物質守恒、化學平衡、壓力和能量守恒方程構成非線性方程組，對其求解即可得到化學平衡組分和燃氣溫度[2]。

目前非線性方程組的求解沒有通用的解析方法，只能通過數值計算方法獲得，非線性方程組的解法長期以來都是數值計算領域的重要研究內容，已經產生了許多有效求解理論和方法，牛頓-拉夫森法是求解非線性方程組的經典算法[3-4]，突出優點是收斂速度快，缺點是為了得到方程組確切解，需要保證給定的初值接近實際解，否則迭代不易收斂。針對燃燒室各種工作條件，比如燃燒室中不同的油氣比、工作壓力以及燃料組成等相關非線性方程組的計算求解，通用的初值很難確定，只能根據具體算例反復嘗試給定[5-7]。在求解非線性方程組時，如何確定初值和迭代是否收斂以及收斂速度的快慢是燃氣溫度計算所關注的問題。本文采用雙變量迭代法和赫夫法，不僅易確定初值，而且2種算法的迭代收斂速度都很快，通過VB編程實現非線性方程組的求解方法，并將算法成功應用于燃燒室出口燃氣溫度測試中。

1 非線性方程組的建立

燃燒室燃燒產生的高溫反應產物是1種多成分的復雜混合物，對于燃燒室中的穩態燃燒過程，通常采用這樣的假設：由于氣體的停留時間較長，燃氣成分的化學平衡能夠建立，可以認為燃燒室出口燃氣處于化學平衡狀態。而碳氫燃料和空氣的燃燒產物平衡成分一般由 CO2、CO、H2O、Ar、O2、N2、H2、OH、NO、O、H、N共12種成分組成，根據這些成分計算燃燒溫度可滿足計算需求。

通過列出物質守恒、壓力、化學平衡和能量守恒方程構建非線性方程組[8-9]。

1.1 物質守恒方程

設定燃料分子式為CHaObNc，根據1 kg原始反應物與燃燒后生成的1 kg燃氣各元素物質的量相等原則列出碳、氫、氧、氮、氬的守恒方程

式中：As=p·Mg，Mg為燃氣的物質的量；f和 ω 分別為燃燒室油氣比和進口空氣的質量含濕量；RO2、SN2、TCO2和DAr分別為O2、N2、CO2和Ar在空氣中的體積分數；Ma、Mf、Mω分別為空氣、燃料和水蒸氣的物質的量；pi為燃氣中各成分分壓；Ni表示1 kg原始反應物中各元素的物質的量。

1.2 壓力方程

式中：p為燃燒室的工作壓力。

1.3 化學平衡方程

根據化學反應式（4）建立以平衡常數表示的化學平衡方程

式中：Kpi為按氣體分壓計算的化學平衡常數；pˉi為pi與標準狀態壓力p0之比，p0=101325 Pa。

1.4 能量守恒方程

碳氫燃料與空氣燃燒的通用化學反應式為

式中：n0為參與燃燒空氣的物質的量，根據已知油氣比可通過式（7）求出；ng為燃燒生成燃氣的物質的量，在計算出燃氣各組分的物質的量百分數后通過任意1種元素的物質守恒方程求出，如可利用碳元素守恒式（8）求出；hH2O為燃燒室進口空氣中的物質的量含濕量。

式中：a、b與CHaObNc中的下標是同一含義，為燃料中氫和氧的原子數；L0為單位質量燃料完全燃燒所需理論空氣質量；rCO2和rCO為燃氣中CO2和CO的物質的量百分數，其中

根據式（6）建立能量守恒方程

式中：Qf為燃料的低熱值；η為燃氣熱態平衡時的燃燒效率，根據燃氣中未燃組分，包括離解組分可以計算其值；d0為定溫燃燒焓差；Hi,T為各組分在溫度T的焓值；T3為燃燒室進口空氣溫度；T4為燃氣溫度。

式（1）、（3）、（5）、（9）構成非線性方程組，其中未知變量包括12種氣體組分分壓pi、As和燃氣溫度T4共14個。

2 數值解法的應用

求解非線性方程組的根本目的是獲得燃氣溫度，引入的雙變量迭代法和赫夫法主要用于求解化學平衡組分，下面具體介紹2種數值解法在燃氣溫度計算中的應用[10-11]。

2.1 雙變量迭代法

雙變量迭代法的基本方程式可以參考文獻[12]，該算法的基本原理是給出pN2和pO2的初值，通過該算法變換后的方程組計算出其他變量和pN2、pO2新的數值。經過大樣本計算驗證一般取pN2=0.780881p、pO2=0.209495p作為迭代初值，可以保證對燃燒室各種工作條件下的迭代計算都有效收斂。針對貧油燃氣，采用雙變量松弛法同時迭代方案計算更為快捷，迭代公式為

當T＞3000 K或者富油時，采用雙變量同時迭代一般不易收斂，此時采用雙變量分別迭代可滿足計算收斂，即先固定pN2值，對pO2進行迭代，然后再對pN2進行迭代，如此反復迭代直至滿足精度要求為止。雙變量分別迭代法與同時迭代法的方程式形式基本相同，不同的是雙變量分別迭代公式采用雙點弦割法，如式（12）所示。雙變量分別迭代方法編程相比于雙變量同時迭代法要稍微復雜。

通常給定精度小于10-9p的情況下，雙變量迭代法的迭代次數在數次到十幾次。目前針對主燃燒室貧油燃氣的工作條件，雙變量同時迭代方案基本能滿足計算要求，而且其計算收斂時間比雙變量分別迭代要快得多。因此一般在計算主燃燒室出口溫度時采用雙變量同時迭代法基本能滿足計算收斂。

2.2 赫夫法

赫夫法是將牛頓法各方程式和變量取對數形式后求解非線性方程組，這樣處理使得燃氣組分的初值可以任意給定，此外該方法對任何原始反應物和反應產物都是通用的，迭代收斂很好，而且適用于富油燃燒[13-15]。

將式（1）整理成

式中：Δi為各變量的修正量；Ri、fi為過程計算量，無實際意義。

根據式（3）整理成

將式（5）整理成

將式（1）、（3）、（5）以及其中的變量先進行對數處理后，對修正量來說處理后的由式（13）、（16）、（17）構成的方程組是線性的，而對方程組的系數矩陣可以采用高斯主元消去法求解，首先需要給定氣體組分分壓和As的初值，可以選擇按照貧油完全燃燒產物的相應數值給定，這樣可以減少迭代次數，但是這樣確定初值增加了計算量；可以任意給定1組初值，但是取值不能為0，一般在計算中取氮氣的分壓初值為0.8 p，其余氣體和As的分壓初值為10-5p，經過大量抽樣計算驗證，這種初值取值能滿足絕大多數計算工況。經高斯主元消去法得到各修正量并代入式（20）和（21）得到新的各分壓值，再將新的數值重新代入系數矩陣中計算，直至相鄰2次計算結果滿足給定的精度要求為止。在給定精度小于 的情況下，赫夫法的迭代次數在數次到二十幾次，該算法的迭代收斂時間要大于雙變量迭代法的。

2.3 溫度迭代計算

利用雙變量迭代法或赫夫法求解化學平衡成分的計算過程實質是在求解由式（1）、（3）、（5）構成的非線性方程組，計算的前提條件是燃燒室工作壓力、油氣比和燃氣溫度已知，但是由于計算過程中變量燃氣溫度是未知的，因此在開始迭代計算時，需要先假定1個燃氣溫度以便迭代計算。將能量守恒方程式（9）整理為

圖1 燃氣溫度的計算流程

式中：fT為溫度T的函數。

在計算出燃氣平衡成分后，將各成分的物質的量百分數代入式（22），如果|fT|在限定的精度范圍內，則結束計算；如果|fT|的值超出限定值，則采用牛頓法迭代公式（23）計算出新的溫度，再重新計算化學平衡組分，直至|fT|滿足精度要求即可。燃氣溫度的計算流程如圖1所示。

3 算例計算

選用算例對采用雙變量同時迭代法、雙變量分別迭代法和赫夫法計算得到的燃氣平衡組分和燃氣溫度進行比較分析。設定燃料的氫碳比為1.923，低熱值為42650 kJ/kg，進口溫度為500℃，燃燒室工作壓力為0.5 MPa，按貧油和富油2種情況油氣比分別取0.05453和0.07496。采用前面所述的數值解法得到的溫度計算結果見表1、2。

從表1、2中結果可得出：

（1）針對同一工作條件，采用數值解法得到的化學平衡組分和燃氣溫度基本一致；

（2）從計算時間上比較，由快到慢分別為雙變量同時迭代法、雙變量分別迭代法和赫夫法，但是在燃燒室溫度場試驗中，測試軟件采集數據的頻率是1 s采集1組數據，因此任意1種數值解法在測試軟件中的應用都滿足測試需求；

表1 數值解法的結果對比（油氣比0.05453）

表2 數值解法的結果對比（油氣比0.07496）

（3）富油條件下雙變量同時迭代法無法完成迭代計算。

在某5頭部扇形燃燒室出口溫度場試驗中，將熱電偶與燃氣分析的測量結果進行比較。燃燒室采用燃料的氫碳比為1.923，燃料低熱值為42650 kJ/kg，進口溫度為400℃，工作壓力為0.55 MPa，油氣比為0.03。燃燒室中間3個頭部的熱電偶與燃氣分析測得的出口溫度場結果對比如圖2所示，具體數值見表3，其中燃氣分析結果包括前面所述采用數值解法得到的溫度結果。

圖2 油氣比0.03時熱電偶與燃氣分析溫度對比

表3 熱電偶與燃氣分析溫度結果對比

從圖2和表3中可見：

（1）熱電偶所測溫度與燃氣分析計算得到的溫度沿著周向變化的趨勢基本一致，應用雙變量迭代法和赫夫法的燃氣分析結果相同；

（2）燃氣分析的測溫結果偏大于熱電偶的測溫結果，最大相差89℃，最小相差0.4℃，平均相差34.9℃。

4 結論

本文利用雙變量迭代法和赫夫法迭代求解用于計算燃氣溫度的非線性方程組，得出如下結論：

（1）雙變量迭代法計算簡單，收斂速度快，初值易確定并能滿足所有條件下的計算收斂，針對不同的燃燒情況可以采用雙變量同時迭代和分別迭代法完成計算。針對貧油燃氣，采用雙變量同時迭代法收斂更快。

（2）赫夫法是對牛頓-拉夫森法的改進優化，變化后的赫夫法各變量初值幾乎可以任意給定，迭代收斂較好并且該方法適用范圍非常廣，缺點是每次迭代都需要重新計算1次系數矩陣，計算量較大。

（3）通過算例驗證可知，所有數值解法得到的計算結果基本一致，但從計算迭代收斂的速度來看，雙變量同時迭代法計算最快，赫夫法計算相對較慢。

（4）針對某5頭部扇形燃燒室中間3個頭部出口溫度的對比結果可見，燃氣分析計算得到的溫度大于熱電偶測得的溫度。