王大銳，程圣清，張 楠

（北京航天動力研究所，北京100076）

0 引言

近年來，耦合詳細化學(xué)反應(yīng)的液體火箭發(fā)動機燃燒過程研究成為熱點，掌握主要推進劑組合的化學(xué)反應(yīng)機理是開展該研究的基礎(chǔ)。目前，甲基肼（MMH）/四氧化二氮（NTO）,液氫（LH2）/液氧（LO2）和煤油（Kerosene）/液氧（LO2）是3種主要的航天推進劑組合。與后兩種推進劑組合相比，MMH/NTO屬于自燃推進劑，其詳細機理很難獲得。此前，對于此類推進劑組合火箭發(fā)動機燃燒CFD計算通常采用總包反應(yīng)代替，如O.Knab將四步總包反應(yīng)加入到發(fā)動機燃燒模型中[1]，對400 N液體姿軌控火箭發(fā)動機燃燒過程進行了數(shù)值模擬。研究后發(fā)現(xiàn)，總包反應(yīng)對于中間反應(yīng)忽略過大，數(shù)值計算結(jié)果與試驗測量結(jié)果出入較大，因此需要構(gòu)建詳細機理對發(fā)動機燃燒過程進行模擬，其中Laurent構(gòu)建了一套含有82種組分和403步基元反應(yīng)的詳細機理（機理未公布）[2-3]，該機理通過與Seamans的爆炸理論在點火延遲時間上計算結(jié)果的比較[4]，驗證了該機理的正確性，但該機理過于復(fù)雜，并不適合耦合在湍流場中進行發(fā)動機燃燒CFD計算。因此，需要將詳細機理進行適當(dāng)簡化，以滿足發(fā)動機燃燒CFD計算要求。

目前對于復(fù)雜燃燒機理的簡化大體可以分為2大類，一類是以去除影響較小的基元反應(yīng)來進行簡化，如敏感性分析法SA[5-6]，準(zhǔn)穩(wěn)態(tài)假定QSSA[7]，主成分分析法PCA[8-9]等；另一類則是以去除不重要的組分來進行簡化，如直接關(guān)系圖法DRG[10]以及多代通量法PFA[11-12]等。以上各簡化方法各有優(yōu)勢，但也都有局限性。

本文提出將PCA和PFA方法聯(lián)合在一起使用，對43種組分和201步基元反應(yīng)的甲基肼/四氧化二氮反應(yīng)機理進行簡化研究，并將簡化機理與國外文獻以及詳細機理在典型工況下點火延遲時間、平衡溫度以及主要平衡產(chǎn)物計算結(jié)果進行對比。最后，對比不同簡化方法對MMH和NTO，CH4和O2以及H2和O2三種推進劑組合化學(xué)反應(yīng)機理的簡化結(jié)果，說明PP法簡化效果更好，適用范圍更廣。

1 主成分聯(lián)合多代通量機理簡化方法

“第一步簡化”先對含有m種組分和p步基元反應(yīng)機理的速度敏感性系數(shù)矩陣，其形式如下所示：

矩陣F由下式計算得到：

其中

Rj=kjrj(c)

式中：fi為組分i的生成速率；kj為第j步反應(yīng)的反應(yīng)速率常數(shù)；υi,j為基元反應(yīng)j中組分i的化學(xué)計量系數(shù)；Rj為第j步反應(yīng)的生成物的生成速率；rj(c)為第j步反應(yīng)中反應(yīng)物濃度的乘積；Fi,j為第j步反應(yīng)中組分i的生成或者消耗速率與組分i的凈生成率的比值。再將FTF對角化，得到矩陣

式中λ1,λ2,...,λp為FTF的特征值，且λ1>λ2>...>λp。根據(jù)主軸定理，有：

式中：λi為主成分Ψi中各基元反應(yīng)的重要性，即λi大所對應(yīng)的特征向量中，元素越大其代表的基元反應(yīng)越重要。此時，給定簡化機理基元反應(yīng)個數(shù)p1，并要求第一步簡化后機理的著火延遲時間和平衡溫度以及重要組分百分比與原機理計算結(jié)果在一定誤差之內(nèi)，具體表達式如下：

