基于非球形顆粒離散單元法分析單模塊裝藥藥粒散布特性

李梓鈺, 余永剛

(南京理工大學 能源與動力工程學院, 江蘇 南京 210094)

0 引言

模塊裝藥點傳火過程是影響火炮內彈道性能的重要環節,而點傳火過程中模塊藥盒破裂、火藥顆粒在藥室內的飛散運動則可能影響壓力波的產生和增長。這種壓力波現象會對內彈道性能的穩定性產生影響,嚴重時可能引起膛內局部壓力異常的升高從而造成彈丸早炸、膛炸等安全事故。

關于模塊裝藥膛內壓力波的研究,許多學者是將火藥顆?？醋鬟B續介質,建立歐拉-歐拉模型對壓力波的影響因素進行試驗和模擬仿真分析。如Ruth等[1]對155 mm XM216 模塊裝藥進行裝藥安全性研究,探究常溫和極端溫度條件下內彈道循環過程中的點火、藥床移動性和壓力波的發展,發現模塊盒體破裂前燃燒室內出現快速增壓,點火循環后期壓力波動較小。陸兵中等[2-3]基于短管試驗炮,考慮模塊運動性和初始裝填分布的影響,對點火及擴散、壓力波動的發展過程進行模擬和預測,發現火焰擴散期間膛內壓力差較小,使得裝藥床移動速度慢、移動距離短;后又考慮藥盒間的隔倉效應及環形間隙的作用等,對膛內壓力波動進行精確模擬,發現相同裝藥量下,模塊藥盒材料的含能密度越低,膛內壓差越小。王育維等[4-5]則從裝藥長度及模塊盒體的能量特性兩個角度對影響小號模塊裝藥壓力波的因素進行分析,依據兩相流內彈道理論,發現隨裝藥長度的增加,壓力越穩定,壓力波越小;模塊藥盒的火藥力越大,則產生的壓力波也越大。除此之外,他們還建立了雙一維多相流內彈道模型,考察不同溫度下模塊藥盒能量大小對壓力波的影響。但隨著計算機技術和相關理論的飛速發展,通過數值模擬從顆粒層次對壓力波的影響因素進行研究成為可能。

在有關能源、化工、制藥等工業生產以及農業、食品等領域中,為研究顆粒的動力學特性及其運動機制,國內外學者建立了多種數值模型對顆粒系統的復雜行為進行模擬,以探究試驗難以獲得的深層次、多尺度運動機理。其中,對于非球形顆粒采用DEM模型結合計算流體力學(CFD)的數值模型已成功應用于各種復雜流-固系統,如噴動床[6-10]、氣力輸送[11]、鼓泡流化床[12-18]等氣-固兩相流動過程。目前對于模塊裝藥藥粒散布的試驗與仿真研究,是將藥粒簡化為球形顆粒進行模擬仿真[19-20],然而其實際形狀更接近于柱狀或桿狀。并且當前研究普遍認為,顆粒的形狀對顆粒行為具有較大影響。因此,在研究中應用球形顆粒假設得到的顆粒動力學特性極有可能與真實的非球形顆粒系統不同。

針對上述情況,本文以單模塊裝藥火藥床為研究對象,采用離散單元法模擬柱狀火藥顆粒的散布過程,對模塊藥盒破裂、柱狀火藥顆粒在藥室內的飛散運動和堆積形態進行數值模擬。根據模擬結果,對單模塊裝藥藥粒散布特性進行深入分析,以期從顆粒層次上對單模塊裝藥點傳火過程中火藥顆粒運動規律有更深入的理解和認識,進而為后續研究多模塊裝藥點傳火過程中藥粒散布形態提供理論基礎和手段。

1 物理模型

