劉伯偉，姜 毅

(北京理工大學(xué) 宇航學(xué)院，北京 100081)

0 引言

燃?xì)庹羝綇椛浔粡V泛用于陸基及潛射戰(zhàn)略導(dǎo)彈的發(fā)射系統(tǒng)中，與其他彈射形式相比，具有結(jié)構(gòu)簡(jiǎn)單、溫度適中、壓力輸出平穩(wěn)、輸出能量可調(diào)等優(yōu)點(diǎn)［1］。該方式通常使用液態(tài)水作為冷卻介質(zhì)，通過(guò)液態(tài)水的熱傳導(dǎo)和汽化吸熱兩種方式對(duì)燃?xì)膺M(jìn)行降溫，其發(fā)射過(guò)程非常復(fù)雜，涉及燃?xì)馍淞?、氣液兩相流、傳熱、汽化傳質(zhì)等多方面內(nèi)容。傳統(tǒng)的零維內(nèi)彈道算法不考慮流體參數(shù)沿管路軸線的變化，且無(wú)法給出流場(chǎng)中液相的流型，也就無(wú)法具體研究冷卻水的汽化過(guò)程及流場(chǎng)中各處的流動(dòng)狀態(tài)，而應(yīng)用計(jì)算流體力學(xué)技術(shù)(CFD)求解氣液兩相流場(chǎng)，可較好地解決這一問(wèn)題。

某型集中注水式燃?xì)庹羝麖椛溲b置采用高能雙基藥作為能量源，并使用液態(tài)水作為冷卻介質(zhì)，本文對(duì)其建立了二維軸對(duì)稱(chēng)兩相凍結(jié)流場(chǎng)計(jì)算模型，使用Mixture多相流模型求解氣液兩相凍結(jié)流場(chǎng)，使用組分輸運(yùn)模型將流場(chǎng)中的氣相工質(zhì)簡(jiǎn)化為包含燃?xì)?、水蒸汽和空?種組分的混合氣體，利用水的汽化模型模擬了水的汽化過(guò)程，利用動(dòng)網(wǎng)格技術(shù)模擬了發(fā)射過(guò)程中的流場(chǎng)變化，研究了汽化效應(yīng)對(duì)流場(chǎng)產(chǎn)生的影響。

1 物理模型

1.1 多相流模型中的組分輸運(yùn)方程

不考慮流場(chǎng)中的化學(xué)反應(yīng)［2］，組分輸運(yùn)模型中第i組分守恒方程統(tǒng)一形式為

與單一組分的流場(chǎng)守恒方程相比，式(1)中增加了組分i的質(zhì)量分?jǐn)?shù)項(xiàng)Yi。如果總的組分?jǐn)?shù)為N，則有

計(jì)算時(shí)，需要求解前N－1種組分的守恒方程，而第N種組分的質(zhì)量分?jǐn)?shù)可由式(2)求得。

當(dāng)把式(1)應(yīng)用于多相流模型中，對(duì)q相的第i種組分，式(1)變?yōu)?/p>

式中 αq為q相的體積分?jǐn)?shù)qjpi為由q相第j組分到p相第i組分的質(zhì)量轉(zhuǎn)移源項(xiàng)piqj依此類(lèi)推。

1.2 真實(shí)氣體模型

在導(dǎo)彈發(fā)射過(guò)程中，流場(chǎng)中除燃?xì)馔?，還存在液態(tài)水汽化而來(lái)的水蒸汽，以及初始時(shí)刻的少量空氣。當(dāng)壓力小于20 MPa、溫度大于1 400 K時(shí)，可將燃?xì)夂涂諝庾骼硐霘怏w處理，但由于水蒸汽離液態(tài)不遠(yuǎn)，將其簡(jiǎn)化為理想氣體具有一定偏差。因此，有必要將水蒸汽作真實(shí)氣體考慮［3］。

對(duì)水蒸汽應(yīng)用Soave-Redlich-Kwong真實(shí)氣體模型［4－6］:

其中

該模型需要提供3個(gè)參數(shù)，即臨界溫度Tc、臨界壓力pc、偏心因子ω。

1.3 水的汽化方程

根據(jù)水的飽和溫度計(jì)算水的汽化率，對(duì)計(jì)算域中各個(gè)網(wǎng)格內(nèi)的氣相和液相流體分別求解。當(dāng)混合相溫度大于水的飽和溫度，水吸收能量汽化為水蒸汽;當(dāng)混合相溫度小于水的飽和溫度時(shí)，水蒸汽釋放能量凝結(jié)為液態(tài)水［7］。

液態(tài)水汽化公式

水蒸汽凝結(jié)公式

某一單元格內(nèi)液態(tài)水的凈汽化率為

水汽化造成的能量變化為

