鋼箱內部爆炸破壞的SPH 數值模擬*

強洪夫，孫新亞，王 廣，黃拳章

（火箭軍工程大學，陜西 西安 710025）

橋梁是公路、鐵路等交通網的咽喉。為了減小恒載，對斜拉橋、懸索橋等大跨度橋廣泛采用鋼箱梁作其結構支撐的主梁。隨著恐怖襲擊的不斷演化，戰略性橋梁成為恐怖主義爆炸破壞的重點目標，因此，有必要對鋼箱結構內部爆炸過程中鋼箱的變形規律進行研究[1-4]。Tang 等[5]運用LS-DYNA 軟件研究了汽車炸彈爆炸沖擊作用下橋塔、橋墩、橋面結構的局部破壞模式；姚術健等[2,6]通過實驗，研究了鋼箱內部爆炸過程中內壁失效變形的模式, 并利用ANSYS 軟件分析了箱體的損傷特性。鋼箱內部爆炸相比外部爆炸而言，爆炸沖擊波在約束空間內多次反射、疊加、匯聚，該物理過程相當復雜，目前的研究大多局限于爆炸過程中鋼箱變形的宏觀現象。光滑粒子流體動力學(smoothed particle hydrodynamics,SPH)作為一種無網格的粒子方法，因其在計算過程中可避免網格重分和算法耦合，被廣泛用于研究爆炸沖擊問題。自Monaghan 等[7]修正的人工黏度解決了SPH 計算結果的非物理震蕩問題后，強洪夫等[8]、Qiang 等[9]和Liu 等[10]利用該方法對爆炸沖擊波沖擊鋼板的過程成功地進行了數值模擬。本文中，擬利用自編的SPH 程序對鋼箱內部爆炸過程進行數值模擬，并與實驗對比，驗證SPH 算法在模擬鋼箱內部爆炸問題上的有效性；分析爆炸過程中鋼箱表面壓力、應力波的傳播規律以及特定點處速度變化趨勢；得到鋼箱表面的von Mises 應力分布；進一步對鋼箱內不同位置爆炸的過程進行數值模擬，得到不同工況下鋼箱的失效變形模式，以期通過研究為鋼箱類結構的爆炸毀傷提供一種評估方法。

1 模型與算法

1.1 控制方程和狀態方程

1.1.1 控制方程

鋼箱的內部爆炸是絕熱過程，外力相對于內部驅動力可忽略不計，即可不考慮重力的影響，在流體控制方程中忽略黏性項。因此，其控制方程組表述為：

式中：ρ、v、p、e、x 分別為密度、速度、壓強、能量、位移。

1.1.2 狀態方程

數值模擬過程中，TNT 炸藥的爆轟氣體壓力通過標準的JWL 狀態方程[11]獲得：

對于鋼箱材料，選用Q235 鋼。爆炸過程中鋼箱的各向同性壓力p 利用Mie-Grüneisen 狀態方程[12]求得:

1.2 本構方程

鋼箱材料模型采用B?rvik 等[13]、Nilson[14]提出的修正Johnson-Cook 本構模型，該模型通過流動應力和失效應變描述金屬材料的大變形、高應變率等動態力學特性。其流動應力的表達式為：

1.3 控制方程的SPH 離散

SPH 算法是模擬流體流動的一種拉格朗日型粒子方法，通過使用一系列任意分布的粒子來求解具有各種邊界條件的積分方程或偏微分方程。SPH 方法通常通過核函數插值實現場變量或場變量梯度的插值，通過粒子近似實現對核函數估計積分表達式的粒子離散[9]。

本文中利用SPH 方法模擬鋼箱內部爆炸過程屬于該方法在具有大變形和高度非均勻的動力學極端情形下的應用。為了解決爆炸域內求解結果的非物理振蕩，且將沖擊面內的動能以熱能的形式耗散，引入人工黏性[7]，同時為了更好地處理爆炸過程中鋼箱大變形、大扭曲等引起的密度和光滑長度變化劇烈的問題，采用Qiang 等[9]提出的完全變光滑長度SPH 方法對鋼箱內部爆炸過程進行數值模擬，結合后的方程組如下：

式中： fi為 修正系數； Wij=W(xi-xj,h)為核函數，一般選用三次樣條函數；Πij為人工黏度；N 為求解區域內粒子總數；α、β 表示空間坐標軸方向； h 為光滑長度，在連續性方程中將其對時間求導，得：

式中：d 為空間維數。

2 鋼箱內部爆炸過程的數值模擬

2.1 數值模型

為了證明SPH 方法在數值模擬鋼箱結構內部爆炸過程時的有效性，本文中設定了與姚術健等[2]的實驗一致的工況條件，如圖1(a)所示，鋼箱模型為立方體形狀，邊長為600 mm，且采用了120 mm 的延伸結構，Q235 鋼板厚度為4 mm；TNT 炸藥為立方體形狀，邊長為40 mm。具體粒子配置如圖1（b）所示，粒子間距為4 mm，炸藥粒子數為1 331，鋼箱粒子數為265 600。

