采用態(tài)－態(tài)模型的熱化學(xué)非平衡噴管流數(shù)值研究

徐 丹， 曾 明， 張 威， 柳 軍

（國(guó)防科學(xué)技術(shù)大學(xué)航天科學(xué)與工程學(xué)院，湖南長(zhǎng)沙 410073）

采用態(tài)－態(tài)模型的熱化學(xué)非平衡噴管流數(shù)值研究

徐 丹， 曾 明， 張 威， 柳 軍

（國(guó)防科學(xué)技術(shù)大學(xué)航天科學(xué)與工程學(xué)院，湖南長(zhǎng)沙 410073）

數(shù)值求解耦合態(tài)－態(tài)模型的N2／N準(zhǔn)一維噴管非平衡流方程，獲得駐室溫度3000 K～10000 K，壓力1 atm～10 atm條件下的流場(chǎng)參數(shù)分布和詳細(xì)的振動(dòng)能級(jí)分布，分析噴管流動(dòng)中的宏觀和微觀特性，并對(duì)這類復(fù)合占優(yōu)的流動(dòng)考察雙溫度模型中假設(shè)的合理性.

態(tài)－態(tài)模型；熱化學(xué)非平衡；振動(dòng)能級(jí)分布；高超聲速噴管流動(dòng)；數(shù)值模擬

0 引言

高超聲速高焓流動(dòng)中常同時(shí)存在熱力和化學(xué)非平衡現(xiàn)象.對(duì)熱力非平衡流動(dòng)，氣體各種模式的內(nèi)能無(wú)法用統(tǒng)一的溫度表征.一般而言，轉(zhuǎn)動(dòng)與平動(dòng)內(nèi)能模式較易達(dá)到平衡，可以用同一溫度描述，記作平動(dòng)／轉(zhuǎn)動(dòng)溫度或平動(dòng)溫度T，且平動(dòng)和轉(zhuǎn)動(dòng)能級(jí)分布為該溫度下的Boltzmann分布.在電子能可以忽略的情況下，熱力非平衡主要表現(xiàn)為振動(dòng)非平衡.有兩類描述流動(dòng)中熱力非平衡現(xiàn)象和熱力－化學(xué)耦合的方法.一類是多溫度模型法，以雙溫度模型法［1］為代表，通過(guò)振動(dòng)能變化率方程描述振動(dòng)非平衡過(guò)程，以不同于平動(dòng)溫度的振動(dòng)溫度Tv表征非平衡振動(dòng)能，將化學(xué)反應(yīng)控制溫度設(shè)為平動(dòng)溫度與振動(dòng)溫度的某種平均值以反映熱力－化學(xué)耦合；另一類是態(tài)－態(tài)模型法［2－3］，直接描述分子在各振動(dòng)能級(jí)間的非平衡躍遷過(guò)程，給出各能級(jí)粒子數(shù)變化率方程，熱力－化學(xué)耦合則通過(guò)分別描述各能級(jí)上分子的化學(xué)反應(yīng)過(guò)程自然體現(xiàn).

雙溫度模型法引入振動(dòng)能變化率方程和流動(dòng)控制方程耦合求解，得到流場(chǎng)各點(diǎn)的非平衡振動(dòng)能后借用平衡振動(dòng)能和溫度的關(guān)系式，折算出振動(dòng)溫度Tv，它與平動(dòng)溫度的差別能直觀地反映流動(dòng)中的振動(dòng)非平衡現(xiàn)象.振動(dòng)溫度進(jìn)入化學(xué)反應(yīng)控制溫度也是一個(gè)簡(jiǎn)便地體現(xiàn)熱力－化學(xué)耦合的方法.雙溫度模型及其后來(lái)拓展的多溫度模型一定程度上反映了流動(dòng)中的熱力非平衡現(xiàn)象和熱力－化學(xué)耦合，在高溫非平衡流的數(shù)值模擬中獲得了廣泛的應(yīng)用，在許多條件下也獲得了良好的數(shù)值結(jié)果，但是也存在局限［4］，例如在復(fù)合區(qū)的效果不如離解區(qū)，無(wú)法合理解釋AVCO實(shí)驗(yàn)數(shù)據(jù)，等.其實(shí)化學(xué)反應(yīng)控制溫度（即化學(xué)反應(yīng)速率系數(shù)Arrhenius公式中的溫度）中包含振動(dòng)溫度的內(nèi)涵是值得進(jìn)一步分析的.Arrhenius速率系數(shù)公式中的溫度參數(shù)一方面是反映粒子碰撞頻率，與隨機(jī)熱運(yùn)動(dòng)（平動(dòng)）速度有關(guān)；另一方面則是反映碰撞能夠引起化學(xué)反應(yīng)的概率，與參與碰撞粒子的內(nèi)能（包括平動(dòng)、轉(zhuǎn)動(dòng)、振動(dòng)、電子能）有關(guān)，這里以溫度反映粒子內(nèi)能的物理內(nèi)涵其實(shí)是基于參與反應(yīng)的粒子能級(jí)分布滿足該溫度下的平衡分布（Boltzmann分布）的.因此從這個(gè)角度看，在化學(xué)反應(yīng)控制溫度中引入振動(dòng)溫度，可以說(shuō)隱含著對(duì)振動(dòng)非平衡流采用了振動(dòng)能級(jí)分布滿足振動(dòng)溫度下Boltzmann分布的假設(shè)，這個(gè)假設(shè)是否合理，需要進(jìn)一步研究.另外，采用何種方式對(duì)平動(dòng)溫度和振動(dòng)溫度取平均以確定化學(xué)反應(yīng)控制溫度主要依據(jù)經(jīng)驗(yàn)，也具有相當(dāng)?shù)牟淮_定性.

