應(yīng)迪通，張彥春

（中國刑事警察學(xué)院，沈陽 110854）

全球恐怖活動頻發(fā)，涉爆案件時有發(fā)生。偵查人員在面對復(fù)雜的爆炸案件現(xiàn)場時，應(yīng)當(dāng)具備優(yōu)秀的爆炸現(xiàn)場分析能力。爆炸現(xiàn)場分析是基于爆炸案件現(xiàn)場、爆炸試驗以及爆炸仿真的理論方法。由于爆炸試驗成本高、周期長、危險度高，而使用計算軟件模擬爆炸過程，可以形象地展示結(jié)構(gòu)的破壞和變形規(guī)律，因此爆炸有限元仿真成為了相關(guān)學(xué)者的研究重心。

爆炸有限元仿真常使用ANSYS/LS-DYNA和LS-PrePost軟件進(jìn)行操作。使用ANSYS/LS-DYNA 8.1進(jìn)行實體建模的思路與AutoCAD相似，先構(gòu)建炸藥的幾何體模型和空氣的立方體模型，通過布爾操作和轉(zhuǎn)換工作平面的多次切割，方能構(gòu)建出完整的炸藥模型[1]。若使用LS-PrePost進(jìn)行建模，過程相對而言會簡潔一些。以建立球形裝藥在無限空域內(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)幾何體后點擊Tria Mesh使幾何體表面殼單元劃分完畢，點下方TetMesh完成實體有限元模型的劃分。隨后合并炸藥模型和空氣模型的共節(jié)點，以滿足流體歐拉算法。上述方法屬于自由網(wǎng)格劃分方法，操作簡單，適用于立方體的炸藥模型，而對于球形、柱形裝藥的有限元模型網(wǎng)格劃分卻不容易得到優(yōu)質(zhì)的網(wǎng)格[2]。

由于方形裝藥的建模操作簡單、網(wǎng)格質(zhì)量高，當(dāng)需要降低建模成本，又對仿真結(jié)果精度要求不高時，實際的仿真工作常常會將球形裝藥或者圓柱形裝藥簡化為立方體模型。然而在大部分爆炸案件中，犯罪嫌疑人往往采用的是近似球形裝藥、近似柱形裝藥，甚至是不規(guī)則形狀的炸藥。在進(jìn)行爆炸仿真時使用的裝藥形狀應(yīng)盡量與實際裝藥形狀保持一致，否則會降低爆炸仿真結(jié)果的準(zhǔn)確度。為解決球形裝藥建模過程中凸顯的問題，本文將在第一部分專門研究球形裝藥建模方法，再以爆炸仿真中裝藥形狀過度簡化的問題為出發(fā)點，基于LS-PrePost軟件研究不同形狀裝藥所產(chǎn)生爆炸空氣沖擊波的差異。該項研究具有一定的理論意義，對涉爆案件的偵破具有重要的實際應(yīng)用價值。

1 球形裝藥建模方法

查閱資料后發(fā)現(xiàn)，基于ANSYS/LS-DYNA和LS-PrePost軟件對球形裝藥進(jìn)行建模的有多種方法，除使用ANSYS/LS-DYNA和LS-PrePost軟件的基礎(chǔ)建模功能以外有以下三種：1）使用關(guān)鍵字*LOAD_BLAST建模。2）使用SPH方法建模。3）使用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建模。

本文通過總結(jié)文獻(xiàn)資料和進(jìn)行實際建模操作的方式研究三種球形裝藥建模方法，分析比較了三種方法在操作方法、計算時間和適用場景應(yīng)用上的優(yōu)缺點。

1.1 使用關(guān)鍵字*LOAD_BLAST建模

LS-DYNA 中儲存的TNT空中爆炸超壓實驗數(shù)據(jù)可以利用關(guān)鍵字按以下方法直接調(diào)用。首先采用*LOAD_BLAST_ENHANCED關(guān)鍵字定義TNT當(dāng)量、爆炸點及爆炸類型，再結(jié)合關(guān)鍵字*LOAD_SEGMENT與*LOAD_BLAST_ SEGMENT_SET定義爆炸波的接觸面。此方法不需要建立空氣與炸藥的數(shù)值分析模型，施加方法簡單，且能滿足大部分爆炸案件現(xiàn)場的要求。

