甲基肼/四氧化二氮反應化學動力學模型構建及分析

巴延濤 侯凌云,* 毛曉芳 汪鳳山

(1清華大學航天航空學院,北京100084;2北京控制工程研究所,北京100190)

1 引言

甲基肼(MMH)和四氧化二氮(NTO)由于具有優良的物理化學性質,被廣泛用作雙組元小推力液體火箭發動機推進劑.研究MMH/NTO反應機理對理解液體火箭發動機內部工作過程,特別是燃燒室內燃燒過程具有重要意義.然而目前國內外對MMH/NTO反應機理的研究還很少,主要原因是MMH/NTO有劇毒,且由于其具有很高的反應活性,即使在極其稀薄的條件下仍能迅速反應,因此制備MMH/NTO預混氣難度很大,所以,通過實驗手段研究MMH/NTO反應動力學難度很大.

MMH/NTO反應過程不同于其他燃料的燃燒過程:MMH和NTO一旦接觸,就會迅速發生反應,反應釋放的熱量使系統溫度升高.1-8在著火之前,MMH和NTO會分別發生熱分解反應.值得注意的是,不同于普通的分解反應,MMH的分解反應會釋放出大量的熱,這部分釋放的熱量也會促進系統自燃著火.9-11在自燃著火發生之后,MMH/NTO反應系統會發生高溫狀態下的一系列反應.目前公開發表的MMH/NTO反應機理研究很少,研究者大多依據量子化學原理,使用計算機模擬的方法討論分析MMH/NTO的反應特性,主要集中在自燃著火階段.Frank等12模擬了MMH/NTO在自燃著火之前的反應過程,指出甲基二氮烯(CH3N＝NH)為該階段主要產物.Ishikawa等13模擬了NO2分子奪去MMH分子上H原子的過程,指出有CH3NHNH,CH2NHNH,H2C＝NNH2和CH3NHN參與的反應應該被加入到MMH/NTO系統的化學反應動力學研究中去,并通過分析部分中間產物的結構穩定性,指出了可能的反應路徑,但沒有給出具體的反應.Catoire等14利用激波管對MMH在不同壓力和溫度范圍內的熱解反應進行了實驗研究,得到了99步的MMH熱解機理和14步簡化機理.關于MMH/NTO自燃著火過程,Catoire等2通過量子化學手段得到了MMH/NTO的82組分、403個基元反應的詳細機理,并通過靈敏度分析方法得到了影響反應物自燃著火過程的關鍵反應,但沒有給出該機理的細節.另外該機理規模過于龐大,無法應用于發動機內流動燃燒過程的數值模擬.聶萬勝15,16和Knab17等研究自燃推進劑MMH/NTO火箭發動機不穩定燃燒時采用了一步總包反應,Xu等18對空間小推力室進行數值模擬時分別采用了一步和四步總包反應機理,這兩種總包機理基于經驗得出,多應用于自然推進劑火箭發動機穩態流動燃燒的數值模擬中,無法描述MMH/NTO的自燃著火過程.因此,構建一個能夠反映MMH/NTO自燃著火過程的簡化反應動力學模型具有重要的意義.

本文在現有研究的基礎上構建一個能夠準確描述MMH/NTO反應過程的簡化反應動力學模型,為實現MMH/NTO化學反應與計算流體力學耦合提供有價值的反應動力學模型.

2 構建

MMH/NTO反應歷程按照先后順序劃分為冷反應、熱分解反應、自燃著火反應以及高溫反應四個階段,2隨著反應進行,這四個階段各自起到不同的作用.另外本文還考慮了燃料MMH的分解、氧化劑NTO的離解、NO2的分解、中間產物CH3的氧化及HONO的分解等反應.

2.1 冷反應

冷反應特指從MMH/NTO充分混合之后到發生自燃點火之前的一段過程,由于在該過程中反應系統的溫度還較低,一般認為是在常溫(298 K)到400 K之間,所以稱之為冷反應階段.

MMH與NO2一經接觸,MMH分子中間N原子上的H原子會被強氧化劑NO2逐步奪去,生成中間產物CH3NNH2、CH3NNH、CH3N2以及 HONO,該階段反應的活化能很低,能夠在較低的溫度下(甚至是273 K以下)自發進行,并放出大量的熱,該過程產生的中間產物HONO不穩定,會繼續分解為小分子.冷反應階段的重要反應主要是R1,R2和R3(參見表1).這三個反應的熱力學參數和動力學參數從文獻11,14,19,20中得到,本文對R1和R3的活化能數值做了適當的調整.Thaxton等21對以下反應:

進行了研究,給出了該反應的活化能數值,為104.2 kJ·mol-1,在此基礎上Catoire等2考慮到NH3分子中的N―H鍵(鍵能447.7 kJ·mol-1)相比MMH分子中間N原子的N―H鍵(368.2 kJ·mol-1)更穩定,鍵能更高,給出了R1反應的活化能(24.7 kJ·mol-1),其僅考慮了鍵能的因素,結果不夠合理.由于NO2分子有未成對的電子,化學活性很高,若認為NO2為活性自由基分子,應用活化能的估算經驗公式,22MMH+NO2反應活化能取為MMH中間N原子上N―H鍵鍵能的6%左右,約為21.3 kJ·mol-1.本文考慮到NH3分子為正四面體結構,其與NO2分子發生碰撞時任意一個N―H鍵均可能發生斷裂,而MMH分子為直鏈結構,中間N原子上存在的N―H鍵僅存在于MMH分子的局部,雖然其相比NH3分子中的N―H鍵來說更脆弱,鍵能更低,但是該鍵與NO2分子發生碰撞并斷裂的幾率較小.本文根據雙分子反應碰撞理論,碰撞是否發生反應的概率與反應活化能和碰撞幾率相關,即應增大文獻2中所取的活化能數值,所以在由經驗公式得到的數值(21.3 kJ·mol-1)的基礎上,又參考了文獻2提出的數值(24.7 kJ·mol-1),線性增大該反應活化能的數值,取為28 kJ·mol-1.由于CH3NNH2和CH3NNH等自由基的反應尚無相關研究,Catoire等2通過類比以下反應:

表1 Mech23反應動力學模型Table 1 Mech23 chemical kinetic model

得到了R2的動力學參數,并認為R3的動力學參數與R1相同,本文遵循該思路,得到R2和R3的動力學參數.

2.2 熱分解反應

MMH具有很高的反應活性,在較高溫度下會發生分解,這就是MMH的熱分解反應.另外,冷反應階段產生的HONO會在高溫條件下分解成小分子.其中MMH的熱分解反應會釋放出大量的熱,對推進劑自燃點火起到重要作用.熱分解反應階段的主要反應為R4和R5(參見表1).其中R4參考Catoire等23的研究,該反應正向進行時能釋放大量的熱,這部分熱量將使得系統溫度大幅升高,導致后續自燃著火現象的發生.由于該機理(機理含有23種組分,因此簡稱為Mech23機理)缺少許多中間組分的反應路徑,導致系統達到平衡時溫度較低,且著火延遲時間變長,為了更好地體現燃料自燃溫升的全過程,本文適當降低了該反應的活化能,以促進該反應正向進行,使系統能夠迅速著火,且達到平衡時溫度較合理,經過嘗試,最終取R4的活化能為文獻23所給數值的約1/3.由于HONO在高溫下不能穩定存在,所以增加反應R5,反應中M稱為第三體,為促進HONO轉化為小分子產物OH和NO提供能量,該反應的動力學參數參照Atkinson等19的相關研究,并做了適當調整,一方面提高HONO的轉化率,另一方面不至于使R5進行得過快,吸收大量的熱,而使得系統溫升緩慢,著火延遲時間變長.

Mech23中涉及的分解反應還有NTO以及NO2的分解(具體反應方程式參見表1中R6和R7).

2.3 自燃著火反應

自燃著火是MMH/NTO反應歷程的關鍵之處,也是MMH/NTO雙組元液體燃料的一大特色優勢.隨著上述兩個反應階段不斷產生和積累熱量,反應系統擁有了足夠的能量來激發著火過程,著火過程經歷的時間很短,在這個過程中整個反應系統的溫度和壓力驟然升高,其中所涉及的重要反應為R1、R3和R4(見表1).

這三個反應對于MMH/NTO自燃著火起著重要作用.首先因為反應物分子MMH和CH3NNH化學穩定性低,而且NO2分子外層軌道上有一個未成對的電子,氧化性很強,所以這兩個反應極易發生,甚至在很低的溫度(0°C以下)都能夠迅速發生,正是由于這兩個反應發生的活化能較低,才會在低溫條件下迅速發生,并引發一系列后續反應,使得MMH/NTO能夠發生自燃著火并建立高溫,也可以說,R1和R3是整個反應過程的起始環節.對于反應R4來說,由于其正反應具有很強的放熱效應,為反應物系統熱量的積累以及著火階段溫度的劇增起到了至關重要的推動作用.