態(tài)－態(tài)模型法描述振動(dòng)非平衡及其與化學(xué)反應(yīng)耦合的方法則更為精細(xì)和直接.它將粒子在能級(jí)間的躍遷與化學(xué)反應(yīng)視為同一性質(zhì)的過(guò)程研究，能夠得到流體組元詳細(xì)能級(jí)分布而不只是總的非平衡振動(dòng)能的變化規(guī)律，掌握了更多的細(xì)節(jié)；它直接研究各能級(jí)上粒子的化學(xué)反應(yīng)，避免了多溫度模型中經(jīng)驗(yàn)性的熱力－化學(xué)耦合方法.不過(guò)，盡管態(tài)－態(tài)模型在理論上展現(xiàn)出諸多優(yōu)點(diǎn)，但各種內(nèi)能模式能級(jí)間的躍遷過(guò)程復(fù)雜，使得計(jì)算中處理的信息量巨大，并且計(jì)算結(jié)果強(qiáng)烈依賴于所使用的微觀模型和動(dòng)力學(xué)給出的相關(guān)數(shù)據(jù)，所以目前還無(wú)法實(shí)現(xiàn)態(tài)－態(tài)模型與多維流場(chǎng)程序的耦合［5］，國(guó)外采用態(tài)－態(tài)模型的研究一般集中在不涉及流動(dòng)的零維問(wèn)題［5］，或是正激波后、邊界層［6－7］和準(zhǔn)一維噴管非平衡流動(dòng)［8－9］問(wèn)題.

國(guó)內(nèi)在態(tài)－態(tài)模型的研究方面剛剛起步，一些學(xué)者對(duì)模型基本原理和實(shí)驗(yàn)方法進(jìn)行了總結(jié)和討論［10］，筆者曾采用態(tài)－態(tài)模型研究零維問(wèn)題［11］，對(duì)等溫等壓或等溫等容的N2／N混合物系統(tǒng)，詳細(xì)研究其熱化學(xué)非平衡過(guò)程，揭示了非平衡過(guò)程中的一些微觀分布信息，加深了對(duì)態(tài)－態(tài)模型和非平衡過(guò)程的理解.

進(jìn)一步將態(tài)－態(tài)模型與流動(dòng)耦合，研究N2／N混合物的準(zhǔn)一維熱化學(xué)非平衡噴管流動(dòng).國(guó)外學(xué)者求解耦合態(tài)－態(tài)模型的準(zhǔn)一維非平衡噴管流控制方程時(shí)是采用空間推進(jìn)方法［8－9］，考慮到非平衡噴管流中喉部非聲速等原因造成的入口條件試湊的困難，本文采用時(shí)間推進(jìn)法.在駐室溫度3 000 K～10 000 K，壓力1 atm～10 atm條件下數(shù)值求解N2／N非平衡噴管流動(dòng)，獲得流場(chǎng)參數(shù)分布特別是詳細(xì)的振動(dòng)能級(jí)分布，分析噴管流動(dòng)中的宏觀和微觀特性變化特點(diǎn).對(duì)這類復(fù)合占優(yōu)的流動(dòng)，考察雙溫度模型中假設(shè)的合理性.