此種建模方法借用了CONWEP模型，它考慮了空氣的可壓縮性與輕質(zhì)的特點，忽略了空氣的剛度與慣性效應(yīng)，因此使用它進(jìn)行空中爆炸分析時不用為空氣介質(zhì)建模，也就是不用考慮爆炸空氣沖擊波在空氣中的傳播過程，僅有結(jié)構(gòu)模型參與爆炸響應(yīng)計算就可以，這樣在較高的精度范圍內(nèi)，相比于ALE方法來說大大節(jié)省了爆炸分析的計算量[3]。此方法適用于遠(yuǎn)距離爆炸場景中玻璃介質(zhì)在爆炸沖擊下的動力響應(yīng)研究，例如張文嵐等使用了CONWEP計算模型方法對爆炸荷載作用下的高速列車車窗玻璃進(jìn)行了動態(tài)響應(yīng)研究[4]。

1.2 使用SPH方法建模

SPH（smoothed particle hydrodynamics）全稱光滑粒子流體動力學(xué)方法，是一種無網(wǎng)格的純Lagrange方法，與有限元法、有限差分法等基于網(wǎng)格的數(shù)值方法不同，它用一系列隨機分布的粒子質(zhì)點來代表整個連續(xù)介質(zhì)流體[5]。這些離散出來的光滑粒子承載著所有的力學(xué)量，也能夠表示材料，并確保了描述材料物理行為的狀態(tài)方程能夠進(jìn)行有效和穩(wěn)定的模擬。不同材料的粒子相對運動便形成界面滑移，且這些粒子能夠按照計算公式任意流動，克服了許多基于網(wǎng)格的方法在求解過程中存在的問題，從而比較自然地模擬高速碰撞等物理現(xiàn)象[6]。在研究炸藥爆炸及其作用規(guī)律的工作中，SPH方法不失為一種科學(xué)有效的方法。

建立SPH模型十分簡單，在LS-PrePost的Mesh界面選擇SPH Generation功能，輸入球心坐標(biāo)、半徑、粒子分布密度即可，不需要建立空氣模型；或者使用Shape Mesher事先建好實體模型，再通過SPH Generation將實體模型轉(zhuǎn)化為SPH模型。具體操作的參數(shù)選擇和設(shè)置如圖1所示。

圖1 建立SPH模型（a：直接建立SPH模型；b：將實體模型轉(zhuǎn)換為SPH模型）Fig.1 SPH modeling (a: direct modeling; b: one converted from an ontic model)

1.3 使用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建模

此方法在傳統(tǒng)方法基礎(chǔ)上有所改善，是一種全新的爆炸載荷施加方式，該方法首先需要定義炸藥、空氣的材料參數(shù)與狀態(tài)方程；然后借助關(guān)鍵字*ALE_MUTI-MATERIAL_GROUP定義一個包含炸藥與空氣的材料組；最后借助關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMTRY進(jìn)行填充。在LS-DYNA關(guān)鍵字手冊中詳細(xì)記載了關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMTRY的用法，現(xiàn)以如何施加球形爆炸載荷為例介紹此方法。

圖2中所建模型為1 m×1 m×1 m的立方體空氣與半徑5.305 cm的球形炸藥，當(dāng)采用傳統(tǒng)建模時需分別單獨建立空氣與炸藥的網(wǎng)格，而采用該方法建模時僅需要建立空氣網(wǎng)格作為背景網(wǎng)格，然后利用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY在指定區(qū)域?qū)⒄ㄋ幉牧咸畛渲两ê玫目諝饩W(wǎng)格中，并且不要求空氣網(wǎng)格和炸藥網(wǎng)格共節(jié)點。操作完成后，在炸藥幾何模型所占據(jù)的空間中，原網(wǎng)格中的空氣被炸藥材料所替代。需要注意的是，無論所選擇的填充形狀是什么，最終的被填充區(qū)域是炸藥模型與空氣模型相交的區(qū)域，即圖中的1/8球體。