2.4 高溫反應

在高溫反應階段,MMH/NTO系統主要發生以下反應:MMH的分解以及MMH和NTO的燃燒.其中MMH的分解反應前已述及,MMH/NTO的氧化過程中冷反應階段涉及的反應在此也不再贅述,這里主要給出兩條重要的CH3的氧化過程.

由于反應系統中NO2的分解產生了O原子,所以加入了CH3與O原子的一系列反應:

又由于NO2與HONO的分解產生了具有較強氧化性的OH自由基和NO自由基,所以考慮在本模型中加入OH、NO與CH3的反應.高溫反應階段中添加了許多重要的小分子(如CH3、OH、O、NO等)的反應路徑,主要是為了消耗前三個反應階段生成的中間產物,使得整個反應系統最終轉化為H2O、N2、H2以及NO等較穩定的小分子,釋放反應中間產物的化學能,使系統趨于穩定.高溫階段的反應均直接取自甲烷空氣氧化反應機理GRI Mech 2.11中,選取反應時遵循文獻24中提到的原則,即對于競爭反應,選取反應活化能較低的反應路徑,以盡可能地降低組分數,減小機理的規模.

綜合以上所述,得到了Mech23的所有反應,如表1所示.

3 驗證及分析

上述構建的Mech23機理基本涵蓋了MMH/NTO反應的整個過程,由于目前尚未有公開發表的對該反應的實驗研究,我們針對文獻中給出的理論分析和計算機模擬結果,采用CHEMKIN中平衡反應器模型和封閉全混均質反應器模型對Mech23簡化機理進行了驗證分析.

對于MMH/NTO反應系統來說,著火延遲時間刻畫了從推進劑噴入燃燒室到產生正推力的延遲時間,該指標具有重要意義.采用構建的簡化機理計算的著火延遲時間與詳細反應機理以及理論分析獲得的數據進行對比分析.

3.1 驗證

3.1.1 著火延遲時間驗證

為了驗證Mech23反應機理的準確性,對其在封閉全混均質反應器中發生的自燃著火過程進行了模擬,同時為模擬地面試驗環境,選取反應工況為初始溫度T0=298 K,初始壓力p0=2.42×104Pa,采用定容過程,得到了著火溫升曲線.

圖1 不同機理著火溫升曲線對比Fig.1 Comparison of temperature rising using different mechanisms

圖1為不同機理著火溫升曲線對比.在圖1中,實線為Mech23機理的溫升曲線,著火延遲時間為3.2 ms;點畫線為Catoire等2利用量子化學方法得到的詳細反應機理的溫升曲線,計算得到著火延遲時間為3.6 ms;虛線是Catoire等僅考慮了點火前產物(CH3N(NH2)NO2和CH3N(NH2)ONO)作用,在詳細機理基礎上保留這兩種組分相關的反應,刪除了與它們無關的反應而得到的機理1的溫升曲線,著火延遲時間為1.5 ms;點線是Catorie等去除了與CH3N(NH2)NO2和CH3N(NH2)ONO有關組分參與的反應提出的機理2的溫升曲線,只考慮了點火前產物(CH3N(NH2)NO2和CH3N(NH2)ONO)2的影響,得到的著火延遲時間為6.3 ms.另外,根據Seamans等25的熱爆炸理論計算得到的著火延遲時間為3.6 ms.

依據Mech23機理得到的著火延遲時間與詳細反應機理2和通過熱爆炸理論25得到的著火延遲時間吻合較好,從圖1中還可以看出Mech23簡化機理明顯優于Catorie提出的簡化機理1和簡化機理2,從而一定程度上驗證了Mech23機理的準確性.

3.1.2 壓力對著火延遲時間影響驗證

圖2主要說明了不同初始溫度和不同反應物比例條件下壓力對著火延遲時間的影響,其中離散點為應用Mech23機理計算得到的值,對比了Seamans等25依據熱爆炸理論得到的著火延遲時間隨壓力的變化曲線.從圖2可以看出對應不同線型的三種工況,隨著反應初始壓力的減小,MMH/NTO著火延遲時間呈負指數規律迅速增大,應用Mech23機理得到的點火延遲時間隨初始壓力的變化規律與Seamans等25得到的結果比較吻合.