1 基本方程與計(jì)算方法

態(tài)－態(tài)模型將位于不同振動(dòng)能級(jí)上的分子視為不同組元，位于振動(dòng)能級(jí)v上的N2記作N2（v）.本文采用有68（0～67）個(gè)振動(dòng)能級(jí)的N2能級(jí)模型，因此N2（v）／N混合物中含68個(gè)分子組元和1個(gè)N原子組元，共69個(gè)組元.態(tài)－態(tài)模型直接描述分子在不同振動(dòng)能級(jí)間的躍遷和各能級(jí)上分子N2（v）與原子N的化學(xué)反應(yīng)，共5類過(guò)程：分子間的振動(dòng)－平動(dòng)能量交換過(guò)程（vTm）、分子與原子間的振動(dòng)－平動(dòng)能量交換過(guò)程（vTa）、振動(dòng)－振動(dòng)能量交換過(guò)程（vv）；原子作為碰撞參與者的離解－復(fù)合反應(yīng)（dra）和分子作為碰撞參與者的離解－復(fù)合反應(yīng)（drm）.具體表示如下［12］

其中v和w代表N2的不同振動(dòng)能級(jí).根據(jù)上述各類躍遷或化學(xué)反應(yīng)的速率方程就可給出N2（v）和N共計(jì)69個(gè)組元的質(zhì)量生成率.各類過(guò)程的速率系數(shù)具體計(jì)算公式和相關(guān)數(shù)據(jù)見文獻(xiàn)［6，12－13］.

對(duì)非平衡流動(dòng)，需要69個(gè)組元連續(xù)方程.加上動(dòng)量、能量方程后，守恒形式的準(zhǔn)一維熱化學(xué)非平衡噴管流控制方程為

其中ρi為N2（v）和N各組元的密度，A為噴管截面面積，E為單位體積的總能，˙wi為的各組元的質(zhì)量生成率，由上述（1）～（5）式給出的各類過(guò)程的速率方程表達(dá)，其中包含各組元密度ρi．

對(duì)控制方程（6）式進(jìn)行差分離散時(shí)，無(wú)粘通量項(xiàng)E和源項(xiàng)W采用全隱格式，右端項(xiàng)F采用顯式離散，隱式求解方法為L(zhǎng)U-SGS方法，具體參見文獻(xiàn)［14］.源項(xiàng)隱式求解時(shí)需要的（?W／?U）矩陣需根據(jù)態(tài)－態(tài)模型各躍遷過(guò)程的速率系數(shù)公式推導(dǎo)計(jì)算［15］.

定解條件中，噴管入口給出總溫總壓條件，流速由內(nèi)點(diǎn)外插得到，假設(shè)駐室至入口第一點(diǎn)流動(dòng)接近平衡因而等熵.入口點(diǎn)的N原子和各振動(dòng)能級(jí)N2（v）的質(zhì)量分?jǐn)?shù)采用駐室溫度壓力下的態(tài)－態(tài)模型平衡結(jié)果給定.噴管超聲速出口邊界上所有物理量由內(nèi)點(diǎn)外插得到.初場(chǎng)的確定首先根據(jù)噴管各點(diǎn)膨脹比分布馬赫數(shù)，以熱化學(xué)平衡條件給出對(duì)應(yīng)馬赫數(shù)的宏觀流動(dòng)參數(shù)，N和N2（v）組元的質(zhì)量分?jǐn)?shù)則取為駐室條件下的態(tài)－態(tài)模型的平衡結(jié)果.

為驗(yàn)證筆者編制的“態(tài)－態(tài)模型”數(shù)值計(jì)算程序，曾對(duì)文獻(xiàn)［12］的算例，進(jìn)行了N2／N系統(tǒng)熱化學(xué)非平衡過(guò)程的數(shù)值模擬，程序結(jié)果和文獻(xiàn)結(jié)果吻合很好［11］.

2 數(shù)值算例與結(jié)果分析

噴管流動(dòng)計(jì)算針對(duì)一軸對(duì)稱的噴管進(jìn)行.噴管半徑沿軸線的分布曲線為拋物線

噴管長(zhǎng)度為1 m，喉部位置在x＝0.5 m處.

噴管駐室溫度T0在3 000 K～10 000 K范圍，駐室壓力p0在1 atm～10 atm范圍，本文選取了多個(gè)算例條件求解N2／N非平衡噴管流場(chǎng)，下面以T0＝8 000 K（N2離解程度較高）和T0＝3 000 K（離解程度很低）兩組算例為代表給出數(shù)值計(jì)算結(jié)果和分析.

