基于五階WENO格式的燃氣在藥床中流動過程二維兩相流研究*

程申申，王 浩，薛 紹，陶如意

（南京理工大學能源與動力工程學院，江蘇 南京 210094）

膛內點傳火過程是一個復雜的高溫高壓傳熱傳質過程，包括點火管火藥燃燒、點火管燃氣流出以及膛內藥床的傳熱燃燒過程，對火炮的藥床燃燒以及后續推進過程有重要的影響。在實際實驗中，現有測試手段很難得到實際狀態下這一階段的實驗現象。因此，開展藥床燃氣流動的內彈道模擬對整個火炮發射過程的研究具有重要意義[1-3]。

近年來針對點火管和膛內火藥床燃燒的研究成果有很多。王浩等[4]通過點傳火特性實驗，對比了不同裝填條件下的傳火管多點壓力曲線及火焰傳播速率，分析了底火焰在傳火管中的傳輸特性，是早期對點傳火進行的實驗研究。王珊珊等[5-6]對點火管和中心管建立了一維兩相流模型進行內彈道兩相流研究，分析了不同結構點火管的點傳火性能。郭夢婷等[7]對中心炸管定容階段進行了二維兩相流研究，探究了中心管式拋撒機構的內彈道流場特性。但是以上研究都著眼于點火管內火藥燃燒和傳播，對于燃氣在膛內流動及火藥燃燒的研究較少。Miura 等[8]使用帶破裂膜片的燃燒裝置進行實驗，并采用兩相流體動力學代碼進行模擬計算，用于探索膛內負壓差產生的機理,并得出結論,壓差是由固體推進劑運動導致的，壓差波動取決于底火噴口的流量和噴口的位置，但實驗和模擬過程未對火焰發展及流場規律進行詳細描述和觀測。Jaramaz 等[9]對整個火炮膛內流動過程進行實驗和一維兩相流研究，但是未考慮徑向效應。程誠等對某制導炮彈二維兩相流內彈道性能進行研究，并采用高階黎曼近似模型求解兩相流方程，詳細研究了帶推進系統的整個內彈道過程中膛內流場特性及其發展[10-11]，但主要研究火藥床燃燒及彈丸推進過程，并未仔細對初始階段點傳火過程進行流場分析。因此，開展中心傳火管裝藥點火過程的理論與試驗研究，對進一步分析火炮發射起始階段膛內的點火過程燃氣流動規律，降低起始壓力波的產生，保證火炮射擊的安全性至關重要。

本文中，擬設計低壓可視化實驗裝置，觀察點火藥燃氣在藥床中的流動過程?？紤]到可視化裝置的承壓能力，采用假藥粒裝填結構，并結合二維兩相流數值模擬來預測這一階段的實際過程。使用五階精度的WENO格式，相比于經典的Mac Cormack 格式與TVD格式，其具有精度高、間斷分辨清晰，而且效率高的優點，是十分有效的算法[12-14]。以期研究結果對進一步開展火炮起始階段燃氣流動及壓力變化規律研究具有應用價值。

1 數理模型

1.1 物理模型

對膛內初始階段點傳火過程建立的計算模型如圖1所示。膛內填充由尼龍材料仿真實藥粒大小制成的假藥粒，中心為點傳火管，內裝小粒黑火藥，左端用端蓋密封，裝有底火以點燃點火管內火藥，右端用鋁合金膜片密封。

圖1 裝藥結構示意圖Fig.1 Schematic of the charge structure

根據上述模型建立真實藥粒運動燃燒和假藥粒運動的兩相流模型，為了簡化內彈道數值模擬的復雜性，需對實際物理過程進行適當假設，基本假設如下：

（1）采用雙流體模型假設，將發射藥顆粒群作為具有連續介質特性的擬流體來處理;

（2）忽略湍流影響和氣固兩相間黏性；

（3）火藥顆粒服從幾何燃燒和指數燃燒定律[2]；

（4）氣相狀態方程采用Nobel-Abel狀態方程[2]；

（5）假設底火內點火藥完全以燃氣形式加入點火管內，點火管采用一維兩相流計算，破膜后源項在膛內當地釋放；

（6）假藥粒不產生燃氣，尺寸與真實藥粒相同。

1.2 膛內軸對稱二維兩相流方程組

對于內彈道膛內氣-固兩相流過程，由于膛內沿身管軸線方向呈對稱形式，因此在內彈道膛內氣-固兩相流二維計算時，往往采用軸對稱模型，守恒形式可以表示如下：

式中：U為時間守恒矢量，F、G分別為徑向和軸向的矢通量，H為軸對稱變換矢量，S為源項。具體可表示如下：