式中：PAB為組分B相對于A的生成相關(guān)系數(shù)；PA和CA分別為組分A的生成通量和消耗通量。類似地，二代和三代相關(guān)系數(shù)表達式如下：

式中Mi1和Mi2為中間組分。

從文獻 [12]上可以知道，理論上代數(shù)越多，簡化結(jié)果越好，但通過對GRI-Mech3.0機理簡化結(jié)果來看，四代計算精度提高并不明顯，所以本文仍采用三代通量計算來進行簡化。最后定義rAB為組分B對組分A的一代、二代及三代相關(guān)系數(shù)相加之和，表達式如下：

判斷組分B對A的重要程度時，將rAB與設(shè)定的閾值ε作比較，如果大于該閾值，則認為該組分對組分A重要，予以保留，否則在簡化機理時應(yīng)予以剔除。當(dāng)全部對A重要的其他組分被找到后，再以這些組分替換A，按照上面的方法繼續(xù)尋找對該組分重要的其他組分，直到最后一次尋找重要組分時沒有找到新的重要組分為止，將此時得到的簡化機理2帶入式 (4)等驗證，驗證后的簡化機理作為PP方法的最終簡化結(jié)果。

2 簡化結(jié)果與分析

本文詳細反應(yīng)機理MN-Detailed包含43種組分201步基元反應(yīng)，動力學(xué)參數(shù)參考文獻 [13]，熱力學(xué)參數(shù)參考Burcat熱力學(xué)庫[14]，燃燒數(shù)值求解由Chemkin軟件完成，簡化機理由自編程序PP-Beta計算得出。

2.1 簡化機理與文獻計算結(jié)果對比

首先對比詳細機理和簡化機理與其他文獻計算結(jié)果。采用Chemkin封閉全混均質(zhì)反應(yīng)器對甲基肼/四氧化二氮自燃著火過程進行模擬，模擬反應(yīng)工況為初始溫度為T=298 K，初始壓力為p=242×104Pa，定容燃燒，得到不同機理著火溫度曲線如下圖1所示。

圖1 不同機理著火過程溫度對比Fig.1 Comparison of temperature in combustion process of different mechanisms

圖1中實線為PP法簡化結(jié)果MN-Reduced機理著火曲線，該機理著火延遲時間為3.33 ms；虛線為MN-Detailed機理著火曲線，著火延遲時間為3.5 ms；點劃線為Catorie利用量子化學(xué)得到的Catorie-Detailed機理[2]，著火延遲時間為3.82 ms；Seamans利用爆炸理論計算結(jié)果為 3.6 ms[4]，MN-Detailed機理和MN-reduced機理著火延遲時間與量子化學(xué)詳細機理和爆炸理論計算的結(jié)果符合較好。

為進一步驗證簡化結(jié)果的正確性，對比了不同壓力、不同混合比、不同初溫下MN-Reduced機理與Seamans根據(jù)爆炸理論著火延遲時間計算的結(jié)果，如圖2所示其中離散點為MN-Reduced機理計算結(jié)果。在溫度和混合比一定情況下，著火時間隨著壓力成負指數(shù)下降，而MN-Reduced機理計算壓力對著火時間影響規(guī)律與Seamans理論計算結(jié)果也符合較好。

圖2 不同壓力下著火延遲時間對比Fig.2 Comparison of ignition delay time at different pressure

2.2 發(fā)動機特定燃燒工況下PP法簡化結(jié)果與原機理計算結(jié)果比較

要將甲基肼/四氧化二氮簡化機理運用到液體火箭發(fā)動機燃燒CFD中，還需要比較簡化機理與原詳細機理在發(fā)動機燃燒典型工況下的計算結(jié)果。液體火箭發(fā)動機燃燒是定壓燃燒過程，所以著火過程采用封閉全混均質(zhì)反應(yīng)器，求解器采用定壓求解。表1中Φ為氧化劑和推進劑質(zhì)量比。

表1 不同壓力下推進劑燃燒初始條件Tab.1 Initial conditions of propellant combination combustion at different pressure