為更好地探究模塊裝藥點傳火過程中火藥顆粒運動規律及其堆積形態,建立了該過程的數理模型,并為試驗設計了可視化點傳火平臺,其結構圖如圖1所示,它由135 mm 口徑的有機玻璃身管、擊發點火裝置、控壓膜片等組成。身管由金屬套筒固定,其上方開有測壓孔,該孔還用于定位;藥室右側為半密閉端,裝有控壓膜片,使藥室壓力控制在安全極限值內。單模塊裝藥點傳火過程可簡單概括為:首先擊發點燃中心傳火管內點火藥包,點火藥包燃燒產生高溫高壓燃氣,隨著壓力不斷升高,達到傳火孔破孔壓力,火藥燃氣和灼熱固體微粒通過對流和輻射換熱點燃發射藥,模塊藥盒發生破裂,盒內未燃完的發射藥會在藥室內飛散運動,最終形成非均勻堆積散布。

圖1 可視化點傳火試驗裝置結構圖Fig.1 Schematic diagram of a visual experimental platform for the ignition and flame-spreading process

為充分考慮氣相和固相間的相互作用以及固相顆粒間的碰撞,本文采用三維非穩態氣-固兩相流模型,并做如下假設:

1)用球元法構建的類圓柱形顆粒代替標準圓柱形顆粒,假設所有模擬藥粒大小和形狀完全相同;

2)模塊藥盒的破口壓力約為10 MPa,且破口時藥盒右端蓋端瞬時全部破開;

3)模塊藥盒與藥室等直徑處理,不考慮藥盒在藥室內的運動;

4)采用Realizablek-ε(k為湍流動能,ε為湍流耗散率)湍流模型描述點火燃氣射流的湍流流動。

2 數學模型

2.1 氣相控制方程

氣相控制方程包括連續性方程、動量守恒方程和能量方程,其表達式如下:

(1)

(2)

式中:ρg和εg分別表示氣體密度和氣體體積分數;ug為氣體速度;p為氣體壓力;τg為氣體黏度;g為重力加速度;S為氣-固兩相間的動量交換量。

2.2 柱狀火藥顆粒的接觸模型及運動方程

本文采用球元法來描述柱狀火藥顆粒,如圖2所示,圖2(a)為標準柱狀顆粒,圖2(b)為3個球元構成的類柱狀顆粒。由球元法構建的非球形顆粒,其受力和運動模型得到了簡化,將顆粒與顆粒間、顆粒與壁面間的碰撞力轉移到球元上計算,從而可以采用成熟的球形顆粒計算公式。但缺點便是引入了人為粗糙度,增加了顆粒摩擦力,并且隨著構建顆粒所需球元數目的增多,顆粒受力計算誤差增加,計算時間增加。根據Markauskas 等提出的計算人為粗糙度的方法,任冰[21]對不同方法構造的正方體顆粒的人為粗糙度進行計算并通過顆粒自由下落試驗和堆積試驗進行驗證,發現隨著構建球元數量的增加(用一種球元構建的情況下),所構建的顆粒與真實顆粒間的誤差就越小?？紤]試驗所用模擬發射藥粒的形狀大小,本文采用3個球元構建柱狀顆粒。

圖2 柱狀火藥顆粒模型示意圖Fig.2 Diagram of the cylindrical propellent pellets

由球元法所構建的柱狀火藥顆粒的碰撞接觸形式可分為六類:球元-球元接觸,球元-顆粒接觸,顆粒-顆粒交叉接觸,顆粒-顆粒平行接觸,球元-壁面接觸,顆粒-壁面接觸,具體形式如圖3所示。

圖3 柱狀顆粒的碰撞機制Fig.3 Collision mechanism of cylindrical particles

圖4為柱狀顆粒碰撞后的受力示意圖,顆粒主要受到重力、碰撞力以及氣-固間作用力,其運動方式主要是平動和轉動兩種方式,由牛頓第二定律描述如下:

圖4 柱狀顆粒受力示意圖Fig.4 Diagram of forces and motion of the cylinder pellet

(3)