式中 Sh為水汽化吸收的能量或水蒸汽凝結(jié)釋放的能量，當(dāng)m·為正，表示當(dāng)前單元格內(nèi)總體表現(xiàn)為液態(tài)水汽化吸熱，流場(chǎng)能量降低，Sh為負(fù)，反之亦同;Qlat為水的汽化潛熱，根據(jù)當(dāng)?shù)貕毫Σ轱柡退c和飽和蒸汽表得到。

2 計(jì)算模型

2.1 計(jì)算域

計(jì)算域采用二維軸對(duì)稱(chēng)模型，由噴管、水室、霧化器、低壓室、彈托等部分組成，如圖1所示。其中，水室中存放著冷卻水。整個(gè)流場(chǎng)均采用四邊形網(wǎng)格劃分，因霧化器開(kāi)孔較細(xì)小，故對(duì)霧化器開(kāi)孔附近的網(wǎng)格進(jìn)行了加密，如圖2所示。

圖1 計(jì)算域示意圖Fig.1 Diagram of computational domain

圖2 霧化器附近網(wǎng)格示意圖Fig.2 Diagram of mesh nearing atomizer

2.2 邊界條件

圖3為邊界條件示意圖。噴管入口處為壓力入口邊界，需要給定入口的總溫、總壓，其中總溫為3 000 K，總壓由發(fā)動(dòng)機(jī)空放試驗(yàn)得到，如圖3(b)所示;計(jì)算域下邊界為軸對(duì)稱(chēng)邊界;其余外邊界為壁面邊界，物面邊界采用無(wú)滑移絕熱壁面邊界條件，近壁面湍流計(jì)算采用標(biāo)準(zhǔn)壁面函數(shù)模型;其中，彈托壁面為運(yùn)動(dòng)邊界，計(jì)算時(shí)需要積分此面上的混合相工質(zhì)靜壓，從而得到導(dǎo)彈某一時(shí)刻的運(yùn)動(dòng)加速度，并結(jié)合動(dòng)網(wǎng)格技術(shù)實(shí)現(xiàn)計(jì)算域的變形。由于彈托壁面附近的網(wǎng)格非常規(guī)整，故使用動(dòng)態(tài)分層法生成新網(wǎng)格。

圖3 邊界條件示意圖Fig.3 Diagram of boundary conditions

2.3 初始條件

計(jì)算初始時(shí)刻由高壓室破膜時(shí)刻開(kāi)始，破膜壓力為 2 MPa，預(yù)加冷卻水量為 2.20 kg。

2.4 工況介紹

計(jì)算分為2種工況。其中，工況1未考慮汽化效應(yīng)，液態(tài)水僅以熱傳導(dǎo)的形式為燃?xì)饨禍?，并沒(méi)有水蒸汽生成，故將氣相工質(zhì)簡(jiǎn)化為單一組分的理想燃?xì)?，液相為不可壓縮的液態(tài)水;工況2引入汽化模型和組分輸運(yùn)模型，液態(tài)水汽化為水蒸汽，將氣相工質(zhì)看作燃?xì)?、水蒸汽和少量空?種組分的混合氣體，通過(guò)引入Soave-Redlich-Kwong模型，將水蒸汽看作真實(shí)氣體。

3 結(jié)果與分析

3.1 氣液兩相流場(chǎng)分布

以工況2為例，圖4為工況2初始時(shí)刻和離筒時(shí)刻的流場(chǎng)示意圖，對(duì)稱(chēng)軸上方為初始時(shí)刻流場(chǎng)示意圖，預(yù)加冷卻水加于水室中，如圖4所示;對(duì)稱(chēng)軸下方為離筒時(shí)刻流場(chǎng)，此時(shí)流場(chǎng)軸向長(zhǎng)度約為初始時(shí)刻的6.6倍。

圖4 發(fā)射前后流場(chǎng)對(duì)比示意圖Fig.4 Diagram of computational domain after and before launching

圖5為破膜后0～10 ms流場(chǎng)中的液態(tài)水和水蒸汽分布圖。圖5(a)中，燃?xì)鈱⒁簯B(tài)水沖向低壓室，經(jīng)過(guò)霧化器后，一方面液態(tài)水可較均勻地分布于整個(gè)低壓室，從而減少燃?xì)鈱?duì)彈托底部的直接沖擊;另一方面，可增大燃?xì)夂鸵簯B(tài)水的交界面，有利于液態(tài)水的汽化。

圖5(b)中，初始時(shí)刻流場(chǎng)中并沒(méi)有水蒸汽，隨著與燃?xì)饨唤缣幍囊簯B(tài)水被加熱，其溫度超過(guò)了飽和溫度，從而在交界面處被汽化為水蒸汽。由于生成的水蒸汽主要分布于燃?xì)馀c液態(tài)水的交界面處，故水蒸汽充當(dāng)了燃?xì)馀c液態(tài)水的能量傳遞中介。