2.1 T0＝8 000 K算例

采用了態(tài)－態(tài)模型進(jìn)行噴管非平衡流計(jì)算得到的宏觀流動(dòng)參數(shù)如壓力、溫度、密度、速度在氣流膨脹過(guò)程中的變化特點(diǎn)與雙溫度模型得到的結(jié)果一致，這里不再給出，下面主要討論組元質(zhì)量分?jǐn)?shù)或N2（v）與N2的粒子數(shù)之比沿軸線的分布，噴管中不同位置的振動(dòng)能級(jí)分布及其偏離Boltzmann分布的特點(diǎn).

2.1.1 參數(shù)沿軸線分布

圖1給出了N質(zhì)量分?jǐn)?shù)CN沿噴管軸線的分布.p0＝1 atm條件下CN約在0.7～0.8之間，5 atm時(shí)約為0.4～0.5，10 atm時(shí)約為0.3～0.4.這是總壓升高抑制駐室N2離解，而膨脹過(guò)程中提高壓力又促進(jìn)N復(fù)合生成N2的緣故.三個(gè)總壓條件下噴管喉道下游都表現(xiàn)出化學(xué)凍結(jié)的現(xiàn)象，但總壓越大，凍結(jié)出現(xiàn)得越晚.

圖2以第0（基態(tài)）和第30振動(dòng)能級(jí)為代表給出了N2（v）質(zhì)量分?jǐn)?shù)CN2（v）沿軸線的分布.各振動(dòng)能級(jí)N2（v）的質(zhì)量分?jǐn)?shù)在入口到x／L＝0.4這段范圍內(nèi)變化都較小，直到喉道附近才發(fā)生明顯變化，最后在喉道下游的不同位置趨于凍結(jié).喉道前溫度變化較小，各能級(jí)的CN2（v）也相對(duì)穩(wěn)定，到喉道附近后溫度快速降低，CN2（v）發(fā)生明顯變化.喉道下游氣流膨脹過(guò)程中溫度繼續(xù)下降，但由于壓力和密度的快速下降，各類躍遷和化學(xué)反應(yīng)速率都急劇下降，而流速又快速增大，因此喉道下游不遠(yuǎn)處就表現(xiàn)出流動(dòng)熱化學(xué)凍結(jié)的特點(diǎn).不同振動(dòng)能級(jí)出現(xiàn)凍結(jié)的位置略有差別，但大致都在x／L＝0.6～0.8之間.總壓越高流動(dòng)凍結(jié)出現(xiàn)得越晚，這是由于壓力增加（密度增大）增大了粒子碰撞頻率，使流動(dòng)的非平衡程度降低.振動(dòng)基態(tài)的CN2（0）在膨脹過(guò)程中隨溫度降低而增加（見圖2（a）），而其它較高振動(dòng)能級(jí)的CN2（v）變化趨勢(shì)相反（見圖2（b））.從CN2（v）凍結(jié)值來(lái)看，總壓越大，基態(tài)及低能級(jí)的凍結(jié)值越高，而高能級(jí)的凍結(jié)值越低，這與總壓升高推遲流動(dòng)凍結(jié)的特點(diǎn)是一致的.

圖1 T0＝8 000 K條件下CN沿噴管軸線分布Fig.1 Distribution of CNalong the axis of nozzle at T0＝8 000 K

圖2 T0＝8 000 K條件下CN2（v）沿噴管軸線分布Fig.2 Distribution of CN2（v）along the axis of nozzle at T0＝8 000 K

圖2反映了總壓條件對(duì)CN2（v）在膨脹過(guò)程中變化規(guī)律的影響，但由于不同總壓條件下N2／N混合物的組成不同，質(zhì)量分?jǐn)?shù)還不能反映不同振動(dòng)能級(jí)上粒子分布.圖3分別給出了基態(tài)和第30振動(dòng)能級(jí)的N2（v）粒子數(shù)在總的N2粒子數(shù)中所占比例沿軸線的分布.由圖3（a）可見，在流動(dòng)凍結(jié)之前，p0＝1 atm條件下N2（0）所占比例最高，p0＝10 atm條件下最低，這對(duì)應(yīng)于總壓較低時(shí)氣流溫度較低，更多的粒子通過(guò)躍遷進(jìn)入低振動(dòng)能級(jí)；從N2（0）比例的凍結(jié)值來(lái)看，則是p0＝10 atm條件下最高，p0＝1 atm條件下最低.這是由于壓力較低條件下能級(jí)分布較早地進(jìn)入凍結(jié)狀態(tài).較高振動(dòng)能級(jí)的粒子比例分布的變化特點(diǎn)及總壓條件的影響則與基態(tài)能級(jí)情況相反，這可由圖3（b）看出.