提供輔助方程來封閉基本方程。膛內點傳火過程源項由點火管提供，點火管為一維線源，從破膜后點火燃氣流入膛內過程視為源項當地釋放。源項與膛內燃燒過程進行耦合計算，具體一維點火管控制方程及流出流量公式可參考文獻[5]，其他如狀態方程、燃燒規律、相間阻力、相間傳熱、顆粒間應力、顆粒表面溫度等輔助方程等參見文獻[2]。

2 數值計算方法

基于有限差分法，對上述模型在時間項上采用三階TVD型龍格庫塔法進行時間推進，空間項采用五階WENO格式進行空間離散。

2.1 五階WENO格式

WENO格式基本思路是將基架點分為多個組，每個組獨立計算某一點的導數逼近，得到多個差分，由于流場中間斷不可能處處存在，所以根據每個模板的光滑程度，設定權重，如果光滑度越高，則權重越大，存在間斷則權重趨近零。

對于向量形式的守恒方程：

采用全局Lax-Friedrichs分裂對空間項F、G進行流通矢量分裂：

利用分裂后的矢通量構造矢通量空間導數項的差分逼近式，以F項為例：

采用Jiang 等[14]的五階精度WENO格式構造矢通量導數逼近式的通量形式為（以正特征值為例）：

hj+1/2為正通量，由3個通量加權平均而得：

加權系數為（k=1,2,3）：

式中：C1=1/10，C2=6/10，C3=3/10，ε0=10?6，β=2，Ik為光滑度量因子。Ik的表達式為：

當特征值為負時，利用對稱性，將上式中所有下標“j+m”換成“j?m”即可。對G項的離散同理。

2.2 時間方向離散

空間導數離散化后的半離散方程為：

式中：Q n為離散后空間導數項，R n為源項。使用三階TVD龍格庫塔方法進行離散：

?t為時間步長，可以根據下式確定：

式中：C為柯朗數，dx為計算網格尺度；λ*為全局最大特征值，由上述算法可計算得出。

2.3 初始條件和邊界條件

膛內裝填條件采用如圖1所示的裝填結構，并結合實驗裝填數據和環境參數，確定模擬過程中的初始條件。膛底和出口處邊界條件采用鏡面反射方法，計算至膜片破裂。由于WENO格式采用多個基架點進行計算，在邊界外設置3個虛擬點可將邊界點轉化為內點參與計算，還可以保持邊界點的計算精度。中心軸線邊界滿足軸對稱關系式。

3 實驗驗證與數值結果分析

3.1 假藥粒實驗

為了驗證所建立模型和數值方法的可靠性，對假藥粒情況進行實驗驗證。采用了針對105 mm 口徑火炮藥室縮比的可視化實驗裝置，如圖2所示。裝置一端端蓋內側連接點火管，外側可以裝底火，用以激發點火點燃點火管內的點火藥；另一端端蓋為中通結構，用以鎖緊閉壓膜片。裝置中間裝有耐高壓玻璃管，管內裝有定量的假藥粒，通過裝置中部預留的矩形窗口可以有效觀測裝置內部的實驗現象。裝置上方設置3 個傳感器接頭，分別對應藥室的膛底、中部和彈底（膜片）處，可以通過壓力傳感器測得這幾處的壓力。

整體實驗裝置如圖2所示，可視化窗口為400 mm×120 mm，高速攝像系統可通過此窗口拍攝實驗過程中火焰發展傳播過程。傳感器距膛底分別為40、247.5和455 mm，可通過數據采集系統獲取各位置壓力。整體實驗裝填點火藥為2#小粒黑火藥，質量為30 g；假藥粒為PA 尼龍，質量為3 kg；膜片為鋁合金，厚度為0.5 mm。

圖2 可視化實驗裝置Fig.2 Visual experimental system

3.2 數值計算與實驗數據比對

以上述實驗裝置進行可視化模擬實驗的內彈道過程二維兩相流動為研究對象，其具體結構可參見圖1中的簡化物理模型，計算參量設置為：點火藥質量，30.0 g；點火藥密度，1650 kg/m3；點火藥燃速指數，0.98；點火藥燃速系數，3.68 mm·s?1·MPa?0.98；火藥力，300 J/g；點火藥平均直徑，1 mm；假藥粒密度，1100 kg/m3；假藥粒質量，3.0 kg；藥室長度，495 mm；藥室直徑，110 mm；點火管長度，25 mm；點火管直徑，19 mm；氣體絕熱指數，1.3；網格長度，5 mm；氣體余容，1 cm3/g。

為了簡化計算，忽略點火管的二維效應，即在傳火孔破膜后將點火管作為線性點火源加入膛內進行兩相耦合計算，彈底與膛底為垂直于軸向的固壁平面，采用反射壁面條件。通過壓力傳感器對膛內各位置采集的壓力與數值模擬獲取壓力進行比對，由于模擬計算只到破膜后截止，因此只取破膜前的壓力進行對比分析。