(4)

式中:mp為柱狀模擬發射藥粒的質量;vp為顆粒的運動速度;Gp為顆粒的重力;Fc為碰撞力;Ff為氣-固間作用力;Ip為柱狀顆粒的轉動慣量;ωp為顆粒的角速度;Tp為柱狀顆粒轉動的合力矩。每個柱狀顆粒所受到的碰撞力是組成它的每個球元的碰撞力的總和,即

(5)

式中:Fce為組成柱狀顆粒某球元j所受的碰撞合力;Ne為組成柱狀顆粒的球元數量,本文取3。構成同一顆粒的球元間無相互作用力,顆粒-顆粒、顆粒-壁面間均采用Hertz-Mindlin無滑動接觸模型[22]。

2.3 曳力模型

氣相和固相間的曳力作用作為相間動量傳遞的主要方式將氣-固運動耦合在一起,本文選取的是Free-stream曳力模型[22],即

Fd=0.5CdρgA|vp-g|vp-g

(6)

式中:A為顆粒的投影面積;vp-g為氣體與顆粒的相對速度;Cd為曳力系數,

(7)

Re為雷諾數。

2.4 氣-固相間耦合

本文中分別通過Fluent軟件和EDEM軟件進行CFD和DEM求解,CFD對流相信息進行跟蹤,DEM對顆粒相信息進行跟蹤。CFD首先由初始條件計算流相參數,然后將仿真控制權交給DEM;DEM根據流場信息和顆粒碰撞求解并更新顆粒的速度和位置等信息,并將此信息和控制權傳回CFD;CFD根據返回的顆粒信息和流場信息計算該時間步長內的氣-固相互作用力,計算當前網格的動量、能量源相等,進而完成一次耦合計算。

3 結果與討論

3.1 DEM模型的驗證

為便于記錄和觀察單模塊裝藥點傳火過程中藥盒破裂、藥粒飛散運動,建立可視化點傳火試驗平臺。如圖5所示,試驗裝置由半密閉爆發器點傳火平臺、壓力傳感器、數據采集及處理設備等組成。試驗過程中壓電傳感器獲得壓力信號,經電荷放大器將電信號放大,再輸入瞬態記錄儀進行轉換,最后通過計算系統進行數據的采集和處理,并且利用高速攝像機對單模塊裝藥點傳火的全過程進行記錄。計算模型結構示意圖及柱狀顆粒平面模型如圖6所示,模塊藥盒位于半密閉爆發器身管的左側,與底火端的距離l0為40 mm,藥盒右側的藥室空間則均為空氣。模塊藥盒與藥室身管壁的實際間隙為5 mm,但為了簡化計算模型,模擬仿真時該間距忽略不計,且不考慮模塊藥盒與中心傳火管的壁厚。為確保試驗的安全性,藥盒內僅裝有少量的真火藥粒(約占藥盒體積的2%),絕大部分為圓柱形的模擬發射藥粒,其尺寸為φ8 mm×13 mm。試驗時,利用機械能量點燃模塊裝藥中心傳火管內附著的點火藥包,點火藥燃燒產生高溫高壓燃氣使得靠近傳火管的發射藥首先被點燃,而后混合燃氣不斷擴散,使得模塊藥盒內真火藥顆粒被點燃。過程中由于模塊藥盒內壓力不斷上升,與鄰接藥室自由空間形成較大壓力梯度,達到破口壓力后,燃氣和假火藥會以較大速度向藥室右側密閉端運動。求解時使用壓力求解器,壓力-速度耦合采用SIMPLE 模式,以實現動量、能量等參數的傳遞。由于氣相與顆粒相耦合問題中,CFD 網格尺寸必須大于顆粒尺寸,故取網格尺寸為 15 mm,共計4 445個網格單元,并經過網格無關性驗證。

圖5 試驗裝置示意圖Fig.5 Schematic diagram of experimental devices