3.2 分析

3.2.1 不同壓力對著火過程的影響

圖2 不同初始壓力下著火延遲時間變化對比Fig.2 Comparison of ignition delay time under different initial pressures

圖3給出了氧燃混合比(NTO/MMH配比)為1.65(質量比)的條件下,推進劑在三種不同恒定壓力條件下的著火溫升曲線.由圖3可以看出,推進劑在標準工況下經歷定壓燃燒過程的著火延遲時間為0.7 ms,之后系統很快達到化學平衡,達到平衡時系統溫度穩定在3000 K左右.另外,隨著燃燒室壓力的下降,推進劑著火延遲時間會延長.經過試算,當壓力下降到0.31 MPa時,在2 ms計算時長內推進劑無法自燃著火,從而得到了著火的低壓極限.毛曉芳等26模擬低壓狀態下小推力器脈沖著火特性時指出,當供應管路壓力為0.4 MPa,對應穩態燃燒室壓力為0.3 MPa時,MMH/NTO推進劑在小脈沖寬度內自燃著火困難,這與本文計算的著火低壓極限吻合.因此若要使該發動機在最小脈沖時間內產生啟動,應保證推進劑供應管路壓力在0.4 MPa以上.

圖3 不同壓力下MMH/NTO著火溫升曲線Fig.3 Temperature rising profile of MMH/NTO under different pressures

對應三種不同的壓力條件,反應物系統達到平衡狀態時,其產物含量如圖4所示,在穩態工況下,對應三個不同的壓力值,系統達到化學平衡時存在CH3NNH、H2、H、OH、H2O、CO2、NO、N2和CO等9種組分,它們的摩爾含量之和均達到99%以上.從圖中可以看出,隨著壓力的升高,CO2、N2、H2等相對穩定的小分子含量增多,使得系統的平衡溫度升高,在較高的溫度下,H2O和NO等分子發生分解,含量下降,CH3NNH作為重要反應(R4,見表1)的正向反應產物,隨著系統平衡溫度的升高,含量減少.

3.2.2 不同氧燃比對著火過程的影響

圖5給出了恒定穩態燃燒室壓力0.9 MPa條件下,推進劑在7種不同氧燃比條件下的著火溫升曲線.在0.9 MPa定壓條件下,推進劑氧燃比從0.4變化到3.0時,著火延遲時間單調增加,這是由于隨著燃料(MMH)相對含量的減少,推進劑由富燃向貧燃的狀態變化,導致著火困難.另外還可以看出,系統的平衡溫度呈現先升高后下降的趨勢,這可以從圖6中得到解釋.

圖4 不同壓力下平衡組分含量圖Fig.4 Contents of equilibrium species under different pressures

圖5 不同NTO/MMH配比下溫升曲線Fig.5 Temperature rising profile under different NTO/MMH mass ratios

圖6 不同NTO/MMH配比下平衡組分含量圖Fig.6 Contents of equilibrium species under different MMH/NTO mass ratios

圖6給出了不同氧燃比下產物組成分析,當氧燃比從0.4升高到3時,H原子含量呈現先增后減的趨勢,從反應機理分析可知,H的生成總伴隨著H2O和CO2等完全燃燒產物的生成,所以H原子濃度越高,表明燃燒越完全,反應系統的平衡溫度隨之升高;NO則呈現出先減后增的趨勢,從反應機理中發現,NO來源于NO2的分解,該反應活化能很高,促使其發生需要消耗非常多的熱量,導致反應系統溫度下降明顯,因此NO濃度較低時反應系統溫度較高;對于非常穩定的N2分子,其化學反應活性很低,較高的N2濃度必然伴隨著較高的系統溫度;在3000 K左右的高溫下,相對于CO2分子來說,CO分子更能夠穩定存在,所以CO也呈現了濃度先增后減的趨勢.

綜上可知,系統溫度隨著氧燃比從小到大會呈現先升高后下降的趨勢,與系統達到平衡時不同的平衡組分濃度相關.

3.2.3 反應靈敏度分析