表1和圖3是對初始溫度一定、混合比一定時，不同壓力下簡化機理與詳細機理對著火延遲時間和平衡溫度計算結(jié)果的對比?？梢钥闯?，簡化機理與詳細機理在不同壓力下著火過程均符合較好，著火延遲相差均在3×10-5s以內(nèi)，平衡溫度相差0.1 K以內(nèi)。

圖3 不同壓力下簡化機理和詳細機理著火過程比較Fig.3 Combustion process comparison of simplified mechanism and detailed mechanism at different pressure

表2和圖4是對初始溫度一定、燃燒壓力一定，不同混合比下簡化機理與詳細機理對著火延遲時間和平衡溫度計算結(jié)果的對比?？梢钥闯?，簡化機理與詳細機理在不同混合比下著火過程均符合較好，著火延遲相差均在2×10-5s以內(nèi)，平衡溫度相差0.1 K以內(nèi)。

表2 不同氧燃比下推進劑燃燒初始條件Tab.2 Initial conditions of propellant combination combustion at different mixing ratio

圖4 不同氧燃比下簡化機理與詳細機理燃燒過程比較Fig.4 Combustion process comparison of reduced mechanism and detailed mechanism at different Φ

圖5 不同壓力下簡化機理與詳細機理H2O摩爾百分數(shù)計算結(jié)果對比Fig.5 Comparison of H2O mole fraction calculation results of simplified mechanism and detailed mechanism at different pressure

圖6 不同壓力下簡化機理與詳細機理N2摩爾百分數(shù)計算結(jié)果對比Fig.6 Comparison of N2mole fraction calculation results of simplified mechanism and detailed mechanism at different pressure

圖5和圖6是質(zhì)量混合比1.65，初始溫度298 K，不同壓力下主要平衡產(chǎn)物H2O和N2簡化機理和詳細機理計算結(jié)果的對比?？梢钥闯觯麄€著火過程中簡化機理計算結(jié)果與詳細機理計算結(jié)果均符合較好。

綜上，利用PP法得到的簡化機理在著火延遲時間、平衡溫度以及平衡組分摩爾百分數(shù)都與詳細機理計算結(jié)果保持很好的一致性，這為耦合MMH/NTO化學(xué)反應(yīng)的發(fā)動機燃燒CFD研究奠定了良好的基礎(chǔ)。

2.3 PP法與主成分分析和多代通量方法的比較

主成分法 (PCA)和多代通量分析法 (PFA)都是比較常用的化學(xué)反應(yīng)機理簡化方法，但各自在使用時均有一定的局限性，本文采用PCA，PFA 及 PP法對MMH/NTO，H2/O2[14]及CH4/O2[15]推進劑化學(xué)反應(yīng)詳細機理進行簡化，對比簡化結(jié)果如表3所示。

表3中，MMH/NTO詳細機理在質(zhì)量混合比1.6～1.8，壓力 1.6～2.0 MPa，初始溫度 298 K 范圍內(nèi)簡化；H2/O2詳細機理在摩爾混合比1～4，壓力0.1～1 MPa，著火溫度為1 000 K范圍內(nèi)簡化；CH4/O2詳細機理在摩爾混合比1～2，壓力0.1～1 MPa,著火溫度為1 300～1 500 K范圍內(nèi)簡化。從表3中可看出：一方面PCA簡化機理結(jié)果基元反應(yīng)數(shù)少于PFA簡化結(jié)果，但組分種類要明顯多于PFA，這是因為PCA通過去除非重要基元反應(yīng)達到簡化詳細機理目的，而PFA則是通過去除非重要組分種類來簡化詳細機理，PP法將它們的特點結(jié)合起來，既能夠去除非重要基元反應(yīng)又能夠去除非重要組分，因而得到更為簡單的簡化機理；另一方面，從H2/O2和CH4/O2簡化結(jié)果來看，PCA、PFA存在一定局限性，只能對一部分推進劑反應(yīng)機理進行簡化，而PP法可以簡化機理的范圍比前兩者都大，具有發(fā)展成一種通用簡化方法的潛力。