圖6 計算模型結構示意圖和柱狀顆粒模型示意圖Fig.6 Geometry of the calculation model and the cylindrical pellet model

圖7為高速攝像機拍攝的藥室內藥粒在1 000 ms時的堆積形態,可看出藥粒呈坡狀集中分布于 350～475 mm區域[20]。為定量描述火藥堆積形態,在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為數值計算得到的終態藥粒分布的正視圖,與試驗所取特征點相對應,對終態藥粒坡狀堆積表面取6個點進行測量。將試驗值與數值計算值相比較,如圖9所示,發現數值模擬與試驗測量值的平均誤差約為4.9%,遠小于試驗估計的測量誤差15%,并且與采用等體積球形顆粒的模擬仿真相比,一定程度上縮小了誤差值,進一步表明顆粒的形狀對顆粒行為具有很大影響,建立的類圓柱形顆粒離散單元模型對模擬單模塊裝藥點傳火過程中的藥粒分布形態具有較高的精度。

圖8 終態藥粒分布正視圖(l0=40 mm)Fig.8 Front view of propellant pellets’ final distribution (l0=40 mm)

圖9 藥粒堆積高度的試驗值與數值模擬計算值Fig.9 Exprimental and numerical values of propellant pellets’ stacking height

3.2 仿真結果與分析

圖10為0.23～10 ms時藥室Z=0 mm截面的氣相壓力云圖。由圖10(a)可知,在藥盒破口前,中心傳火管內附著的點火藥包被點燃,經傳火孔向藥盒內輸送高溫高壓燃氣,模塊藥盒壓力快速上升,約0.23 ms即由1大氣壓強上升至10 MPa;同時火藥燃氣也通過傳火管端口向藥室輸入,但由于右側藥室相比模塊藥盒體積更大,故與模塊藥盒內壓力相比壓力上升速度較緩。當藥盒內壓力上升至約 10 MPa 時,假設藥盒右端蓋全部破開,由于此時模塊藥盒內顆粒密集,相間阻力大,火藥燃氣難以滲透,裝藥和藥室自由空間鄰接處形成較大壓力梯度,如圖10(b)所示,使得燃氣和藥粒混合物從模塊藥盒右端蓋泄出。如圖10(c)～圖10(f)所示,當點火藥和少量發射藥粒燃燒完全后,藥室內壓差逐漸減小,最大壓力不再升高。

為更好地理解和研究單模塊裝藥點傳火過程中的藥粒運動特性,圖11給出了仿真過程中0.23～600 ms 內藥粒速度分布及飛散形態。如圖11(a)所示,點火藥燃燒產生的高壓燃氣經由傳火孔向模塊藥盒內輸入,在燃氣射流作用下盒內顆粒獲得速度,且越靠近傳火管壁顆粒速度越大。0.23 ms 時藥盒右端蓋瞬時完全破裂,此時端蓋兩側最大壓差約為 9.6 MPa,大部分藥粒和氣體混合在一起加速向藥室右側運動。1.5 ms 時顆粒加速到最大值74.8 m/s,并且顆粒由于速度不同在藥室內沿身管軸向呈環形分層分布。如圖11(c)所示,6 ms時加速最快即速度分層最右側的部分藥粒最先與右側密閉端固壁產生碰撞,速度迅速減小并且部分顆粒獲得反向速度。此后,由于點火藥和發射藥的混合高壓燃氣不再生成,藥室內的壓力梯度逐漸減小,藥粒向右的運動趨勢不變,但速度在碰撞過程中減小。如圖11(h)所示,火藥顆粒600 ms時最大速度已小于0.5 m/s,且運動顆粒僅占2%左右,藥粒整體趨于靜止,故將該時刻近似作為藥粒在藥室內飛散運動的終止時間,此時的藥粒分布形態視為最終形態。