為了弄清哪幾個反應在著火階段發揮重要作用,對Mech23機理的著火延遲時間進行了靈敏度分析研究,使用以下工況:反應物配比為MMH+2.5 NTO(1.0275 N2O4+1.4725 NO2),初始溫度T0=298 K,初始壓力p0=2.42×104Pa,采用定容封閉全混反應器模型.

圖7 Mech23機理著火延遲時間靈敏度Fig.7 Sensitivity to ignition delay time of Mech23

靈敏度分析的方法是,分別將每個反應的速率常數k乘以2,其他反應不變,然后在相同計算工況下得到新的點火延遲時間τ.在這里,著火延遲時間的靈敏度定義成(τ-τ0)/τ0,其中τ0是Mech23機理相同計算工況下得到的推進劑點火延遲時間.圖7給出了點火延遲時間的靈敏度譜圖,在該工況下反應R1、R3和R4對著火延遲時間有較高的靈敏度,其他反應對著火延遲時間的靈敏度在1%以下,負值表示增大該反應的反應速率可以縮短著火延遲時間,加快著火.R1對著火延遲時間具有最高的靈敏度,原因是該反應是整個MMH/NTO反應過程的起始反應,為后續所有反應提供了反應物和初始能量,因此該反應發生的速率和方向對系統的著火延遲時間具有決定性的作用.反應R3生成大量的CH3N2自由基,且該反應是CH3N2自由基的唯一來源,從Mech23機理可以看出,CH3N2自由基是CH3的唯一來源,且能分解生成穩定小分子N2,對反應物系統升溫有重要影響,因此對系統著火延遲時間較敏感.對于反應R4,前已述及,該反應會釋放出大量的熱,促使了推進劑發生自燃著火,此處的敏感性分析也得到了相同的結論.

4 結論

本文構建了一個包含23種組分和20個基元反應MMH/NTO簡化反應機理,對其進行了驗證,證明該機理能夠準確地描述MMH/NTO燃料自燃著火的整個過程.在不同壓力和氧燃比以及敏感性分析中得到以下結論.

(1)當反應系統初始壓力從0.5 MPa增加到1.3 MPa,MMH/NTO點火延遲時間減少,平衡溫度升高;當壓力下降到0.31 MPa時,在2 ms計算時長內推進劑無法自燃著火,從而得到了著火的低壓極限;

(2)當氧燃比從0.4增加到3時,MMH/NTO點火延遲時間延長,平衡溫度先升高后降低;

(3)反應R1、R3和R4對著火延遲時間有較高的靈敏度.

(1)Yuan,T.;Chen,C.;Huang,B.The Comparison of the Hot-Fire and Cold-Flow Observations of NTO/MMH Impinging Combustion.45thAIAAAerospace Sciences Meeting,Nevada,January 8-11,2007.

(2) Catoire,L.;Chaumeix,N.;Pichon,S.;Paillard,C.J.Propul.Power2006,20(1),120.

(3) Catoire,L.;Swihart,M.T.J.Propul.Power2002,18(6),1242.doi:10.2514/2.6059

(4)Leca,C.;Boh,P.MON and MMH Pressure Surges for a Simplified Propellant Feed System.Third International Conference on Spacecraft Propulsion,Cannes,France,October 10-13,2000;ESA:Noordwijk,Netherlands,2000.

(5) Solomon,Y.;DeFini,S.J.;Pourpoint,T.L.;Anderson,W.E.Gelled MMH Hypergolic Droplet Investigation.47thAIAA/ASME/SAE/ASEE Joint Propulsion Conference.San Diego,California,US,August 3,2011.

(6) Cai,G.Q.;Yu,Q.S.;Dong,N.;Wu,N.C.Acta Phys.-Chim.Sin.1991,7(3),333.[蔡國強,俞慶森,董 南,吳念慈.物理化學學報,1991,7(3),333.]doi:10.3866/PKU.WHXB19910314

(7) Kondrikov,B.N.;Smirnov,S.P.;Minakin,A.V.;Doherty,R.M.Propell.Explos.Pyrot.2004,29(1),27.doi:10.1002/prep.200400027

(8) Osmont,A.;Catoire,L;Klapotke,T.M.;Vaghjiani,G.L.;Swihart,M.T.Propell.Explos.Pyrot.2008,33(3),209.doi:10.1002/prep.200700213

(9) Sun,H.Y.;Law,C.K.J.Phys.Chem.A2007,111(19),3748.doi:10.1021/jp067591l