圖5 前10 ms液態(tài)水和水蒸汽分布Fig.5 Distribution of water and vapor in 0～ 10 ms

空放試驗(yàn)中，利用高速攝影儀拍攝了霧化器后方流場(chǎng)。點(diǎn)火后，首先沖過(guò)霧化器的是水室中預(yù)加的液態(tài)水，隨后燃?xì)膺M(jìn)入流場(chǎng)，圖6為從中截取的液態(tài)水分布圖。由圖6可知，穿過(guò)霧化器后，原先集中于水室中的液態(tài)水變?yōu)樯⒉加诤蠓搅鲌?chǎng)，大部分液態(tài)水分布于流場(chǎng)的軸線附近，只有少部分液態(tài)水由霧化器側(cè)面的開(kāi)孔排出。圖5(a)的仿真結(jié)果與高速攝影中的液態(tài)水分布非常一致，證明Mixture模型能較好地模擬流場(chǎng)中液相的分布，驗(yàn)證了仿真結(jié)果的可靠性。

圖7為工況1和工況2中液態(tài)水的質(zhì)量變化曲線，為加以區(qū)分，工況1加實(shí)心方塊，工況2加空心方塊，下文類(lèi)同。工況1因?yàn)槲醇悠Ｐ?，水的質(zhì)量不變。工況2初始時(shí)刻水的汽化較慢，這是因?yàn)槌跏紩r(shí)刻燃?xì)馀c液態(tài)水的交界面較小;0.2 s后，隨著燃?xì)馀c液態(tài)水的混合越來(lái)越充分，液態(tài)水的汽化速度逐漸加快;0.48 s后，由于大部分液態(tài)水已經(jīng)汽化，燃?xì)馀c液態(tài)水的交界面開(kāi)始減小，水的汽化速度也開(kāi)始減慢;最終液態(tài)水剩余 0.208 kg，與預(yù)加的 2.20 kg相比，汽化率達(dá)到94.5%，汽化效果較好。

圖6 液態(tài)水分布的高速攝影圖片F(xiàn)ig.6 High speed photograph of water distribution

圖7 液態(tài)水的質(zhì)量Fig.7 Mass of water

3.2 汽化效應(yīng)對(duì)內(nèi)彈道性能的影響

圖8為低壓室中以質(zhì)量平均的混合工質(zhì)平均靜溫曲線。2種工況的溫度均在離筒時(shí)刻達(dá)到最大值，但工況1的溫度曲線在離筒前單調(diào)上升，且在離筒時(shí)刻達(dá)到最大值1 113 K;工況2的溫度曲線先上升，并在0.23 s達(dá)到極大值546 K，隨著汽化速度的加快，溫度轉(zhuǎn)而降低，0.45 s后剩余液態(tài)水較少，汽化速度減慢，溫度再次上升，并在離筒前達(dá)到最大值589 K?？深A(yù)見(jiàn)，適當(dāng)增加預(yù)加水量，可使工況2的溫度進(jìn)一步降低。美國(guó)MX導(dǎo)彈采用燃?xì)庹羝綇椛浜?，使工質(zhì)溫度降低到477～533 K［8］，工況2與此溫度范圍吻合，而工況1則偏差過(guò)大。由此可見(jiàn)，與只考慮熱傳導(dǎo)降溫相比，考慮汽化效應(yīng)能大幅降低流場(chǎng)的溫度。

圖9為彈托底部平均靜壓曲線，2種工況的壓力曲線變化趨勢(shì)基本相同，且均在0.21 s達(dá)到最大值，但工況2的壓力變化比工況1更加劇烈，上升與下降速度均快于工況1。0.48 s后，由于液態(tài)水的汽化減慢，流場(chǎng)溫度上升，工況2的壓力轉(zhuǎn)為上升，并使其離筒壓力高于工況1。與圖8的低壓室平均靜溫相比，汽化效應(yīng)對(duì)壓力的影響較小。這是因?yàn)椴豢紤]汽化效應(yīng)時(shí)，雖然流場(chǎng)中的氣相工質(zhì)沒(méi)有增加液態(tài)水汽化而來(lái)的水蒸汽，但同時(shí)也未考慮汽化所吸收的能量損失，造成仿真流場(chǎng)中的溫度升高，根據(jù)理想氣體狀態(tài)方程，溫度升高造成壓力相應(yīng)升高，故而壓力損失并不明顯。

圖8 低壓室平均靜溫Fig.8 Average static temperature of low-pressure chamber

圖9 彈托底部平均靜壓Fig.9 Average static pressure nearing sabot

