應(yīng)迪通,張彥春
(中國(guó)刑事警察學(xué)院,沈陽(yáng) 110854)
全球恐怖活動(dòng)頻發(fā),涉爆案件時(shí)有發(fā)生。偵查人員在面對(duì)復(fù)雜的爆炸案件現(xiàn)場(chǎng)時(shí),應(yīng)當(dāng)具備優(yōu)秀的爆炸現(xiàn)場(chǎng)分析能力。爆炸現(xiàn)場(chǎng)分析是基于爆炸案件現(xiàn)場(chǎng)、爆炸試驗(yàn)以及爆炸仿真的理論方法。由于爆炸試驗(yàn)成本高、周期長(zhǎng)、危險(xiǎn)度高,而使用計(jì)算軟件模擬爆炸過(guò)程,可以形象地展示結(jié)構(gòu)的破壞和變形規(guī)律,因此爆炸有限元仿真成為了相關(guān)學(xué)者的研究重心。
爆炸有限元仿真常使用ANSYS/LS-DYNA和LS-PrePost軟件進(jìn)行操作。使用ANSYS/LS-DYNA 8.1進(jìn)行實(shí)體建模的思路與AutoCAD相似,先構(gòu)建炸藥的幾何體模型和空氣的立方體模型,通過(guò)布爾操作和轉(zhuǎn)換工作平面的多次切割,方能構(gòu)建出完整的炸藥模型[1]。若使用LS-PrePost進(jìn)行建模,過(guò)程相對(duì)而言會(huì)簡(jiǎn)潔一些。以建立球形裝藥在無(wú)限空域內(nèi)爆炸的1/8模型為例,選擇Solid欄中的Box功能,輸入?yún)?shù)建立立方體模型,接著使用Sphere功能直接建立1/8球體,最后使用Boolean功能切割模型即可。方形、柱形裝藥的幾何模型同樣適用此方法。劃分網(wǎng)格常用的做法是使用Tetrahedron Mesher功能,在mesh model欄目中選擇Pick skin Geometry,設(shè)定好合適的Edge Size,選定目標(biāo)幾何體后點(diǎn)擊Tria Mesh使幾何體表面殼單元?jiǎng)澐滞戤?,點(diǎn)下方TetMesh完成實(shí)體有限元模型的劃分。隨后合并炸藥模型和空氣模型的共節(jié)點(diǎn),以滿足流體歐拉算法。上述方法屬于自由網(wǎng)格劃分方法,操作簡(jiǎn)單,適用于立方體的炸藥模型,而對(duì)于球形、柱形裝藥的有限元模型網(wǎng)格劃分卻不容易得到優(yōu)質(zhì)的網(wǎng)格[2]。
由于方形裝藥的建模操作簡(jiǎn)單、網(wǎng)格質(zhì)量高,當(dāng)需要降低建模成本,又對(duì)仿真結(jié)果精度要求不高時(shí),實(shí)際的仿真工作常常會(huì)將球形裝藥或者圓柱形裝藥簡(jiǎn)化為立方體模型。然而在大部分爆炸案件中,犯罪嫌疑人往往采用的是近似球形裝藥、近似柱形裝藥,甚至是不規(guī)則形狀的炸藥。在進(jìn)行爆炸仿真時(shí)使用的裝藥形狀應(yīng)盡量與實(shí)際裝藥形狀保持一致,否則會(huì)降低爆炸仿真結(jié)果的準(zhǔn)確度。為解決球形裝藥建模過(guò)程中凸顯的問(wèn)題,本文將在第一部分專門(mén)研究球形裝藥建模方法,再以爆炸仿真中裝藥形狀過(guò)度簡(jiǎn)化的問(wèn)題為出發(fā)點(diǎn),基于LS-PrePost軟件研究不同形狀裝藥所產(chǎn)生爆炸空氣沖擊波的差異。該項(xiàng)研究具有一定的理論意義,對(duì)涉爆案件的偵破具有重要的實(shí)際應(yīng)用價(jià)值。
查閱資料后發(fā)現(xiàn),基于ANSYS/LS-DYNA和LS-PrePost軟件對(duì)球形裝藥進(jìn)行建模的有多種方法,除使用ANSYS/LS-DYNA和LS-PrePost軟件的基礎(chǔ)建模功能以外有以下三種:1)使用關(guān)鍵字*LOAD_BLAST建模。2)使用SPH方法建模。3)使用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建模。
本文通過(guò)總結(jié)文獻(xiàn)資料和進(jìn)行實(shí)際建模操作的方式研究三種球形裝藥建模方法,分析比較了三種方法在操作方法、計(jì)算時(shí)間和適用場(chǎng)景應(yīng)用上的優(yōu)缺點(diǎn)。
LS-DYNA 中儲(chǔ)存的TNT空中爆炸超壓實(shí)驗(yàn)數(shù)據(jù)可以利用關(guān)鍵字按以下方法直接調(diào)用。首先采用*LOAD_BLAST_ENHANCED關(guān)鍵字定義TNT當(dāng)量、爆炸點(diǎn)及爆炸類型,再結(jié)合關(guān)鍵字*LOAD_SEGMENT與*LOAD_BLAST_ SEGMENT_SET定義爆炸波的接觸面。此方法不需要建立空氣與炸藥的數(shù)值分析模型,施加方法簡(jiǎn)單,且能滿足大部分爆炸案件現(xiàn)場(chǎng)的要求。
此種建模方法借用了CONWEP模型,它考慮了空氣的可壓縮性與輕質(zhì)的特點(diǎn),忽略了空氣的剛度與慣性效應(yīng),因此使用它進(jìn)行空中爆炸分析時(shí)不用為空氣介質(zhì)建模,也就是不用考慮爆炸空氣沖擊波在空氣中的傳播過(guò)程,僅有結(jié)構(gòu)模型參與爆炸響應(yīng)計(jì)算就可以,這樣在較高的精度范圍內(nèi),相比于ALE方法來(lái)說(shuō)大大節(jié)省了爆炸分析的計(jì)算量[3]。此方法適用于遠(yuǎn)距離爆炸場(chǎng)景中玻璃介質(zhì)在爆炸沖擊下的動(dòng)力響應(yīng)研究,例如張文嵐等使用了CONWEP計(jì)算模型方法對(duì)爆炸荷載作用下的高速列車(chē)車(chē)窗玻璃進(jìn)行了動(dòng)態(tài)響應(yīng)研究[4]。
SPH(smoothed particle hydrodynamics)全稱光滑粒子流體動(dòng)力學(xué)方法,是一種無(wú)網(wǎng)格的純Lagrange方法,與有限元法、有限差分法等基于網(wǎng)格的數(shù)值方法不同,它用一系列隨機(jī)分布的粒子質(zhì)點(diǎn)來(lái)代表整個(gè)連續(xù)介質(zhì)流體[5]。這些離散出來(lái)的光滑粒子承載著所有的力學(xué)量,也能夠表示材料,并確保了描述材料物理行為的狀態(tài)方程能夠進(jìn)行有效和穩(wěn)定的模擬。不同材料的粒子相對(duì)運(yùn)動(dòng)便形成界面滑移,且這些粒子能夠按照計(jì)算公式任意流動(dòng),克服了許多基于網(wǎng)格的方法在求解過(guò)程中存在的問(wèn)題,從而比較自然地模擬高速碰撞等物理現(xiàn)象[6]。在研究炸藥爆炸及其作用規(guī)律的工作中,SPH方法不失為一種科學(xué)有效的方法。
建立SPH模型十分簡(jiǎn)單,在LS-PrePost的Mesh界面選擇SPH Generation功能,輸入球心坐標(biāo)、半徑、粒子分布密度即可,不需要建立空氣模型;或者使用Shape Mesher事先建好實(shí)體模型,再通過(guò)SPH Generation將實(shí)體模型轉(zhuǎn)化為SPH模型。具體操作的參數(shù)選擇和設(shè)置如圖1所示。
圖1 建立SPH模型(a:直接建立SPH模型;b:將實(shí)體模型轉(zhuǎn)換為SPH模型)Fig.1 SPH modeling (a: direct modeling; b: one converted from an ontic model)
此方法在傳統(tǒng)方法基礎(chǔ)上有所改善,是一種全新的爆炸載荷施加方式,該方法首先需要定義炸藥、空氣的材料參數(shù)與狀態(tài)方程;然后借助關(guān)鍵字*ALE_MUTI-MATERIAL_GROUP定義一個(gè)包含炸藥與空氣的材料組;最后借助關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMTRY進(jìn)行填充。在LS-DYNA關(guān)鍵字手冊(cè)中詳細(xì)記載了關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMTRY的用法,現(xiàn)以如何施加球形爆炸載荷為例介紹此方法。
圖2中所建模型為1 m×1 m×1 m的立方體空氣與半徑5.305 cm的球形炸藥,當(dāng)采用傳統(tǒng)建模時(shí)需分別單獨(dú)建立空氣與炸藥的網(wǎng)格,而采用該方法建模時(shí)僅需要建立空氣網(wǎng)格作為背景網(wǎng)格,然后利用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY在指定區(qū)域?qū)⒄ㄋ幉牧咸畛渲两ê玫目諝饩W(wǎng)格中,并且不要求空氣網(wǎng)格和炸藥網(wǎng)格共節(jié)點(diǎn)。操作完成后,在炸藥幾何模型所占據(jù)的空間中,原網(wǎng)格中的空氣被炸藥材料所替代。需要注意的是,無(wú)論所選擇的填充形狀是什么,最終的被填充區(qū)域是炸藥模型與空氣模型相交的區(qū)域,即圖中的1/8球體。
關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY的定義如圖2所示,共3行。第1行的4個(gè)關(guān)鍵字定義空氣背景,其中FMSID表示空氣的PART號(hào),F(xiàn)MIDITYP表示PART類型,BAMMG表示空氣材料在材料組中的編號(hào),NTRACE表示填充采樣點(diǎn)數(shù)量,常用默認(rèn)值。第2行的3個(gè)關(guān)鍵字定義炸藥,其中CONTTYP表示炸藥的形狀,“6”為球形(通過(guò)選取不同的編號(hào)可模擬不同的形狀,例如通過(guò)編號(hào)“5”可以定義長(zhǎng)方體或立方體裝藥),F(xiàn)ILLOPT表示填充方式,“0”表示內(nèi)部填充,F(xiàn)AMMG表示炸藥材料在材料組中的編號(hào)。第3行定義炸藥的具體位置和形體參數(shù),其中 XC、YC、ZC表示炸藥球心坐標(biāo),RADIUS表示球體半徑,即裝藥量的大??;參數(shù)設(shè)置完畢后點(diǎn)擊右側(cè)Insert按鈕,并點(diǎn)擊上方Accept保存。圖中亮黃色網(wǎng)狀球體即為建立的球形裝藥示意圖。
圖2 使用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建模Fig.2 Modeling with the keyword *INITIAL_VOLUME_FRACTION_GEOMETRY-led way
使用關(guān)鍵字*LOAD_BLAST建模具有局限性,其可供選擇的爆炸類型只有表面接觸爆炸和空中爆炸,裝藥的形狀只有球形或半球形,且不能對(duì)表面形狀復(fù)雜或有遮擋的考察對(duì)象進(jìn)行動(dòng)力分析。需要注意的是,在近爆情況下,該方法所調(diào)用的經(jīng)驗(yàn)公式精確度并不高。因此在滿足以下三種爆炸條件的情況下使用關(guān)鍵字*LOAD_BLAST建模會(huì)更有優(yōu)勢(shì):1)球形裝藥;2)遠(yuǎn)距離爆炸;3)爆炸波接觸表面幾何結(jié)構(gòu)簡(jiǎn)單。而使用SPH方法的優(yōu)點(diǎn)是無(wú)需網(wǎng)格,沒(méi)有網(wǎng)格畸變問(wèn)題,形成的爆轟效果相對(duì)而言更加逼真,可以形象地模擬炸藥爆轟及產(chǎn)物膨脹飛散過(guò)程。在使用SPH方法時(shí),由于FEM算法和SPH算法各自存在缺陷,通常會(huì)結(jié)合Lagrange和SPH算法,即SPH-FEM耦合法。一部分學(xué)者[7-9]認(rèn)為該方法可避免FEM算法的網(wǎng)格畸變問(wèn)題,同時(shí)在重點(diǎn)部位發(fā)揮SPH算法優(yōu)勢(shì),可以大大減少計(jì)算時(shí)間。
使用關(guān)鍵字* INITIAL_VOLUME_FRACTION_GEOMETRY建模則不需要單獨(dú)建立炸藥單元,僅僅通過(guò)關(guān)鍵字修改相應(yīng)參數(shù),就可以建立不同形狀的模型、快速改變裝藥量,操作靈活,減少了重復(fù)建模的繁瑣過(guò)程,也能避免復(fù)雜模型劃分網(wǎng)格的困難。
在各種爆炸場(chǎng)景中,不同的物理參數(shù)直接對(duì)爆炸空氣沖擊波的超壓峰值造成影響,例如裝藥形狀、炸藥密度、起爆點(diǎn)的空間位置、炸藥的濕度、環(huán)境溫度等參數(shù)。本文嘗試將裝藥形狀參數(shù)與爆炸空氣沖擊波響應(yīng)值聯(lián)系起來(lái),使用上文提到的不同形狀裝藥建模方法進(jìn)行爆炸仿真實(shí)驗(yàn),研究不同形狀裝藥爆炸所形成的空氣沖擊波超壓??紤]到使用關(guān)鍵字* INITIAL_ VOLUME_FRACTION_GEOMETRY可以在不需要重新建模的情況下快速建立不同形狀裝藥的有限元模型,本文將使用該方法建立爆炸仿真模型,研究不同形狀裝藥在無(wú)限空氣域中的爆炸。
為了模擬真實(shí)爆炸案件中常見(jiàn)的近距離爆炸,本次爆炸仿真的測(cè)試點(diǎn)在0~1 m之間。爆炸仿真實(shí)驗(yàn)采用cm-g-μs單位制,使用LS-PrePost進(jìn)行前后處理,使用ANSYSLS-DYNA求解器進(jìn)行計(jì)算。TNT炸藥當(dāng)量定為1 000 g,密度是1.6 kg/cm3。通過(guò)關(guān)鍵字分別建立球形、柱形和方形裝藥。其中球形裝藥半徑為5.305 cm;為了使柱形裝藥盡量接近球體,設(shè)定圓柱體長(zhǎng)徑比為1,直徑和高均為9.260 cm;方形裝藥邊長(zhǎng)為8.550 cm。
三種裝藥均設(shè)置為中心點(diǎn)起爆。使用Shape Mesher功能建立100 cm×100 cm×100 cm的立方體空氣模型,單元?jiǎng)澐殖叽鐬? cm;采用SOLID164單元和歐拉算法。劃分立方體單元共1 000 000個(gè);三種工況的炸藥模型分別通過(guò)關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建立??紤]到此次爆炸仿真的對(duì)稱性,取八分之一實(shí)體進(jìn)行建模,從而減少計(jì)算量。
2.2.1 TNT本構(gòu)模型和狀態(tài)方程參數(shù)
TNT炸藥采用目前比較通用的HIGH_EXPLOSIVE_BURN模型以及JWL狀態(tài)方程加以描述,相關(guān)參數(shù)如表1所示[10]。
表1 炸藥材料參數(shù)Table 1 Parameters for explosive as material
2.2.2 空氣本構(gòu)模型和狀態(tài)方程參數(shù)
空氣采用NULL材料模型以及LINEAR_POLYNOMIAL狀態(tài)方程加以描述,在數(shù)值分析中假設(shè)標(biāo)準(zhǔn)大氣條件,即平均海平面和標(biāo)準(zhǔn)溫度為15 ℃,大氣壓力為101.325 kPa[11],相關(guān)參數(shù)如表2所示[10]。
表2 空氣材料參數(shù)Table 2 Parameters of air as material
無(wú)限空氣域中的爆炸是炸藥在空氣中爆炸的一種理想模型,可以使用此模型的仿真結(jié)果與空氣中炸藥爆炸的過(guò)程做對(duì)比。炸藥在空氣中爆炸會(huì)瞬時(shí)轉(zhuǎn)變?yōu)楦邷馗邏旱谋óa(chǎn)物,爆炸產(chǎn)物在空氣中進(jìn)行膨脹,其結(jié)果是在爆炸產(chǎn)物內(nèi)形成稀疏波。同時(shí),爆炸產(chǎn)物強(qiáng)烈壓縮空氣,在空氣中形成爆炸空氣沖擊波[12]。
通過(guò)仿真手段可以將空中爆炸過(guò)程準(zhǔn)確形象地展現(xiàn)出來(lái),圖3是1 kg球形TNT爆炸后在不同時(shí)刻的壓強(qiáng)等值面效果圖,以壓力云圖的方式表現(xiàn)了爆炸空氣沖擊波的傳播過(guò)程。如圖所示,爆炸空氣沖擊波從起爆中心開(kāi)始,以球形曲面的傳播方式向外圍傳播。爆炸空氣沖擊波伴隨著傳播距離的增加能量逐步耗散,在爆炸傳播的初始階段,超壓峰值急速衰減,而后衰減過(guò)程較為平緩。其原因是在總能量不變的情況下,波陣面的表面積急速擴(kuò)大,通過(guò)單位面積的能量就急劇減??;另一方面是空氣受到了絕熱壓縮使得溫度瞬間升高,產(chǎn)生了不可逆的能量損失,消耗了部分能量[13]。
圖3 球形裝藥爆炸空氣沖擊波在不同時(shí)刻的壓強(qiáng)等值面圖Fig.3 The moment-different pressure-equivalent surfaces from the exploding air-shock waves with spherical charging
為驗(yàn)證此基礎(chǔ)模型是否符合爆炸空氣沖擊波在自由空氣中的傳播規(guī)律,本文將基礎(chǔ)模型的仿真結(jié)果與經(jīng)驗(yàn)公式預(yù)測(cè)值進(jìn)行對(duì)比。常用來(lái)預(yù)測(cè)爆炸空氣沖擊波超壓值的經(jīng)驗(yàn)公式有兩種:TM5-1300手冊(cè)數(shù)據(jù)擬合公式和Henrych J公式,這兩種公式在國(guó)際上受到廣泛認(rèn)可,近爆區(qū)的準(zhǔn)確性相對(duì)較高。
表3和圖4顯示了基礎(chǔ)模型仿真結(jié)果與經(jīng)驗(yàn)公式計(jì)算結(jié)果的對(duì)比情況,圖4中仿真的結(jié)果和經(jīng)驗(yàn)公式的整體趨勢(shì)較為一致,結(jié)合圖3中爆炸空氣沖擊波陣面的傳播形態(tài),可以得知此模型符合爆炸空氣沖擊波的傳播規(guī)律。
圖4 爆炸仿真結(jié)果與經(jīng)驗(yàn)公式計(jì)算結(jié)果對(duì)比Fig.4 Comparison of the results from simulative explosion with those from the empirical equation’s calculation
表3 爆炸仿真結(jié)果與經(jīng)驗(yàn)公式計(jì)算結(jié)果對(duì)比Table 3 Comparison of results from numerical simulative explosion with those from empirical equation’s calculation
從圖4中還可以看出,仿真結(jié)果相對(duì)偏小,尤其是在比例距離小于0.5 m·kg-1/3時(shí),仿真結(jié)果小于經(jīng)驗(yàn)公式的計(jì)算值。其原因之一可能是爆炸的強(qiáng)瞬時(shí)性以及實(shí)驗(yàn)裝置的采樣頻率等受到限制,使得各種實(shí)驗(yàn)結(jié)果的精確數(shù)據(jù)很難在實(shí)驗(yàn)中被測(cè)得,所獲取的數(shù)據(jù)存在很大的離散性,原因之二是仿真建模難以完全還原采集實(shí)驗(yàn)數(shù)據(jù)的真實(shí)爆炸場(chǎng)景,設(shè)備之間或者地面墻體間的反射波會(huì)增強(qiáng)爆炸空氣沖擊波的威力。
由表3和圖4可見(jiàn),爆炸仿真的超壓峰值數(shù)據(jù)隨著比例距離的增大而越來(lái)越接近經(jīng)驗(yàn)公式的數(shù)據(jù),且變化趨勢(shì)與Henrych J公式、TM5-1300手冊(cè)數(shù)據(jù)擬合公式均一致,特別是數(shù)值仿真結(jié)果與Henrych J公式結(jié)果非常接近,由此說(shuō)明,本文建立的爆炸數(shù)值仿真模型是有效可靠的。
在經(jīng)過(guò)驗(yàn)證的有限元模型基礎(chǔ)上,采用控制變量的方法改變爆源形狀參數(shù),使用ANSYSLS-DYNA求解器依次對(duì)球形、柱形、方形裝藥的有限元爆炸模型進(jìn)行求解,得到不同比例距離下的超壓峰值以及相對(duì)應(yīng)的超壓峰值到達(dá)時(shí)間,選取部分有代表性的爆炸仿真結(jié)果,作圖5來(lái)展示三種不同裝藥方式在各比例距離處的超壓峰值對(duì)比情況;分析柱形和方形裝藥相對(duì)于球形裝藥的差值,結(jié)果如表4所示。
圖5 三種裝藥方式在不同比例距離處的超壓峰值對(duì)比圖Fig.5 Comparison among the three charging ways about their peak overpressures at different scaling distances
表4 三種裝藥方式分別在不同比例距離處的超壓峰值及其差值Table 4 Peak overpressures and differences at diverse scaling distances from the indicated three charging ways
爆炸仿真結(jié)果表明:1 kg TNT空中爆炸時(shí),選擇0.1、0.4、0.6、1.0 m·kg-1/3距離處作為四個(gè)觀測(cè)點(diǎn),柱形裝藥與球形裝藥的超壓峰值比值分別為:1.44、1.51、1.496、1.44,差值分別為9、1.73、0.93、0.50 MPa;而方形裝藥與球形裝藥的超壓峰值比值分別為:1.66、2.08、2.29、2.89,差值分別為13.40、3.63、2.41、2.14 MPa。可見(jiàn),在比例距離小于等于1.0 m·kg-1/3時(shí),柱形裝藥的超壓峰值與球形裝藥更接近,而方形裝藥超壓峰值與球形裝藥差異更大;隨著比例距離的增加,柱形裝藥的峰值超壓與球形裝藥的比值變化不大,而方形裝藥的峰值超壓與球形裝藥的比值逐漸增加,變化較大。
本文基于LS-PrePost軟件,研究了使用關(guān)鍵字*LOAD_BLAST建模、SPH法建模以及關(guān)鍵字*INITIAL_VOLUME_ FRACTION_GEOMTRY建模等三種常見(jiàn)建模方法的優(yōu)缺點(diǎn)和適用場(chǎng)景,利用關(guān)鍵字*INITIAL_VOLUME_FRACTION _GEOMTRY建模方法,建立了球形裝藥有限元模型,研究得到了1 kg TNT在不同裝藥形狀下,爆炸產(chǎn)生的空氣沖擊波超壓的不同。研究表明,不同裝藥形狀對(duì)于爆炸產(chǎn)生的沖擊波超壓的影響是不容忽視的。在研究實(shí)際爆炸案件的裝藥及其爆炸破壞作用時(shí),應(yīng)盡可能采用最接近實(shí)際情況的裝藥形狀,在通過(guò)爆炸仿真方法研究問(wèn)題時(shí),也要選擇合適的建模方法和裝藥條件,才能保證分析結(jié)果的準(zhǔn)確可靠。