關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY的定義如圖2所示，共3行。第1行的4個關(guān)鍵字定義空氣背景，其中FMSID表示空氣的PART號，F(xiàn)MIDITYP表示PART類型，BAMMG表示空氣材料在材料組中的編號，NTRACE表示填充采樣點數(shù)量，常用默認(rèn)值。第2行的3個關(guān)鍵字定義炸藥，其中CONTTYP表示炸藥的形狀，“6”為球形（通過選取不同的編號可模擬不同的形狀，例如通過編號“5”可以定義長方體或立方體裝藥），F(xiàn)ILLOPT表示填充方式，“0”表示內(nèi)部填充，F(xiàn)AMMG表示炸藥材料在材料組中的編號。第3行定義炸藥的具體位置和形體參數(shù)，其中 XC、YC、ZC表示炸藥球心坐標(biāo)，RADIUS表示球體半徑，即裝藥量的大小；參數(shù)設(shè)置完畢后點擊右側(cè)Insert按鈕，并點擊上方Accept保存。圖中亮黃色網(wǎng)狀球體即為建立的球形裝藥示意圖。

圖2 使用關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建模Fig.2 Modeling with the keyword *INITIAL_VOLUME_FRACTION_GEOMETRY-led way

1.4 不同建模方法優(yōu)缺點

使用關(guān)鍵字*LOAD_BLAST建模具有局限性，其可供選擇的爆炸類型只有表面接觸爆炸和空中爆炸，裝藥的形狀只有球形或半球形，且不能對表面形狀復(fù)雜或有遮擋的考察對象進(jìn)行動力分析。需要注意的是，在近爆情況下，該方法所調(diào)用的經(jīng)驗公式精確度并不高。因此在滿足以下三種爆炸條件的情況下使用關(guān)鍵字*LOAD_BLAST建模會更有優(yōu)勢：1）球形裝藥；2）遠(yuǎn)距離爆炸；3）爆炸波接觸表面幾何結(jié)構(gòu)簡單。而使用SPH方法的優(yōu)點是無需網(wǎng)格，沒有網(wǎng)格畸變問題，形成的爆轟效果相對而言更加逼真，可以形象地模擬炸藥爆轟及產(chǎn)物膨脹飛散過程。在使用SPH方法時，由于FEM算法和SPH算法各自存在缺陷，通常會結(jié)合Lagrange和SPH算法，即SPH-FEM耦合法。一部分學(xué)者[7-9]認(rèn)為該方法可避免FEM算法的網(wǎng)格畸變問題，同時在重點部位發(fā)揮SPH算法優(yōu)勢，可以大大減少計算時間。

使用關(guān)鍵字* INITIAL_VOLUME_FRACTION_GEOMETRY建模則不需要單獨建立炸藥單元，僅僅通過關(guān)鍵字修改相應(yīng)參數(shù)，就可以建立不同形狀的模型、快速改變裝藥量，操作靈活，減少了重復(fù)建模的繁瑣過程，也能避免復(fù)雜模型劃分網(wǎng)格的困難。

2 不同形狀裝藥形成的爆炸空氣沖擊波

在各種爆炸場景中，不同的物理參數(shù)直接對爆炸空氣沖擊波的超壓峰值造成影響，例如裝藥形狀、炸藥密度、起爆點的空間位置、炸藥的濕度、環(huán)境溫度等參數(shù)。本文嘗試將裝藥形狀參數(shù)與爆炸空氣沖擊波響應(yīng)值聯(lián)系起來，使用上文提到的不同形狀裝藥建模方法進(jìn)行爆炸仿真實驗，研究不同形狀裝藥爆炸所形成的空氣沖擊波超壓。考慮到使用關(guān)鍵字* INITIAL_ VOLUME_FRACTION_GEOMETRY可以在不需要重新建模的情況下快速建立不同形狀裝藥的有限元模型，本文將使用該方法建立爆炸仿真模型，研究不同形狀裝藥在無限空氣域中的爆炸。

