李梓鈺, 余永剛
(南京理工大學 能源與動力工程學院, 江蘇 南京 210094)
模塊裝藥點傳火過程是影響火炮內(nèi)彈道性能的重要環(huán)節(jié),而點傳火過程中模塊藥盒破裂、火藥顆粒在藥室內(nèi)的飛散運動則可能影響壓力波的產(chǎn)生和增長。這種壓力波現(xiàn)象會對內(nèi)彈道性能的穩(wěn)定性產(chǎn)生影響,嚴重時可能引起膛內(nèi)局部壓力異常的升高從而造成彈丸早炸、膛炸等安全事故。
關于模塊裝藥膛內(nèi)壓力波的研究,許多學者是將火藥顆粒看作連續(xù)介質(zhì),建立歐拉-歐拉模型對壓力波的影響因素進行試驗和模擬仿真分析。如Ruth等[1]對155 mm XM216 模塊裝藥進行裝藥安全性研究,探究常溫和極端溫度條件下內(nèi)彈道循環(huán)過程中的點火、藥床移動性和壓力波的發(fā)展,發(fā)現(xiàn)模塊盒體破裂前燃燒室內(nèi)出現(xiàn)快速增壓,點火循環(huán)后期壓力波動較小。陸兵中等[2-3]基于短管試驗炮,考慮模塊運動性和初始裝填分布的影響,對點火及擴散、壓力波動的發(fā)展過程進行模擬和預測,發(fā)現(xiàn)火焰擴散期間膛內(nèi)壓力差較小,使得裝藥床移動速度慢、移動距離短;后又考慮藥盒間的隔倉效應及環(huán)形間隙的作用等,對膛內(nèi)壓力波動進行精確模擬,發(fā)現(xiàn)相同裝藥量下,模塊藥盒材料的含能密度越低,膛內(nèi)壓差越小。王育維等[4-5]則從裝藥長度及模塊盒體的能量特性兩個角度對影響小號模塊裝藥壓力波的因素進行分析,依據(jù)兩相流內(nèi)彈道理論,發(fā)現(xiàn)隨裝藥長度的增加,壓力越穩(wěn)定,壓力波越小;模塊藥盒的火藥力越大,則產(chǎn)生的壓力波也越大。除此之外,他們還建立了雙一維多相流內(nèi)彈道模型,考察不同溫度下模塊藥盒能量大小對壓力波的影響。但隨著計算機技術和相關理論的飛速發(fā)展,通過數(shù)值模擬從顆粒層次對壓力波的影響因素進行研究成為可能。
在有關能源、化工、制藥等工業(yè)生產(chǎn)以及農(nóng)業(yè)、食品等領域中,為研究顆粒的動力學特性及其運動機制,國內(nèi)外學者建立了多種數(shù)值模型對顆粒系統(tǒng)的復雜行為進行模擬,以探究試驗難以獲得的深層次、多尺度運動機理。其中,對于非球形顆粒采用DEM模型結(jié)合計算流體力學(CFD)的數(shù)值模型已成功應用于各種復雜流-固系統(tǒng),如噴動床[6-10]、氣力輸送[11]、鼓泡流化床[12-18]等氣-固兩相流動過程。目前對于模塊裝藥藥粒散布的試驗與仿真研究,是將藥粒簡化為球形顆粒進行模擬仿真[19-20],然而其實際形狀更接近于柱狀或桿狀。并且當前研究普遍認為,顆粒的形狀對顆粒行為具有較大影響。因此,在研究中應用球形顆粒假設得到的顆粒動力學特性極有可能與真實的非球形顆粒系統(tǒng)不同。
針對上述情況,本文以單模塊裝藥火藥床為研究對象,采用離散單元法模擬柱狀火藥顆粒的散布過程,對模塊藥盒破裂、柱狀火藥顆粒在藥室內(nèi)的飛散運動和堆積形態(tài)進行數(shù)值模擬。根據(jù)模擬結(jié)果,對單模塊裝藥藥粒散布特性進行深入分析,以期從顆粒層次上對單模塊裝藥點傳火過程中火藥顆粒運動規(guī)律有更深入的理解和認識,進而為后續(xù)研究多模塊裝藥點傳火過程中藥粒散布形態(tài)提供理論基礎和手段。
為更好地探究模塊裝藥點傳火過程中火藥顆粒運動規(guī)律及其堆積形態(tài),建立了該過程的數(shù)理模型,并為試驗設計了可視化點傳火平臺,其結(jié)構圖如圖1所示,它由135 mm 口徑的有機玻璃身管、擊發(fā)點火裝置、控壓膜片等組成。身管由金屬套筒固定,其上方開有測壓孔,該孔還用于定位;藥室右側(cè)為半密閉端,裝有控壓膜片,使藥室壓力控制在安全極限值內(nèi)。單模塊裝藥點傳火過程可簡單概括為:首先擊發(fā)點燃中心傳火管內(nèi)點火藥包,點火藥包燃燒產(chǎn)生高溫高壓燃氣,隨著壓力不斷升高,達到傳火孔破孔壓力,火藥燃氣和灼熱固體微粒通過對流和輻射換熱點燃發(fā)射藥,模塊藥盒發(fā)生破裂,盒內(nèi)未燃完的發(fā)射藥會在藥室內(nèi)飛散運動,最終形成非均勻堆積散布。
圖1 可視化點傳火試驗裝置結(jié)構圖Fig.1 Schematic diagram of a visual experimental platform for the ignition and flame-spreading process
為充分考慮氣相和固相間的相互作用以及固相顆粒間的碰撞,本文采用三維非穩(wěn)態(tài)氣-固兩相流模型,并做如下假設:
1)用球元法構建的類圓柱形顆粒代替標準圓柱形顆粒,假設所有模擬藥粒大小和形狀完全相同;
2)模塊藥盒的破口壓力約為10 MPa,且破口時藥盒右端蓋端瞬時全部破開;
3)模塊藥盒與藥室等直徑處理,不考慮藥盒在藥室內(nèi)的運動;
4)采用Realizablek-ε(k為湍流動能,ε為湍流耗散率)湍流模型描述點火燃氣射流的湍流流動。
氣相控制方程包括連續(xù)性方程、動量守恒方程和能量方程,其表達式如下:
(1)
(2)
式中:ρg和εg分別表示氣體密度和氣體體積分數(shù);ug為氣體速度;p為氣體壓力;τg為氣體黏度;g為重力加速度;S為氣-固兩相間的動量交換量。
本文采用球元法來描述柱狀火藥顆粒,如圖2所示,圖2(a)為標準柱狀顆粒,圖2(b)為3個球元構成的類柱狀顆粒。由球元法構建的非球形顆粒,其受力和運動模型得到了簡化,將顆粒與顆粒間、顆粒與壁面間的碰撞力轉(zhuǎn)移到球元上計算,從而可以采用成熟的球形顆粒計算公式。但缺點便是引入了人為粗糙度,增加了顆粒摩擦力,并且隨著構建顆粒所需球元數(shù)目的增多,顆粒受力計算誤差增加,計算時間增加。根據(jù)Markauskas 等提出的計算人為粗糙度的方法,任冰[21]對不同方法構造的正方體顆粒的人為粗糙度進行計算并通過顆粒自由下落試驗和堆積試驗進行驗證,發(fā)現(xiàn)隨著構建球元數(shù)量的增加(用一種球元構建的情況下),所構建的顆粒與真實顆粒間的誤差就越小??紤]試驗所用模擬發(fā)射藥粒的形狀大小,本文采用3個球元構建柱狀顆粒。
圖2 柱狀火藥顆粒模型示意圖Fig.2 Diagram of the cylindrical propellent pellets
由球元法所構建的柱狀火藥顆粒的碰撞接觸形式可分為六類:球元-球元接觸,球元-顆粒接觸,顆粒-顆粒交叉接觸,顆粒-顆粒平行接觸,球元-壁面接觸,顆粒-壁面接觸,具體形式如圖3所示。
圖3 柱狀顆粒的碰撞機制Fig.3 Collision mechanism of cylindrical particles
圖4為柱狀顆粒碰撞后的受力示意圖,顆粒主要受到重力、碰撞力以及氣-固間作用力,其運動方式主要是平動和轉(zhuǎn)動兩種方式,由牛頓第二定律描述如下:
圖4 柱狀顆粒受力示意圖Fig.4 Diagram of forces and motion of the cylinder pellet
(3)
(4)
式中:mp為柱狀模擬發(fā)射藥粒的質(zhì)量;vp為顆粒的運動速度;Gp為顆粒的重力;Fc為碰撞力;Ff為氣-固間作用力;Ip為柱狀顆粒的轉(zhuǎn)動慣量;ωp為顆粒的角速度;Tp為柱狀顆粒轉(zhuǎn)動的合力矩。每個柱狀顆粒所受到的碰撞力是組成它的每個球元的碰撞力的總和,即
(5)
式中:Fce為組成柱狀顆粒某球元j所受的碰撞合力;Ne為組成柱狀顆粒的球元數(shù)量,本文取3。構成同一顆粒的球元間無相互作用力,顆粒-顆粒、顆粒-壁面間均采用Hertz-Mindlin無滑動接觸模型[22]。
氣相和固相間的曳力作用作為相間動量傳遞的主要方式將氣-固運動耦合在一起,本文選取的是Free-stream曳力模型[22],即
Fd=0.5CdρgA|vp-g|vp-g
(6)
式中:A為顆粒的投影面積;vp-g為氣體與顆粒的相對速度;Cd為曳力系數(shù),
(7)
Re為雷諾數(shù)。
本文中分別通過Fluent軟件和EDEM軟件進行CFD和DEM求解,CFD對流相信息進行跟蹤,DEM對顆粒相信息進行跟蹤。CFD首先由初始條件計算流相參數(shù),然后將仿真控制權交給DEM;DEM根據(jù)流場信息和顆粒碰撞求解并更新顆粒的速度和位置等信息,并將此信息和控制權傳回CFD;CFD根據(jù)返回的顆粒信息和流場信息計算該時間步長內(nèi)的氣-固相互作用力,計算當前網(wǎng)格的動量、能量源相等,進而完成一次耦合計算。
為便于記錄和觀察單模塊裝藥點傳火過程中藥盒破裂、藥粒飛散運動,建立可視化點傳火試驗平臺。如圖5所示,試驗裝置由半密閉爆發(fā)器點傳火平臺、壓力傳感器、數(shù)據(jù)采集及處理設備等組成。試驗過程中壓電傳感器獲得壓力信號,經(jīng)電荷放大器將電信號放大,再輸入瞬態(tài)記錄儀進行轉(zhuǎn)換,最后通過計算系統(tǒng)進行數(shù)據(jù)的采集和處理,并且利用高速攝像機對單模塊裝藥點傳火的全過程進行記錄。計算模型結(jié)構示意圖及柱狀顆粒平面模型如圖6所示,模塊藥盒位于半密閉爆發(fā)器身管的左側(cè),與底火端的距離l0為40 mm,藥盒右側(cè)的藥室空間則均為空氣。模塊藥盒與藥室身管壁的實際間隙為5 mm,但為了簡化計算模型,模擬仿真時該間距忽略不計,且不考慮模塊藥盒與中心傳火管的壁厚。為確保試驗的安全性,藥盒內(nèi)僅裝有少量的真火藥粒(約占藥盒體積的2%),絕大部分為圓柱形的模擬發(fā)射藥粒,其尺寸為φ8 mm×13 mm。試驗時,利用機械能量點燃模塊裝藥中心傳火管內(nèi)附著的點火藥包,點火藥燃燒產(chǎn)生高溫高壓燃氣使得靠近傳火管的發(fā)射藥首先被點燃,而后混合燃氣不斷擴散,使得模塊藥盒內(nèi)真火藥顆粒被點燃。過程中由于模塊藥盒內(nèi)壓力不斷上升,與鄰接藥室自由空間形成較大壓力梯度,達到破口壓力后,燃氣和假火藥會以較大速度向藥室右側(cè)密閉端運動。求解時使用壓力求解器,壓力-速度耦合采用SIMPLE 模式,以實現(xiàn)動量、能量等參數(shù)的傳遞。由于氣相與顆粒相耦合問題中,CFD 網(wǎng)格尺寸必須大于顆粒尺寸,故取網(wǎng)格尺寸為 15 mm,共計4 445個網(wǎng)格單元,并經(jīng)過網(wǎng)格無關性驗證。
圖5 試驗裝置示意圖Fig.5 Schematic diagram of experimental devices
圖6 計算模型結(jié)構示意圖和柱狀顆粒模型示意圖Fig.6 Geometry of the calculation model and the cylindrical pellet model
圖7為高速攝像機拍攝的藥室內(nèi)藥粒在1 000 ms時的堆積形態(tài),可看出藥粒呈坡狀集中分布于 350~475 mm區(qū)域[20]。為定量描述火藥堆積形態(tài),在X分別取350 mm、375 mm、400 mm、425 mm、450 mm、475 mm,Y=0 mm處分別測量藥粒的高度。
圖7 1 000 ms時刻模塊裝藥點傳火試驗序列圖[20]Fig.7 Sequence diagram of the modular charge’s ignition and flame-spreading experiment at 1 000 ms[20]
圖8為數(shù)值計算得到的終態(tài)藥粒分布的正視圖,與試驗所取特征點相對應,對終態(tài)藥粒坡狀堆積表面取6個點進行測量。將試驗值與數(shù)值計算值相比較,如圖9所示,發(fā)現(xiàn)數(shù)值模擬與試驗測量值的平均誤差約為4.9%,遠小于試驗估計的測量誤差15%,并且與采用等體積球形顆粒的模擬仿真相比,一定程度上縮小了誤差值,進一步表明顆粒的形狀對顆粒行為具有很大影響,建立的類圓柱形顆粒離散單元模型對模擬單模塊裝藥點傳火過程中的藥粒分布形態(tài)具有較高的精度。
圖8 終態(tài)藥粒分布正視圖(l0=40 mm)Fig.8 Front view of propellant pellets’ final distribution (l0=40 mm)
圖9 藥粒堆積高度的試驗值與數(shù)值模擬計算值Fig.9 Exprimental and numerical values of propellant pellets’ stacking height
圖10為0.23~10 ms時藥室Z=0 mm截面的氣相壓力云圖。由圖10(a)可知,在藥盒破口前,中心傳火管內(nèi)附著的點火藥包被點燃,經(jīng)傳火孔向藥盒內(nèi)輸送高溫高壓燃氣,模塊藥盒壓力快速上升,約0.23 ms即由1大氣壓強上升至10 MPa;同時火藥燃氣也通過傳火管端口向藥室輸入,但由于右側(cè)藥室相比模塊藥盒體積更大,故與模塊藥盒內(nèi)壓力相比壓力上升速度較緩。當藥盒內(nèi)壓力上升至約 10 MPa 時,假設藥盒右端蓋全部破開,由于此時模塊藥盒內(nèi)顆粒密集,相間阻力大,火藥燃氣難以滲透,裝藥和藥室自由空間鄰接處形成較大壓力梯度,如圖10(b)所示,使得燃氣和藥?;旌衔飶哪K藥盒右端蓋泄出。如圖10(c)~圖10(f)所示,當點火藥和少量發(fā)射藥粒燃燒完全后,藥室內(nèi)壓差逐漸減小,最大壓力不再升高。
為更好地理解和研究單模塊裝藥點傳火過程中的藥粒運動特性,圖11給出了仿真過程中0.23~600 ms 內(nèi)藥粒速度分布及飛散形態(tài)。如圖11(a)所示,點火藥燃燒產(chǎn)生的高壓燃氣經(jīng)由傳火孔向模塊藥盒內(nèi)輸入,在燃氣射流作用下盒內(nèi)顆粒獲得速度,且越靠近傳火管壁顆粒速度越大。0.23 ms 時藥盒右端蓋瞬時完全破裂,此時端蓋兩側(cè)最大壓差約為 9.6 MPa,大部分藥粒和氣體混合在一起加速向藥室右側(cè)運動。1.5 ms 時顆粒加速到最大值74.8 m/s,并且顆粒由于速度不同在藥室內(nèi)沿身管軸向呈環(huán)形分層分布。如圖11(c)所示,6 ms時加速最快即速度分層最右側(cè)的部分藥粒最先與右側(cè)密閉端固壁產(chǎn)生碰撞,速度迅速減小并且部分顆粒獲得反向速度。此后,由于點火藥和發(fā)射藥的混合高壓燃氣不再生成,藥室內(nèi)的壓力梯度逐漸減小,藥粒向右的運動趨勢不變,但速度在碰撞過程中減小。如圖11(h)所示,火藥顆粒600 ms時最大速度已小于0.5 m/s,且運動顆粒僅占2%左右,藥粒整體趨于靜止,故將該時刻近似作為藥粒在藥室內(nèi)飛散運動的終止時間,此時的藥粒分布形態(tài)視為最終形態(tài)。
圖11 不同時刻藥室內(nèi)的藥粒速度分布及飛散形態(tài)Fig.11 Velocity and distribution of propellant pellets in the chamber at different moments
為研究模塊裝填的初始位置對藥粒分布形態(tài)的影響,如圖12所示,將模塊初始裝填位置分別設置為距底火端10 mm、40 mm、70 mm和100 mm,數(shù)值模擬以上4種工況下單模塊裝藥藥盒破裂后藥粒飛散運動。由圖12可以看出,4種工況下藥室內(nèi)終態(tài)藥粒從底火端向密閉端依次呈水平堆積和坡狀堆積形態(tài)。為更好地對藥粒堆積形態(tài)進行量化比較,引入了水平堆積軸向長度L和坡狀堆積坡度角α兩個特征參數(shù),L表示沿藥室身管軸向水平堆積起點到坡狀堆積起點的距離,α則表示坡狀堆積斜坡表面與水平面的夾角。分別對以上4種工況下終態(tài)藥粒堆積表面進行L和α的測量,詳細數(shù)據(jù)如表1所示。
圖12 不同裝填位置下的終態(tài)藥粒堆積分布Fig.12 Distribution of pellets under different working conditions
表1 藥粒堆積特征參數(shù)
通過比較這兩種特征參數(shù)在4種工況下的變化,可以看出隨著模塊藥盒與藥室底火端的距離增大,水平堆積軸向長度減小,坡狀堆積變平緩,坡度角減小。由圖13可以看出,水平堆積軸向長度L隨著模塊藥盒裝填位置與底火端距離l0的增大線性縮短,其擬合公式為
圖13 水平堆積軸向長度與模塊藥盒裝填軸向位置的關系Fig.13 Relationship between the axial length of the horizontal stacking of propellant pellets and the horizontal loading position of the modular cartridge
L=-1.011 3l0+303.473 3 (R2=0.995 61)
上述結(jié)果與陳安等[20]基于球形顆粒建立的模型所獲得的顆粒堆積變化趨勢一致,但與之相比,采用柱狀顆粒模型仿真獲得的藥粒堆積特性參數(shù)變化更小。之所以出現(xiàn)上述情況,是因為模塊藥盒右端蓋破裂瞬間,顆粒受到氣相壓強的作用,因為等體積球形顆粒垂直于壓力梯度的投影大于大部分柱形顆粒(與壓力梯度方向垂直或近似垂直的較少顆粒除外),球形顆粒獲得更大的壓力,故兩種顆粒在相同質(zhì)量下,球形顆粒將獲得更大的加速度,柱狀顆粒與之相比到達藥室密閉端的速度更小,由圖11(b)和圖11(c)可知,顆粒最大速度約為74.73 m/s,到達右側(cè)密閉端的速度約為64 m/s,而球形顆粒的最大速度可達93 m/s[20]。故在后面較小壓強梯度下加速的顆粒距離藥室右側(cè)較遠距離時速度就不再增加,使得終態(tài)顆粒分布較均勻,坡狀堆積的坡度角較小。而隨著模塊藥盒的裝填位置右移,水平堆積軸向長度縮短,主要原因是模塊藥盒破裂時僅右端蓋全部破開,顆粒在氣相壓力的作用下向右運動,水平堆積起始位置向右移動;顆粒堆積坡度減小則是因為隨著模塊藥盒與藥室底火端的距離增大,顆粒到達藥室右壁面的距離減小、速度增大,無論是顆粒-壁面還是顆粒-顆粒發(fā)生碰撞后,將獲得更大的反向速度,藥粒將反向運動更遠的距離后落下,散布更均勻,使得坡狀堆積角變小,但由于距離變化較小對顆粒速影響比較小,這種變化也不是很明顯。
本文針對單模塊裝藥點傳火過程中藥盒破裂、藥粒在藥室中的飛散和堆積問題,建立了三維非穩(wěn)態(tài)氣-固兩相流模型并進行了數(shù)值模擬。得出如下主要結(jié)論:
1)藥粒在藥室中的終態(tài)分布呈現(xiàn)非均勻狀態(tài),從底火端到密閉端依次呈水平堆積和坡狀堆積,與試驗結(jié)果基本吻合,驗證了該模型的合理性。
2)采用球元法建立了柱狀火藥顆粒的離散元模型,與基于等球形顆粒建立的模型相比,模擬仿真所得的結(jié)果更加精確,更接近實測結(jié)果,表明顆粒形狀對顆粒飛散運動行為有較大影響。
3)單模塊裝藥初始裝填位置與藥室底火端距離由10 mm增大到100 mm時,水平堆積軸向長度呈線性縮短,由293.3 mm縮短為203.8 mm,同時坡狀堆積坡度角由22.01°降為19.82°,堆積逐漸變平緩。