圖 1 鋼箱和TNT 炸藥試件模型圖Fig.1 Model diagrams for the steel box and TNT dynamite specimen

為保證箱體各壁面受力均衡，TNT 炸藥位于鋼箱內部正中心，且各邊與箱體對應邊平行，如圖2 所示。計算中，TNT 炸藥密度為1 500 kg/m3, 采用點起爆方式，起爆點設置在炸藥上表面中心點（-73.9，5.6，50.7），爆速為6 930 m/s，初始內能為4 290 kJ/kg，狀態方程選用JWL 狀態方程；鋼箱密度為7 800 kg/m3,使用Johnson-Cook(J-C)本構模型，具體材料參數A=369 MPa, B=684 MPa, n=0.73, C=0.083, m=1.7,Tm=775 K,=0.1 s-1, CV=875 J/(kg·K), D1=0.13, D2=0.13, D3=-1.5, D4=0.011, D5=0。在組成鋼箱的各板之間施加焊接錨固作用，由文獻[15]可知焊縫金屬強度不應低于母材。為簡化計算，忽略焊縫附近區域母板材料性能的變化，將焊縫與母材視為同一結構，即認為鋼箱是一體的。

圖 2 鋼箱內部爆炸算例模型Fig.2 The case model for the explosion inside the steel box

2.2 鋼箱表面變形情況

數值模擬中，核函數選用三次樣條函數，光滑長度取1.5 倍粒子間距，時間積分采用蛙跳格式，時間步長取0.1 μs。鋼箱內部爆炸的數值模擬中，鋼箱表面的變形情況如圖3 所示。起始時刻，鋼箱表面無變形，1.97 ms 時中心處產生明顯的向外鼓起，并向四周傳播變形；3.34 ms 時變形首次到達表面邊沿處，且表面角隅處發生變形；12.03 ms時，表面變形達到最大，并保持不變。

圖4 給出了鋼箱內部爆炸實驗和數值模擬的對比，圖中A 處鋼板發生卷邊破壞，B 處鋼板向外鼓起，C 處鋼板發生扭曲變形，通過對比發現，數值模擬和實驗所得鋼箱壁面失效變形模式吻合較好。圖5 給出了鋼箱內部爆炸過程中箱體壁面撓度的變化情況。從圖5 可以看出，85 μs 時鋼箱壁面開始變形，撓度值隨時間逐漸增大，12.03 ms 后基本保持不變。數值模擬得到的鋼箱壁面最終撓度值為31.80 mm，而實驗值為32.93 mm，實驗和數值模擬的相對誤差為3.5%。產生這一誤差的主要原因可能是：(1)炸藥理論設計的質量為96 g，由于密度不均勻等原因實驗中炸藥的真實質量卻為98.443 g；(2)本文數值模擬過程中Q235 鋼材料密度是均勻的，沒有考慮制作中的工藝誤差。綜上所述，可得SPH方法能夠數值模擬鋼箱內部爆炸過程中鋼箱壁面的變形情況。

圖 3 鋼箱內部爆炸過程中不同時刻箱體內壁的變形情況Fig.3 Its inner wall deformation at different times for the steel box subjected to inner explosion loading

圖 4 鋼箱壁面破壞的數值模擬和實驗對比Fig.4 Comparison of steel box wall damage between numerical simulation and experiment

圖 5 鋼箱壁面撓度變化趨勢Fig.5 Deflection in the center of the steel box wall varying with time

2.3 鋼箱表面壓力和von Mises 應力分布

圖 6 鋼箱內部爆炸過程中，壁面von Mises 應力分布Fig.6 Von Mises stress inside the steel wall during the explosion in the steel box

圖 7 鋼箱內部爆炸過程中，壁面的壓力分布Fig.7 Pressure distribution of the wall during the explosion in the steel box

裝藥在完全封閉的鋼箱內部爆炸時，爆炸沖擊波在鋼箱四周壁面上發生多次反射，沖擊波強度逐漸降低。圖6 給出了鋼箱內部爆炸時，鋼箱表面的von Mises 應力分布；圖7 給出了鋼箱內部爆炸過程中，鋼箱表面壓力分布情況。從圖6 和圖7 可得：鋼箱內部炸藥爆炸后，很快產生爆炸沖擊波并迅速（85.0 μs）到達鋼箱表面，爆心投影點附近出現壓強和應力最大值；隨后沖擊波沿著鋼箱的內表面向板與板的連接處（角隅位置）運動，同時爆心投影點處超壓峰值逐漸降低，436.1 μs 時各板向角隅運動的沖擊波在角隅處匯聚，形成較強的沖擊波壓力，并出現應力集中和應力增強現象；隨后各板反射的沖擊波與角隅處匯聚加強了的沖擊波向箱室中心匯聚，使箱室中心形成較高超壓，箱室中心超壓又會向四周傳播，此過程與起爆后傳播過程類似，鋼箱表面承受多次類似的沖擊載荷作用后，鋼箱表面壓強和von Mises 應力逐漸消失，且箱內維持較長時間的余壓。