2.1 構(gòu)建空中爆炸有限元模型

為了模擬真實爆炸案件中常見的近距離爆炸，本次爆炸仿真的測試點在0～1 m之間。爆炸仿真實驗采用cm-g-μs單位制，使用LS-PrePost進(jìn)行前后處理，使用ANSYSLS-DYNA求解器進(jìn)行計算。TNT炸藥當(dāng)量定為1 000 g，密度是1.6 kg/cm3。通過關(guān)鍵字分別建立球形、柱形和方形裝藥。其中球形裝藥半徑為5.305 cm；為了使柱形裝藥盡量接近球體，設(shè)定圓柱體長徑比為1，直徑和高均為9.260 cm；方形裝藥邊長為8.550 cm。

三種裝藥均設(shè)置為中心點起爆。使用Shape Mesher功能建立100 cm×100 cm×100 cm的立方體空氣模型，單元劃分尺寸為1 cm；采用SOLID164單元和歐拉算法。劃分立方體單元共1 000 000個；三種工況的炸藥模型分別通過關(guān)鍵字*INITIAL_VOLUME_FRACTION_GEOMETRY建立。考慮到此次爆炸仿真的對稱性，取八分之一實體進(jìn)行建模，從而減少計算量。

2.2 材料參數(shù)選擇

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

2.3 無限空氣域中爆炸過程的仿真

無限空氣域中的爆炸是炸藥在空氣中爆炸的一種理想模型，可以使用此模型的仿真結(jié)果與空氣中炸藥爆炸的過程做對比。炸藥在空氣中爆炸會瞬時轉(zhuǎn)變?yōu)楦邷馗邏旱谋óa(chǎn)物，爆炸產(chǎn)物在空氣中進(jìn)行膨脹，其結(jié)果是在爆炸產(chǎn)物內(nèi)形成稀疏波。同時，爆炸產(chǎn)物強烈壓縮空氣，在空氣中形成爆炸空氣沖擊波[12]。

通過仿真手段可以將空中爆炸過程準(zhǔn)確形象地展現(xiàn)出來，圖3是1 kg球形TNT爆炸后在不同時刻的壓強等值面效果圖，以壓力云圖的方式表現(xiàn)了爆炸空氣沖擊波的傳播過程。如圖所示，爆炸空氣沖擊波從起爆中心開始，以球形曲面的傳播方式向外圍傳播。爆炸空氣沖擊波伴隨著傳播距離的增加能量逐步耗散，在爆炸傳播的初始階段，超壓峰值急速衰減，而后衰減過程較為平緩。其原因是在總能量不變的情況下，波陣面的表面積急速擴大，通過單位面積的能量就急劇減小；另一方面是空氣受到了絕熱壓縮使得溫度瞬間升高，產(chǎn)生了不可逆的能量損失，消耗了部分能量[13]。

圖3 球形裝藥爆炸空氣沖擊波在不同時刻的壓強等值面圖Fig.3 The moment-different pressure-equivalent surfaces from the exploding air-shock waves with spherical charging

2.4 有限元模型驗證

為驗證此基礎(chǔ)模型是否符合爆炸空氣沖擊波在自由空氣中的傳播規(guī)律，本文將基礎(chǔ)模型的仿真結(jié)果與經(jīng)驗公式預(yù)測值進(jìn)行對比。常用來預(yù)測爆炸空氣沖擊波超壓值的經(jīng)驗公式有兩種：TM5-1300手冊數(shù)據(jù)擬合公式和Henrych J公式，這兩種公式在國際上受到廣泛認(rèn)可，近爆區(qū)的準(zhǔn)確性相對較高。

表3和圖4顯示了基礎(chǔ)模型仿真結(jié)果與經(jīng)驗公式計算結(jié)果的對比情況，圖4中仿真的結(jié)果和經(jīng)驗公式的整體趨勢較為一致，結(jié)合圖3中爆炸空氣沖擊波陣面的傳播形態(tài)，可以得知此模型符合爆炸空氣沖擊波的傳播規(guī)律。