圖11 不同時刻藥室內的藥粒速度分布及飛散形態Fig.11 Velocity and distribution of propellant pellets in the chamber at different moments

3.3 模塊裝填的位置對藥粒分布形態的影響

為研究模塊裝填的初始位置對藥粒分布形態的影響,如圖12所示,將模塊初始裝填位置分別設置為距底火端10 mm、40 mm、70 mm和100 mm,數值模擬以上4種工況下單模塊裝藥藥盒破裂后藥粒飛散運動。由圖12可以看出,4種工況下藥室內終態藥粒從底火端向密閉端依次呈水平堆積和坡狀堆積形態。為更好地對藥粒堆積形態進行量化比較,引入了水平堆積軸向長度L和坡狀堆積坡度角α兩個特征參數,L表示沿藥室身管軸向水平堆積起點到坡狀堆積起點的距離,α則表示坡狀堆積斜坡表面與水平面的夾角。分別對以上4種工況下終態藥粒堆積表面進行L和α的測量,詳細數據如表1所示。

圖12 不同裝填位置下的終態藥粒堆積分布Fig.12 Distribution of pellets under different working conditions

表1 藥粒堆積特征參數

通過比較這兩種特征參數在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)

上述結果與陳安等[20]基于球形顆粒建立的模型所獲得的顆粒堆積變化趨勢一致,但與之相比,采用柱狀顆粒模型仿真獲得的藥粒堆積特性參數變化更小。之所以出現上述情況,是因為模塊藥盒右端蓋破裂瞬間,顆粒受到氣相壓強的作用,因為等體積球形顆粒垂直于壓力梯度的投影大于大部分柱形顆粒(與壓力梯度方向垂直或近似垂直的較少顆粒除外),球形顆粒獲得更大的壓力,故兩種顆粒在相同質量下,球形顆粒將獲得更大的加速度,柱狀顆粒與之相比到達藥室密閉端的速度更小,由圖11(b)和圖11(c)可知,顆粒最大速度約為74.73 m/s,到達右側密閉端的速度約為64 m/s,而球形顆粒的最大速度可達93 m/s[20]。故在后面較小壓強梯度下加速的顆粒距離藥室右側較遠距離時速度就不再增加,使得終態顆粒分布較均勻,坡狀堆積的坡度角較小。而隨著模塊藥盒的裝填位置右移,水平堆積軸向長度縮短,主要原因是模塊藥盒破裂時僅右端蓋全部破開,顆粒在氣相壓力的作用下向右運動,水平堆積起始位置向右移動;顆粒堆積坡度減小則是因為隨著模塊藥盒與藥室底火端的距離增大,顆粒到達藥室右壁面的距離減小、速度增大,無論是顆粒-壁面還是顆粒-顆粒發生碰撞后,將獲得更大的反向速度,藥粒將反向運動更遠的距離后落下,散布更均勻,使得坡狀堆積角變小,但由于距離變化較小對顆粒速影響比較小,這種變化也不是很明顯。

4 結論

本文針對單模塊裝藥點傳火過程中藥盒破裂、藥粒在藥室中的飛散和堆積問題,建立了三維非穩態氣-固兩相流模型并進行了數值模擬。得出如下主要結論:

1)藥粒在藥室中的終態分布呈現非均勻狀態,從底火端到密閉端依次呈水平堆積和坡狀堆積,與試驗結果基本吻合,驗證了該模型的合理性。

2)采用球元法建立了柱狀火藥顆粒的離散元模型,與基于等球形顆粒建立的模型相比,模擬仿真所得的結果更加精確,更接近實測結果,表明顆粒形狀對顆粒飛散運動行為有較大影響。

3)單模塊裝藥初始裝填位置與藥室底火端距離由10 mm增大到100 mm時,水平堆積軸向長度呈線性縮短,由293.3 mm縮短為203.8 mm,同時坡狀堆積坡度角由22.01°降為19.82°,堆積逐漸變平緩。