圖3 T0＝8 000 K條件下N2（v）占N2比例沿噴管軸線分布Fig.3 Distribution of ratio of N2（v）in N2along the axis of nozzle at T0＝8 000 K

2.1.2 噴管不同位置的振動(dòng)能級(jí)分布

圖4 T0＝8 000 K條件下CN2（v）分布與Boltzmann分布Fig.4 Distributions of CN2（v）and Boltzmann distribution at T0＝8 000 K

圖4分別給出了p0＝1 atm、5 atm條件下，噴管軸線位置為x／L＝0.2、0.5、0.8處的各振動(dòng)能級(jí)上的質(zhì)量分?jǐn)?shù).圖中實(shí)線代表態(tài)－態(tài)模型計(jì)算得到的實(shí)際振動(dòng)能級(jí)分布，虛線代表該位置處流場(chǎng)平動(dòng)溫度下的Boltzmann分布.可見在x／L＝0.2處，兩分布基本重合，說(shuō)明此處流動(dòng)仍接近振動(dòng)平衡，這是由于流動(dòng)速度還很低，溫度壓力和駐室相比也沒有明顯下降，流動(dòng)過(guò)程中分子間可進(jìn)行足夠多的碰撞以使振動(dòng)能級(jí)分布達(dá)到當(dāng)?shù)販囟认碌钠胶鈶B(tài).在x／L＝0.5（喉道）處，p0＝1 atm條件下的兩種分布出現(xiàn)差別，高能級(jí)上N2（v）的質(zhì)量分?jǐn)?shù)要高于當(dāng)?shù)販囟认碌腂oltzmann分布值，這表明氣流在膨脹降溫和加速過(guò)程中，粒子并沒有及時(shí)完成由高能級(jí)向低能級(jí)的躍遷，流動(dòng)出現(xiàn)振動(dòng)非平衡；而在p0＝5 atm的條件下，兩種分布仍基本重合，這是由于壓力增加提高了碰撞頻率，推遲了噴管流動(dòng)中非平衡的出現(xiàn).在x／L＝0.8處，實(shí)際振動(dòng)能級(jí)分布與當(dāng)?shù)販囟认翨oltzmann分布的差別在p0＝1 atm（圖4（a））、5 atm（圖4（b））和10 atm（圖不再給出）條件下都非常明顯，除基態(tài)外，各能級(jí)N2（v）質(zhì)量分?jǐn)?shù)都要高于Boltzmann分布值，并且總壓越低，差別越大.這是因?yàn)楹淼老掠螝饬鳒囟群蛪毫Γ芏龋┑慕档团c速度的增加都非常迅速，粒子間碰撞頻率明顯降低，非平衡程度進(jìn)一步加劇，流動(dòng)甚至趨近于化學(xué)和熱力學(xué)凍結(jié).

圖4明顯地反映出噴管流動(dòng)中的振動(dòng)非平衡現(xiàn)象.雙溫度或多溫度模型通過(guò)引入表征非平衡振動(dòng)能的振動(dòng)溫度描述振動(dòng)非平衡，隱含假設(shè)振動(dòng)能級(jí)分布滿足振動(dòng)溫度下的Boltzmann分布.這里對(duì)該假設(shè)的合理性做進(jìn)一步分析.首先根據(jù)本文采用了態(tài)－態(tài)模型計(jì)算得到的細(xì)致振動(dòng)能級(jí)分布和各能級(jí)能值求得總的非平衡振動(dòng)能，再采用和雙溫度模型中類似的辦法折算出振動(dòng)溫度，然后比較本文得到的振動(dòng)能級(jí)分布和該振動(dòng)溫度下的Boltzmann分布.