圖3所示為中間測點壓力隨時間變化的實驗結果與模擬結果的對比曲線?？梢钥闯觯M結果與實驗結果符合較好，全局平均誤差為5.35%，最高峰值處（即2.5 ms后區域）誤差僅為1.42%，說明所建立的模型和采用的數值方法合理有效。

圖3 計算壓力曲線與實驗對比圖Fig.3 Comparison of pressure-time curves between numerical simulation and experiment

整個膛內點傳火過程從點火管火藥點燃開始，高溫高壓燃氣經由點火管小孔流入膛內，點燃主裝藥。燃氣在點火管燃燒破膜以及在膛內流動需要一定的時間，而壓力傳感器測點在管壁處，所以初始階段壓力基本沒有變化，從點火開始到1 ms時間段內，均保持在一個大氣壓左右；在隨后的1～2 ms時間段內，壓力快速上升，這是由于假藥粒裝填幾乎占滿整個藥室，藥室透氣性很差，藥室自由空間很小，膛內壓力上升也較快；壓力會在到達峰值時，會產生波動現象，這種波動現象主要由于假藥粒裝填密度過高，裝藥床透氣性差，且假藥粒不會燃燒導致燃氣流動性變差，局部壓力會產生振蕩。

3.3 膛內兩相流動過程分析

通過高速攝像拍攝可視化實驗裝置內火藥燃氣流動過程，與計算得到的壓力場進行對比，可以清晰反映火藥燃氣在膛內的發展狀況。由于可視化裝置只能觀察到膛內部分區域火焰傳播過程，而火焰反映了高溫燃氣的流動過程，進而反映膛內的壓力場的變化過程，因此對計算得到的壓力云圖進行相對應位置截取對比，如圖4所示。

圖4 數值模擬壓力云圖與高速攝影圖像的對比Fig.4 Comparison between pressure clouds by numerical simulation and images by high-speed photography

整個點傳火可分為3 個階段；初始破孔階段，在壓力云圖中可以看出隨著點火管依次破孔，小孔破孔位置出現不同程度壓力梯度，而在高速攝影中，由于點火管被假藥粒包圍，很難直接觀測到初始階段的破孔現象；隨著點火管高溫燃氣不斷向膛內噴射，高溫燃氣在膛內向膛壁和藥室兩端傳播，在壓力云圖中可以觀察到其壓力梯度沿徑向和軸向都有明顯的變化，相對應的高速攝影圖像中，火焰傳播也呈現與壓力云圖同樣的特征；到達2.5 ms時，點火管對應的膛內區域燃氣流動在徑向上基本一致，即沿徑向的二維效應基本消失，而在遠離點火管的區域仍存在一定的徑向效應，從點火管開孔處向膛底和彈底形成明顯的軸向壓力梯度。

為有效呈現膛內壓力隨時間的發展變化，在徑向上取r=20 mm 處觀察壓力沿軸向的變化，如圖5所示。圖5(a)展示點火破膜初期階段（0.5～1.2 ms）膛內壓力沿軸向的變化，這一階段膛內徑向效應明顯，壓力從破孔區域不斷上升，且隨著點火管依次破孔，壓力出現了相對應的峰值，并逐漸向周圍發展。隨著點火管燃氣不斷向膛內噴射，點火管對應的膛內區域壓力基本趨于穩定，如圖5(b)所示：2.5 ms后中心位置壓力約18 MPa，徑向壓力梯度基本消失，此時在軸向壓力梯度作用下，燃氣向藥室兩端傳播，膛底和彈底的壓力不斷升高，沿軸向壓力變化區域也不斷增大。

圖5 r=20 mm 處沿軸向方向的壓力分布Fig.5 Pressure distributions along the axis at r=20 mm

圖6和圖7描述了r=20 mm 處膛內氣相的徑向速度和軸向速度沿藥室軸向的變化規律。從圖6可以看出，點火管破孔初始階段，燃氣不斷從點火管內流向膛內，此時點火管附近區域燃氣的徑向速度要遠大于其軸向速度，沿徑向傳播也比較快。初始破孔階段，由于膛內壓力遠低于管內壓力，氣體噴射出來速度極大（由于壓力差導致由點火管流入膛內氣相速度很大），隨著燃氣不斷流入，徑向壓力差變小，徑向的燃氣速度也逐漸降低。

圖6 0.5～1.0 ms 不同時刻r=20 mm 處沿軸向方向的氣相速度分布Fig.6 Gas phase velocity distributions along theaxis at r=20 mm and different times from 0.5 msto 1.0 ms