表3 不同簡化方法對3種推進劑組合簡化結(jié)果對比Tab.3 Comparison of simplified results of three combined propellants using different simplification methods

3 結(jié)論

本文運用PP方法對MMH/NTO推進劑組合43種組分201步基元詳細機理進行了簡化，簡化機理與國外文獻以及詳細機理在多種典型工況下著火延遲時間、平衡溫度以及平衡組分摩爾百分數(shù)計算結(jié)果都符合較好，說明簡化機理可以有效代替詳細機理進行燃燒計算，并為進行耦合化學(xué)反應(yīng)的液體姿軌控火箭發(fā)動機燃燒CFD研究奠定了良好基礎(chǔ)；PP方法相比于PCA和PFA方法具備可同時減少基元反應(yīng)個數(shù)和組分種類的特點，簡化效果好，簡化機理適用性也更廣。

[1]KNAB O,PRECLIK D,ESTUBLIER D.Flow field prediction within liquid film cooled combustion chambers of storable bipropellant rocket engines,AIAA 98-3370[R].USA：AIAA,1998.

[2]CATOIRE L,CHAUMEIX N,PAILLARD C.Chemical kinetic model for monomethylhydrazine/nitrogen tetroxide gas-phase combustion and hypergolic ignition[J].Journal of Propulsion and Power,2004,20(1)：87-92.

[3]CATOIRE L,SWIHART M T.Thermochemistry of species produced from monomethylhydrazine in propulsion andspace-relatedapplications[J].Journal ofPropulsion and Power,2002,18(6)：1242-1252.

[4]SEAMANS T F,VANPEE M,AGOSTA V D.Development of a fundamental Model of hypergolic ignition in space-ambient engines[J].AIAA Journal,1967,l5(9)：1616-1624.

[5]TURANYI T,BERCES T,VAIDA S.Reaction rate analysis of complex kinetic system[J].International Journal of Chemical Kinetics,1989,21(2)：83-99.

[6]喬瑜,徐明厚.基于敏感性分析的甲烷反應(yīng)機理優(yōu)化簡化[J].華中科技大學(xué)學(xué)報(自然科學(xué)版),2007,35(5)：85-87.

[7]PETERS N,WILLIAMS F A.The asymptotic structure of stoichiometric methane-airflame[J].Combustion and flame,1987,68(2)：185-207.

[8]LU T F,CHUNG K L.A directed relation graph method for mechanism reduction[J].Proceedings of the Combustion Institute,2005,30：1333-1341.

[9]VAJDA S,VALKO P,TURANYI T.Principle component analysis of kinetic model[J].International Journal of Chem-Kinet,1985,17：55-81.

[10]董剛,邱榕,蔣勇,張和平.利用主成分分析法簡化甲烷/空氣層流預(yù)混火焰的反應(yīng)機理[J].火災(zāi)科學(xué),2004,13(3)：158-162.

[11]SUN W T,CHEN Z,GOU X L,et al.A path flux analysis method for the reduction of detailed chemical kinetic mechanisms[J].Combustion and Flame,2010,157(7)：1298-1307.

[12]茍小龍,王衛(wèi),桂瑩,施萬玲.一種多代路徑通量分析化學(xué)機理簡化方法[J].推進技術(shù),2012,33(3)：412-417.

[13]ANDERSON W R.,MCQUAID M J,MICHAEL J,et al.A detailed,finite-rate,chemical kinetics mechanism for monomethylhydrazine-red fuming nitric acid systems,ARL-TR88-502010[R].[S.l.]：ARL-TR,1988.

[14]MATSUO A,FUJII K,FUJIWARA T.Flow features of shock-induced combustion around projectile traveling at hypervelocities[J].AIAA Journal,1995,33(6)：1056-1063.

[15]YANG B,POPE S B.An investigation of the accuracy of manifold methods and splitting schemes in the computational implementation of combustion chemistry[J].Combustion and Flame,1998,112：16-32.

[16]聶萬勝,豐松江.液體姿軌控火箭發(fā)動機燃燒動力學(xué)模型與數(shù)值模擬[M].北京：國防工業(yè)出版社,2011.