圖5給出了T0＝8 000 K，p0＝1 atm、5 atm和10 atm條件下，平動(dòng)溫度和振動(dòng)溫度沿噴管軸線的分布.可見喉道前和喉道后部分區(qū)域振動(dòng)溫度與平動(dòng)溫度都是基本重和的（高總壓條件下重合區(qū)域更大），在氣流繼續(xù)流向下游的過(guò)程中二者分離，平動(dòng)溫度仍然有明顯下降，而振動(dòng)溫度下降緩慢并趨于凍結(jié)不變.由此可知引入振動(dòng)溫度確實(shí)能在一定程度上方便地反映出振動(dòng)非平衡現(xiàn)象，但是否能夠反映出細(xì)致的非平衡能級(jí)分布還需繼續(xù)考察.

圖6對(duì)比了T0＝8 000 K，p0＝1 atm、5 atm條件下采用了態(tài)－態(tài)模型算得的實(shí)際振動(dòng)能級(jí)分布和振動(dòng)溫度下的Boltzmann分布.可見二者有明顯偏離，與圖4非常類似.當(dāng)然，在x／L＝0.8處，實(shí)際振動(dòng)能級(jí)分布與振動(dòng)溫度下Boltzmann分布的差別要比與平動(dòng)溫度下Boltzmann分布的差別稍小.這和喉道前及喉道后近下游平動(dòng)溫度與振動(dòng)溫度重合而在遠(yuǎn)下游處振動(dòng)溫度凍結(jié)在高于平動(dòng)溫度值的特點(diǎn)（見圖5）是一致的.該現(xiàn)象表明，對(duì)實(shí)際的非平衡流動(dòng)，即使引入振動(dòng)溫度表征非平衡振動(dòng)能，也不能簡(jiǎn)單地用振動(dòng)溫度下的Boltzmann分布描述細(xì)致的能級(jí)分布.

圖5 T0＝8 000 K條件下平動(dòng)溫度和振動(dòng)溫度沿噴管軸線分布Fig.5 Distributions of translational and vibrational temperature along the axis of nozzle at T0＝8 000 K

圖6 T0＝8 000 K條件下CN2（v）分布與振動(dòng)溫度下Boltzmann分布Fig.6 Distributions of CN2（v）and Boltzmann distribution at vibrational temperature，T0＝8 000 K

2.2 T0＝3 000 K算例

該算例條件下N2的離解非常微弱.計(jì)算中發(fā)現(xiàn)N質(zhì)量分?jǐn)?shù)在喉道后出現(xiàn)振蕩，并略有上升（見圖7）.這可能是因?yàn)樗憷龡l件下N為微量組元，而68個(gè)振動(dòng)能級(jí)上的N2（v）的離解或其逆反應(yīng)均可能影響到它的質(zhì)量分?jǐn)?shù)，容易引起數(shù)值振蕩.

根據(jù)這個(gè)情況，對(duì)該算例條件另外進(jìn)行了只計(jì)及N2的vTm和vv過(guò)程的噴管非平衡流場(chǎng)計(jì)算，分析化學(xué)反應(yīng)和vTa過(guò)程對(duì)宏觀流動(dòng)特性和微觀粒子能級(jí)分布的影響.

圖8對(duì)比了計(jì)及所有振動(dòng)躍遷和化學(xué)反應(yīng)與只計(jì)及N2的vTm和vv躍遷過(guò)程算得的溫度沿噴管軸線分布，其中振動(dòng)溫度是根據(jù)非平衡振動(dòng)能折算得出.可見化學(xué)反應(yīng)和vTa過(guò)程對(duì)溫度分布的影響非常有限，除入口附近無(wú)化學(xué)反應(yīng)情況下算得的值略高外，其它位置是否計(jì)及化學(xué)反應(yīng)得到的結(jié)果基本重合.流速、壓力、密度特性也基本不受化學(xué)反應(yīng)和vTa過(guò)程影響.

圖7 T0＝3 000 K、p0＝1 atm條件下CN沿軸線分布Fig.7 Distribution of CNalong the axis of nozzle at T0＝3 000 K，p0＝1 atm

圖8 T0＝3 000 K、p0＝1 atm條件下溫度沿軸線分布（包含或不含化學(xué)反應(yīng)）Fig.8 Distribution of temperature along the axis with or without chemical reaction at T0＝3 000 K，p0＝1 atm