(10) Nusca,M.J.;Michaels,R.S.Progress in the Development of a Computational Model for theArmy′sImpinging Stream Vortex Engine.40thAIAA/ASME/SAE/ASEE Joint Propulsion Conference.Fort Lauderdale,Florida,US,July 11-14,2004.

(11) Ross,D.S.;Kirshen,N.A.;Hendry,D.G.Study of the Basic Kinetics of Decomposition of MMH and MHF and the Effects of Impurities on their Stability.National Techmical Information Service,June,1970.

(12)Frank,I.;Hammerl,A.;Klaptke,T.M.;Nonnenberg,C.;Zewen,H.Propell.Explos.Pyrot.2005,30(1),44.doi:10.1002/prep.200400084

(13) Ishikawa,Y.;McQuaid,M.J.J.Mol.Struct.-Theochem2007,818(1-3),119.doi:10.1016/j.theochem.2007.05.014

(14) Catoire,L.;Bassin,X.;Dupre,G.;Paillard,C.Experimental Study and Kinetic Modeling of the Thermal Decomposition of Gaseous Monomethylhydrazine.Application to Detonation Sensitivity.InShock Waves,15th International Colloquium on the Dynamics of Explosions and Reactive Systems;Univ.Colorado,Boulder,Co,July 30-August 4,1995;Springer Verlag:New York,September 1996;pp 139-146.

(15) Nie,W.S.Hypergolic Propellant Rocket Engine Combustion Stability Studies.Ph.D.Dissertation,National University of Defense Technology,Changsha,1998.[聶萬勝.自燃推進劑火箭發動機燃燒穩定性研究[D].長沙:國防科技大學,1998.]

(16) Zhuang,F.C.;Zhang,Z.G.;Nie,W.S.;Zou,Q.Journal of Propulsion Technology2001,22(2),155.[莊逢辰,張中光,聶萬勝,鄒 勤.推進技術,2001,22(2),155.]

(17) Knab,O.;Preclik,D.;Estublier,D.Flow Field Prediction within Liquid Cooled Combustion Chambers of Storable Bi-propellant Rocket Engines.34thAIAA/ASME/SAE/ASEE Joint Propulsion Conference.OH,US,July 13,1998.

(18)Xu,K.M.;Cai,G.B.;Zhou Y.T.Flow Field Numerical Study within Thrust Chamber ofAerospacesmall Liquid Propellant Rocket Engine.InAIAA57th InternationalAstronautical Congress;American Institute ofAeronautics andAstronautics Inc.:IAC,2006,October 2;pp 6510-6519.

(19)Atkinson,R.;Baulch,D.L.;Cox,R.A.;Hampson,R.F.,Jr.;Kerr,J.A.;Rossi,M.J.;Troe,J.J.Phys.Chem.Ref.Data1997,26(6),1329.doi:10.1063/1.556010

(20) Catoire,L.;Chaumeix,N.;Paillard,C.J.Propul.Power2004,20(1),87.doi:10.2514/1.9234

(21)Thaxton,A.G.;Hsu,C.C.;Lin,M.C.Int.J.Chem.Kinet.1997,29(4),245.doi:10.1002/(SICI)1097-4601(1997)29:4<245::AID-KIN2>3.0.CO;2-U

(22) Fu,X.C.;Chen,R.H.Physical Chemistry(Vol.2);People′s Education Press:Beijing,1982;p 213.[傅獻彩,陳瑞華.物理化學(下冊).北京:人民教育出版社,1982:213.]

(23) Catoire,L.;Ludwig,T.;Dupre,G.;Paillard,C.Symposium(International)on Combustion1998,No.2,2359.

(24)Zhai,G.H.;Wang,H.;Yang,H.F.;Ran,X.Q.;Wang,Y.B.;Wen,Z.Y.Acta Phys.-Chim.Sin.2001,17(4),348.[翟高紅,王 惠,楊海峰,冉新權,王育彬,文振翼.物理化學學報,2001,17(4),348.]doi:10.3866/PKU.WHXB20010414

(25) Seamans,T.F.;Vanpee,M.;Agosta,V.D.AIAA J.1967,5(9),1616.doi:10.2514/3.4259

(26)Mao,X.F.;Tang,Z.Y.;Cai,G.B.Journal of Aerospace Power2013,28(3),550.[毛曉芳,唐振宇,蔡國飆.航空動力學報,2013,28(3),550.]