圖10為彈托底部平均靜溫曲線，圖11為0.12 s兩工況的流場(chǎng)溫度云圖，對(duì)稱(chēng)軸之上為工況1，對(duì)稱(chēng)軸之下為工況2。0.12 s后，2種工況的高溫燃?xì)舛紝⒌蛪菏抑械囊簯B(tài)水排向兩側(cè)，對(duì)彈托產(chǎn)生了一定的直接沖擊，但工況1的沖擊作用更加明顯，最高溫度達(dá)到了987 K;而工況2由于考慮了汽化效應(yīng)，在燃?xì)夂鸵簯B(tài)水的交界面處，由于汽化吸熱對(duì)燃?xì)獾慕禍刈饔幂^明顯，同時(shí)汽化生成的水蒸汽比熱容要高于燃?xì)猓瑢?duì)降溫也有一定效果，工況2的最高溫度僅為764 K。

圖12與圖13分別為導(dǎo)彈加速度和速度曲線，結(jié)合圖9，由于彈托底部壓力的變化趨勢(shì)基本一致，故2種工況的加速度與速度曲線趨勢(shì)也相近，但由于工況2彈托底部壓力略高于工況1，故工況2的最大加速度與離筒速度均高于工況1。

圖10 彈托底部平均靜溫Fig.10 Average static temperature nearing sabot

圖11 0.12 s流場(chǎng)溫度云圖Fig.11 Temperature contours at 0.12 s

圖12 導(dǎo)彈加速度Fig.12 Acceleration of missile

圖13 導(dǎo)彈速度Fig.13 Velocity of missile

表1列出了兩種工況下部分內(nèi)彈道參數(shù)與對(duì)應(yīng)的試驗(yàn)值，考慮汽化效應(yīng)的工況2與試驗(yàn)值更加一致，誤差均在6%以內(nèi)，而沒(méi)有考慮汽化效應(yīng)的工況1，則誤差較大。

表1 內(nèi)彈道參數(shù)的對(duì)比Table 1 Comparison of interior ballistic parameters

4 結(jié)論

(1)與不考慮汽化效應(yīng)相比，考慮汽化效應(yīng)使仿真流場(chǎng)的最大平均溫度降低了524 K，彈托底部的最大平均溫度降低了223 K，二者降低的幅度分別達(dá)到47.1%和22.6%。雖然不考慮汽化效應(yīng)也能得到接近試驗(yàn)結(jié)果的壓力與速度曲線，但誤差相對(duì)較大，尤其溫度誤差過(guò)大，難以反映流場(chǎng)中的真實(shí)情況，從而無(wú)法對(duì)彈射裝置的設(shè)計(jì)與改進(jìn)提供可靠性指導(dǎo)。所以，在該類(lèi)問(wèn)題的仿真中必須考慮汽化效應(yīng)。

(2)結(jié)合Mixture模型與液態(tài)水汽化模型能較好地模擬發(fā)射過(guò)程中液態(tài)水的流型與汽化過(guò)程。

(3)Soave-Redlich-Kwong真實(shí)氣體模型能準(zhǔn)確地反映發(fā)射過(guò)程中水蒸汽的狀態(tài)變化，從而給出準(zhǔn)確的內(nèi)彈道結(jié)果。

［1］芮守禎，邢玉明.幾種導(dǎo)彈彈射動(dòng)力系統(tǒng)內(nèi)彈道性能比較［J］.北京航空航天大學(xué)學(xué)報(bào)，2009，35(6):766-770.

［2］袁曾鳳.火箭導(dǎo)彈彈射內(nèi)彈道學(xué)［M］.北京:北京工業(yè)學(xué)院出版社，1987.

［3］趙險(xiǎn)峰，王俊杰.潛地彈道導(dǎo)彈發(fā)射系統(tǒng)內(nèi)彈道學(xué)［M］.哈爾濱:哈爾濱工程大學(xué)出版社，2001.

［4］康玉東，孫冰.燃?xì)夥瞧胶饬髟偕鋮s流動(dòng)傳熱數(shù)值模擬［J］.推進(jìn)技術(shù)，2011，32(1):119-124.

［5］佘守賢.真實(shí)氣體與氣液相變［M］.北京:人民教育出版社，1963.

［6］Soave G.Equilibrium constants from a modified redlich-kwong equation of state［J］.Chemical Engineering Science，1972，27(6).

［7］LEE W H.Pressure iteration scheme for two-phase flow modeling［J］.Ieee/engineering in Medicine and Biology Society Annual Conference，1980(1).

［8］李廣裕.戰(zhàn)略導(dǎo)彈彈射技術(shù)的發(fā)展［J］.國(guó)外導(dǎo)彈與航天運(yùn)載器，1990(7):38-49.