圖7 1.5～3.2 ms不同時刻r=20 mm 處氣相速度沿軸向的分布Fig.7 Gas phase velocity distributions along theaxis at r=20 mm and different times from 1.5 msto 3.2 ms

在點火1.5 ms后，大量燃氣由點火管破孔處向兩端流動，并逐漸堆積形成高壓區，在1.5～2.0 ms內，燃氣逐漸向膛底和彈底運動，氣相軸向速度漸漸升高，如圖7(a)所示；待燃氣運動到膛底或彈底時遇到阻力后又氣相軸向速度慢慢降低至反向運動，如圖7(b)所示，由于高密度裝填，在彈底形成了藥粒壅塞，所以氣相速度在接近彈底處已開始有反向的波動；隨著燃氣軸向運動加劇，燃氣的徑向速度在點火管對應區域基本為零，在2.0 ms前燃氣向膛底和彈底方向傳播，由于燃氣壓力較低，初始傳播時受顆粒擾動較大，徑向速度波動劇烈，如圖7(c)所示；2.2 ms后，由于燃氣充盈膛底和中部區域，燃氣徑向速度除在膛底有小范圍波動外，其他區域燃氣徑向速度降為零，燃氣的徑向擾動逐漸向彈底傳播，如圖7(d)所示。

火藥顆粒運動對初始階段膛內點傳火研究很重要，為進一步研究藥床中藥粒運動規律，圖8～9展示了r=20 mm 膛內固相速度沿藥室軸向發展變化。從圖8可以看出，破孔初始階段火藥顆粒軸向和徑向速度非常低，最高速度不超過0.5 m/s，說明膛內的火藥顆粒基本保持不動，呈現一種多孔介質狀態。從實驗實際拍攝的高速影像來看，固相顆粒確實未發生大范圍運動。

圖8 0.5～1.2 ms不同時刻r=20 mm 處固相速度沿軸向的分布Fig.8 Solid phase velocity distributionsalong the axis at r=20 mm and different times from 0.5 ms to 1.2 ms

隨著火藥燃氣在藥床內不斷流動，勢必會帶動火藥顆粒運動。在高密度裝填的藥床中，火藥顆粒可運動的范圍十分有限，且模擬和實驗采用的假藥粒不會燃燒且變形較小，所以顆粒間相互作用力十分大，導致藥床內的固相速度非常低。如圖9(a)所示，火藥顆粒徑向速度隨時間發展逐漸為零，從破孔位置沿軸向至彈底方向有一定的徑向速度，至膛底處徑向速度也逐漸趨向于零，與徑向氣相速度相對應，說明氣相速度也是影響顆粒運動的重要原因。圖9(b)中火藥顆粒軸向速度也驗證了這一點，由于燃氣軸向壓力梯度大，氣相軸向速度高，所以導致火藥顆粒向膛底和彈底兩側運動速度較高。在真實裝藥情況下，會造成膛底火藥粒擠壓破碎，形成小范圍高壓區域，繼而推動藥床向彈底運動。

圖9 2.5～3.2 ms不同時刻r=20 mm 處固相速度沿軸向的分布Fig.9 Solid phase velocity distributionsalong the axis at r=20 mm and different times from 2.5 ms to 3.2 ms

4 結 論

以火炮膛內初期點傳火性能研究為背景，建立了軸對稱二維兩相流數理模型，對假藥粒進行相應合理化假設，在空間上采用五階WENO格式進行空間離散，時間上采用三階TVD型龍格庫塔法進行時間推進，并進行實驗與數值模擬對比研究。數值分析了點火初期膛內壓力的發展及氣固兩相速度，具體結論如下。

（1）可視化實驗裝置與高精度數值方法，有效模擬了初始階段點火燃氣在藥床中的傳播規律，通過與可視化實驗平臺測試數據和現象進行對比，驗證了數理模型和數值方法的有效性。

（2）點火初期膛內徑向效應明顯，隨著傳火不斷進行，點火管對應膛內區域徑向效應逐漸消失，燃氣的軸向傳播開始占據主導地位，膛內壓力呈現徑向一致、軸向梯度分布的特征。

（3）通過數值模擬與實驗觀測可以得出，在高密度裝填情況下，初始階段藥床藥?；静粫l生運動，最高速度為0.5 m/s，后續固相速度的發展與氣相速度有緊密的關系；氣相速度在初始階段沿徑向傳播較快，后續在軸向壓力梯度作用下，軸向速度開始占據主導，在膛底和中部區域氣相速度降低為零，徑向速度擾動也逐漸向彈底傳播，在到達彈底之前由于藥粒壅塞現象導致阻力變大，氣相速度會出現振蕩。