圖4 爆炸仿真結(jié)果與經(jīng)驗公式計算結(jié)果對比Fig.4 Comparison of the results from simulative explosion with those from the empirical equation’s calculation

表3 爆炸仿真結(jié)果與經(jīng)驗公式計算結(jié)果對比Table 3 Comparison of results from numerical simulative explosion with those from empirical equation’s calculation

從圖4中還可以看出，仿真結(jié)果相對偏小，尤其是在比例距離小于0.5 m·kg-1/3時，仿真結(jié)果小于經(jīng)驗公式的計算值。其原因之一可能是爆炸的強瞬時性以及實驗裝置的采樣頻率等受到限制，使得各種實驗結(jié)果的精確數(shù)據(jù)很難在實驗中被測得，所獲取的數(shù)據(jù)存在很大的離散性，原因之二是仿真建模難以完全還原采集實驗數(shù)據(jù)的真實爆炸場景，設(shè)備之間或者地面墻體間的反射波會增強爆炸空氣沖擊波的威力。

由表3和圖4可見，爆炸仿真的超壓峰值數(shù)據(jù)隨著比例距離的增大而越來越接近經(jīng)驗公式的數(shù)據(jù)，且變化趨勢與Henrych J公式、TM5-1300手冊數(shù)據(jù)擬合公式均一致，特別是數(shù)值仿真結(jié)果與Henrych J公式結(jié)果非常接近，由此說明，本文建立的爆炸數(shù)值仿真模型是有效可靠的。

2.5 爆炸仿真結(jié)果及分析

在經(jīng)過驗證的有限元模型基礎(chǔ)上，采用控制變量的方法改變爆源形狀參數(shù)，使用ANSYSLS-DYNA求解器依次對球形、柱形、方形裝藥的有限元爆炸模型進(jìn)行求解，得到不同比例距離下的超壓峰值以及相對應(yīng)的超壓峰值到達(dá)時間，選取部分有代表性的爆炸仿真結(jié)果，作圖5來展示三種不同裝藥方式在各比例距離處的超壓峰值對比情況；分析柱形和方形裝藥相對于球形裝藥的差值，結(jié)果如表4所示。

圖5 三種裝藥方式在不同比例距離處的超壓峰值對比圖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空中爆炸時，選擇0.1、0.4、0.6、1.0 m·kg-1/3距離處作為四個觀測點，柱形裝藥與球形裝藥的超壓峰值比值分別為：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?？梢?，在比例距離小于等于1.0 m·kg-1/3時，柱形裝藥的超壓峰值與球形裝藥更接近，而方形裝藥超壓峰值與球形裝藥差異更大；隨著比例距離的增加，柱形裝藥的峰值超壓與球形裝藥的比值變化不大，而方形裝藥的峰值超壓與球形裝藥的比值逐漸增加，變化較大。

3 結(jié)論

本文基于LS-PrePost軟件，研究了使用關(guān)鍵字*LOAD_BLAST建模、SPH法建模以及關(guān)鍵字*INITIAL_VOLUME_ FRACTION_GEOMTRY建模等三種常見建模方法的優(yōu)缺點和適用場景，利用關(guān)鍵字*INITIAL_VOLUME_FRACTION _GEOMTRY建模方法，建立了球形裝藥有限元模型，研究得到了1 kg TNT在不同裝藥形狀下，爆炸產(chǎn)生的空氣沖擊波超壓的不同。研究表明，不同裝藥形狀對于爆炸產(chǎn)生的沖擊波超壓的影響是不容忽視的。在研究實際爆炸案件的裝藥及其爆炸破壞作用時，應(yīng)盡可能采用最接近實際情況的裝藥形狀，在通過爆炸仿真方法研究問題時，也要選擇合適的建模方法和裝藥條件，才能保證分析結(jié)果的準(zhǔn)確可靠。