2.4 鋼箱表面中心點處速度和壓強

圖8 給出了鋼箱內部爆炸過程中，鋼箱表面中心點處速度的變化趨勢；圖9 給出了鋼箱內部爆炸過程中，鋼箱表面中心點處壓強的變化趨勢。從圖8 可以看出：101.2 μs 時（如點A 所示），表面中心點的速度跳躍式增加至22.3 m/s，其后速度值又迅速減小，最終減小至0，且速度值在B 區域內出現上下微小浮動。從圖9 可以看出：85.0 μs 時，鋼箱表面中心點處的壓強迅速增大，達到峰值487 MPa； 131.6 μs 時表面壓強達到第二次峰值394 MPa；隨后壓強值上下擾動并逐漸減弱，出現擾動的主要原因是爆炸沖擊波在鋼箱表面上的來回反射；2 000.0 μs 時壓強值逐漸減弱直至消失。由文獻[16]可得：壓力和速度時程曲線形狀與實驗及經典沖擊波理論相吻合，即初次峰值較大，然后發生多次反射，速度和壓強逐漸減弱。

圖 8 鋼箱壁面中心點處的速度變化趨勢Fig.8 Velocity varying with time at the wall center of the steel box

圖 9 鋼箱壁面中心點處的壓強變化趨勢Fig.9 Pressure varying with time at the wall center of the steel box

3 鋼箱內部不同位置爆炸失效模式

為了進一步探討同一爆炸當量下，鋼箱內部不同位置爆炸時，鋼箱結構的破壞程度。本文設定了以下3 種計算工況（如圖10 所示）：工況1，炸藥位于鋼箱底部正中心；工況2，炸藥在底部中心正上方60 mm 處；工況3，炸藥位于底板與側板之間的角隅處，距離各板60 mm。

圖11 給出了3 種工況下鋼箱表面的失效變形模式。圖11(a)為工況1 鋼箱的變形情況，鋼箱底板局部特征破壞明顯，產生圓形的沖切破壞，斷裂口處相對比較平滑，且底板產生裂紋（A）；同時底板與側板的角隅處產生了長條棗核形狀的裂口（B）。圖11(b)為工況2 鋼箱的變形情況，鋼箱底板劇烈破壞，產生不規則的撕裂口，斷裂口呈鋸齒狀，且底板沒有產生比較規則的裂紋（A）；同時底板與側板的角隅處也產生了長條狀棗核形狀的裂口，裂口周圍向外凸起（B）。圖11(c)為工況3 鋼箱的變形情況，距離炸藥較近的角隅處發生嚴重破壞，斷裂口周圍鋼板發生扭曲變形（A）；同時發現，在遠離炸藥的3 個鋼板的角隅處發生了撕裂破壞（B）。結合各工況下鋼箱失效變形的情況，可得炸藥在鋼箱角隅處爆炸時破壞最為嚴重。

圖 10 鋼箱中炸藥位置示意圖(a=60 mm)Fig.10 The explosive location diagram in the steel box(a=60 mm)

圖 11 3 種工況下，鋼箱表面局部破壞示意圖Fig.11 Surface partial destruction diagrams of the steel box under three different conditions

4 結 論

采用光滑粒子流體動力學(smoothed particle hydrodynamics, SPH)方法對鋼箱內部爆炸的過程進行數值模擬，研究了爆炸過程中鋼箱表面的變形，分析了相同當量炸藥在鋼箱內部不同位置爆炸時鋼箱表面的失效模式。

（1）通過比較實驗與SPH 方法數值模擬獲得的結果，數值模擬與實驗所得的失效變形模式吻合較好，實驗和數值模擬得到的鋼箱表面撓度值相對誤差為3.5%，驗證了本文SPH 方法在模擬鋼箱內部爆炸問題上的有效性；

（2）采用SPH 方法數值模擬鋼箱內部爆炸的過程，得出與實驗一致的變形結果。同時也發現，在鋼箱內部爆炸過程中，鋼箱角隅處出現超壓和應力集中現象，這與實驗吻合較好。

（3）進一步研究了相同當量炸藥在鋼箱內部不同位置爆炸時鋼箱失效破壞情況，發現炸藥在鋼箱角隅處爆炸時鋼箱破壞最嚴重。

以上研究結果可為工事結構的防護設計提供參考。