但是，是否計(jì)及化學(xué)反應(yīng)和vTa過(guò)程獲得的粒子能級(jí)分布結(jié)果有明顯差異.圖9給出了噴管3個(gè)位置上的振動(dòng)能級(jí)分布情況.不考慮化學(xué)反應(yīng)時(shí)，振動(dòng)能級(jí)分布與當(dāng)?shù)仄絼?dòng)溫度下的Boltzmann分布在喉道前及喉道附近完全重合，在下游偏離Boltzmann分布的程度也比考慮化學(xué)反應(yīng)的情況下弱得多.這主要是因?yàn)镹與N2高振動(dòng)能級(jí)間的躍遷速率較高，而高振動(dòng)能級(jí)的質(zhì)量分?jǐn)?shù)本身就很小，存在化學(xué)反應(yīng)時(shí)N的復(fù)合會(huì)對(duì)其產(chǎn)生較明顯的影響；而沒有化學(xué)反應(yīng)時(shí)，vTm和vv過(guò)程會(huì)使能級(jí)分布較快趨近于平衡態(tài).同時(shí)可以看到此時(shí)高振動(dòng)能級(jí)的質(zhì)量分?jǐn)?shù)很小，在10－10以下，所以雖然考慮化學(xué)反應(yīng)與否得到的高振動(dòng)能級(jí)分布有較明顯差異，但基本不影響宏觀流動(dòng)特性.

圖9 T0＝3 000 K、p0＝1 atm條件下CN2（v）分布與Boltzmann分布Fig.9 Distributions of CN2（v）and Boltzmann distribution at T0＝3 000 K，p0＝1 atm

3 結(jié)論

用耦合態(tài)－態(tài)模型與流動(dòng)控制方程數(shù)值研究準(zhǔn)一維熱化學(xué)非平衡噴管流，在駐室溫度3 000 K～10 000 K，壓力1 atm～10 atm條件下數(shù)值求解了N2／N非平衡噴管流動(dòng)，分析流動(dòng)宏觀和微觀特性，得到如下結(jié)論：

1）采用時(shí)間推進(jìn)法求解耦合態(tài)－態(tài)模型的準(zhǔn)一維非平衡噴管流，可有效解決非平衡噴管流求解時(shí)由于喉道非聲速等原因造成的入口條件試湊困難問(wèn)題.

2）噴管流動(dòng)中喉道后振動(dòng)非平衡現(xiàn)象顯著.雖然引入振動(dòng)溫度能夠表征非平衡振動(dòng)能，卻無(wú)法細(xì)致描述各振動(dòng)能級(jí)上的粒子分布，真實(shí)的粒子分布與振動(dòng)溫度下的Boltzmann分布仍有明顯差別，復(fù)合占優(yōu)的膨脹流中粒子偏離平衡分布的程度比離解占優(yōu)的激波后流動(dòng)中更強(qiáng).膨脹過(guò)程中不同振動(dòng)能級(jí)分布的變化特點(diǎn)不同，細(xì)致的微觀能級(jí)分布規(guī)律必須通過(guò)采用態(tài)－態(tài)模型獲得.

3）在較低的總溫條件（如T0＜3 000 K）下N2的離解非常微弱，化學(xué)反應(yīng)對(duì)流場(chǎng)宏觀特性的影響可以忽略，因此流動(dòng)計(jì)算時(shí)是否考慮化學(xué)反應(yīng)均可得到合理的流場(chǎng)宏觀參數(shù)分布.但是否計(jì)及化學(xué)反應(yīng)得到的高振動(dòng)能級(jí)分布偏離平衡分布的程度有明顯不同，本文初步分析認(rèn)為不計(jì)及化學(xué)反應(yīng)得到的結(jié)果更加合理，而計(jì)及微量化學(xué)反應(yīng)的計(jì)算結(jié)果具有較大的數(shù)值不確定性.該問(wèn)題在未來(lái)的研究中可以進(jìn)一步深入.

4）態(tài)－態(tài)模型較雙溫度或多溫度模型更具理論優(yōu)勢(shì)，但由于其計(jì)算量巨大，目前用于多維非平衡流動(dòng)模擬還有困難.耦合態(tài)－態(tài)模型進(jìn)行準(zhǔn)一維噴管流研究，一方面更深入了解了重要非平衡區(qū)域微觀過(guò)程對(duì)宏觀流動(dòng)特性的影響，可為改進(jìn)現(xiàn)有的宏觀熱化學(xué)模型提供思路，另一方面也為未來(lái)將態(tài)－態(tài)模型用于二維或三維非平衡流動(dòng)研究做了一定技術(shù)儲(chǔ)備.

［1］Park C.Problems of rate chemistry in the flight regimes of aeroassisted orbital transfer vehicles［J］.Progress in Astronautics and Aeronautics，1985，96：511－537.

［2］Capitelli M，Armenise I，Gorse C.State-to-state approach in the kinetics of air components under re-entry conditions［J］. Journal of Thermophysics and Heat Transfer，1997，11（4）：570－578.

［3］Silva M L，Guerra V，Loureiro J.State-resolved dissociation rate for extremely nonequilibrium atmospheric entries［J］.Journal of Thermophysics and Heat Transfer，2007，21（1）：40－49.

［4］Park C.The limits of two-temperature model［R］.AIAA 2010－911，2010.

［5］Colonna G，Armenise I，Bruno D，et al.Reduction of state-to-state kinetics to macroscopic models in hypersonic flows［J］．Journal of Thermophysics and Heat Transfer，2006，20（3）：477－486.

［6］Capitelli M，Armenise I，Colonna G，Koudriavtsev N，et al.Nonequilibrium vibrational kinetics during hypersonic flow of a solid body in nitrogen and its influence on the surface heat flux［J］.Plasma Chemistry and Plasma Processing，1995，15（3）：501－527.

［7］Armenise I，Capitelli M，Gorse C.Nitrogen nonequilibrium vibrational distribution and non-Arrhenius dissociation constants in hypersonic boundary layers［J］.Journal of Thermophysics and Heat Transfer，1998，12（1）：45－51.

［8］Colonna G，Tuttafesta M，Capitelle M，et al.NO formation in one-dimensional nozzle air flow with state-to-state nonequilibrium vibrational kinetics［R］.AIAA 98－2951，1998.

［9］Colonna G，Tuttafesta M，Capitelle M，et al.Non-Arrhenius no formation rate in one-dimensional nozzle air flow［J］.Journal of Thermophysics and Heat Transfer，1999，13（3）：372－375.

［10］戴東旭，楊學(xué)明.基元化學(xué)反應(yīng)態(tài)態(tài)動(dòng)力學(xué)研究［J］.化學(xué)進(jìn)展，2007，19（11）：1633－1645.

［11］徐丹，張威，曾明，等.態(tài)－態(tài)模型下N2／N系統(tǒng)的熱化學(xué)非平衡過(guò)程研究［J］.空氣動(dòng)力學(xué)學(xué)報(bào)，2014，32（3）：280－288.

［12］Colonna G，Pietanza L D，Capitelle M.Recombination-assisted nitrogen dissociation rates under nonequilibrium conditions ［J］.Journal of Thermophysics and Heat Transfer，2008，22（3）：399－406.

［13］Esposito F，Armenise I，Capitelli M.N-N2State to state vibrational-relaxation and dissociation rates based on quasiclassicalcalculations［J］.Chemical Physics，2006，331：1－8.

［14］柳軍.熱化學(xué)非平衡流及其輻射現(xiàn)象的實(shí)驗(yàn)和數(shù)值計(jì)算研究［D］.長(zhǎng)沙：國(guó)防科學(xué)技術(shù)大學(xué)，2004.

［15］徐丹.高溫?zé)峄瘜W(xué)非平衡模型研究［D］.長(zhǎng)沙：國(guó)防科學(xué)技術(shù)大學(xué)，2012.

Numerical Study of Thermochemical Nonequilibrium Nozzle Flow in State-to-State Model

XU Dan，ZENG Ming，ZHANG Wei，LIU Jun
（College of Aerospace Science and Engineering，National University of Defense Technology，Changsha 410073，China）

A state-to-state model is used to simulate thermochemical nonequilibrium nozzle flow of N2／N mixture.Quasi-onedimensional nozzle flow at reservoir temperature from 3 000 K to 8 000 K and pressure from 1 atm to 10 atm is simulated numerically to obtain flow properties，especially detail vibrational population distribution.Evolution of macro-and micro-properties of nonequilbrium flow during expanding is analyzed.Rationality of two-temperature or multi-temperature approach is investigated for recombination dominating flows.

state-to-state model；thermochemical nonequilibrium；vibrational population distribution；hypersonic nozzle flow；numerical simulation

date：2013－10－09；Revised date：2014－02－21

V211.22

A

2013－10－09；

2014－02－21

國(guó)家自然科學(xué)基金（11102231）資助項(xiàng)目

徐丹（1987－），男，遼寧鳳城，博士生，主要從事計(jì)算流體力學(xué)、高溫氣體動(dòng)力學(xué)研究，E-mail：449761494＠qq.com

1001-246X（2014）05-0531-08