爆炸沖擊下試驗臺架數值模擬方法

莫昌鋒，周云波，付條奇，張 明，劉粟濤，孫曉旺

(南京理工大學 機械工程學院， 南京 210094)

1 引言

現代戰爭中軍用防護型車輛在執行任務的時經常面臨著地雷和簡易爆炸裝置(improvised explosive device，IED)的威脅，其產生的爆炸沖擊波不僅使車輛結構發生破壞，同時還會對車內乘員造成嚴重的生理損傷[1]。隨著計算機數值仿真技術的成熟，有限元分析已經成為研究軍用防護型車輛爆炸沖擊響應的重要工具。由于地雷和IED發生爆炸時會在短時間內產生巨大的能量，因此在車輛防護領域中，研究準確高效的爆炸沖擊數值模擬算法是這些年來的熱點和難點問題。

LS-DYNA使用顯示時間積分的高度非線性瞬態動態進行有限元分析，可以很好地模擬爆炸沖擊過程。國內外學者在應用LS-DYNA進行爆炸沖擊數值仿真上做了大量的研究工作，如畢程程等[2]研究了初始體積分數法在爆炸模擬中的應用，并將該方法結合ALE流固耦合算法在爆炸算例中進行可行性驗證，得出了初始體積分數法定義的炸藥在爆炸模擬很好地保證計算精度。曾愛等[3]對防爆油箱進行爆炸沖擊響應數值模擬研究，分析了爆炸環境下燃油箱和燃油的動態響應過程。Thanh等[4]在LS-DYNA中分別應用了經驗爆炸算法、ALE算法和拉格朗日-歐拉耦合算法分析爆炸載荷下鋼板的受力情況。Slavik[5]提出了一種新的爆炸仿真方法，在拉格朗日結構周圍的空氣進行ALE建模，而對于遠場爆炸則只需在ALE空氣域的外表面與經驗爆炸方程結合即可，該數值仿真方法擴展了LS-DYNA分析爆炸沖擊新場景的能力。Samiee等[6]利用LS-DYNA的顯示有限元分析，研究了聚脲對爆炸載荷下雙層聚脲鋼板性能的影響。

近年來，在對爆炸沖擊問題分析的數值模擬上，由于ALE算法結合了拉格朗日算法和歐拉算法的優點，能很好地體現流體的運動規律，同時又能解決網格畸變的問題，因此ALE算法用來模擬爆炸沖擊過程已經較為成熟且精度較高[7]。在進行爆炸分析時，ALE算法需要建立大量的空氣和土壤單元，這導致ALE的計算時間比較長，造成計算效率不高。為了提高ALE算法的計算效率，LS-DYNA在ALE算法的基礎上開發了一種新的算法，利用結構化網格來加快計算時間，這種算法稱為結構化ALE算法(S-ALE)。王明振等[8]應用S-ALE流固耦合算法對飛機水上迫降進行數值分析，結果表明相對于傳統的ALE耦合算法，S-ALE流固耦合算法較好地避免了流體的滲漏。Kurtolu[9]分別應用ALE算法和S-ALE算法模擬了鋼結構艙體在爆炸沖擊下的響應，研究表明采用S-ALE求解器時占用內存更小，并且在求解相同的模型下需要的時間更短。

本文中主要研究內容是用傳統ALE算法和S-ALE算法對爆炸沖擊試驗臺架進行數值模擬，比較這2種算法下鋼板的最大殘余變形和內能。最后進行臺架實爆試驗，驗證并比較這2種算法的精確度，通過仿真和試驗結果對比分析，S-ALE算法與試驗結果更接近，計算時間更短。

2 爆炸仿真算法研究和有限元模型

2.1 爆炸仿真算法

傳統的ALE算法和S-ALE算法的區別在于描述流場作用區域時是采用非結構化網格還是結構化網格。采用ALE算法分析時，LS-DYNA的ALE求解器會先讀取網格元素和節點信息，然后根據這些信息輸入來建立流場作用區域的非結構化網格，這個過程將會消耗大量的讀取時間。S-ALE算法則將流場作用區域通過三維坐標來定義，LS-DYNA的S-ALE求解器根據三維坐標定義的幾何信息直接生成結構化三維網格。結構化網格的特點就是網格內的節點與相鄰的點連接關系是固定不變的，并且隱含在生成的網格中[10]，這避免了求解模型時讀取網格元素和節點信息的過程，加快數值計算的速度，并且使得模型占用的內存較少。同時結構化網格描述單元和節點的連接情況時簡單明了，簡化了ALE耦合時的搜索算法。

ALE算法和S-ALE算法在處理爆炸沖擊問題時，都需要考慮爆炸流場與目標結構之間的耦合作用。基于多物質單元耦合模擬爆炸沖擊過程是一種國內外公認成熟的仿真算法，該算法能很好地預測實爆試驗的結果[11]。在仿真過程中，爆炸沖擊波通過土壤、空氣等介質傳播，最后作用到目標結構，這個過程中流固耦合的控制主要通過罰函數來實現[12]。罰函數算法是檢測到接觸面發生穿透的時候，通過節點施加邊界力，這個邊界力定義為：

fi=-kidi

(1)

其中：fi為邊界力；di是接觸面間的穿透距離；ki是接觸剛度。如圖1所示，在接觸邊界上，罰函數法被認為每一個耦合點創建一個虛擬彈簧，接觸剛度的定義為：

(2)

其中：fsi是接觸面耦合比例因子；Kv、Al、Vl分別是體積模量、流體網格的斷面面積和流體單元的體積。

2.2 有限元模型的建立

為了驗證某試樣鋼板的抗爆炸性能，設計了埋爆炸藥沖擊試驗臺架，見圖2。試驗臺架由底部支架、鋼板固定支架以及配重塊組成，材料都為Q235，其材料參數見表1所示。試樣鋼板的總體尺寸為1 500 mm×1 500 mm×20 mm，四周用高強度螺栓固定在試驗臺架上。根據爆炸試驗要求，采用TNT圓柱形炸藥，高徑比為1/3，當量為8 kg，炸點位于試樣鋼板的正中間下方，炸藥填埋深度100 mm。

圖1 流固耦合邊界力示意圖

圖2 試驗臺架示意圖

表1 Q235的材料參數

在對臺架以及試樣鋼板進行劃分網格時，網格的基本尺寸為10 mm，除了配重塊為實體單元，其余部分均為殼單元。采用ALE流固耦合的方法模擬爆炸沖擊波作用于試樣鋼板，因此需要建立空氣和土壤網格。空氣網格區域只需要包含試驗臺架即可，空氣和網格用實體單元建模，總尺寸為2 862 mm×2 600 mm×1 790 mm，土壤的高度390 mm，有限元模型見圖3所示。

圖3 試驗臺架有限元模型示意圖

空氣區域采用Mat-Null材料模型與線性多項式狀態方程進行定義[13]:

P=C1+C2μ+C3μ2+(C4+C5μ+C6μ2)E

(3)

其中：P表示氣體壓力；Ci(i=1、2、3、4、5、6)為多項式方程系數；E為單位體積初始內能。這些參數對應的數值分別為：C0=-0.1；C1=C2=C3=C6=0；C4=C5=0.4；E0=0.253 mJ/mm3；V0=1[14]。土壤為超塑性材料模型，通過Mat_Soil_And_Foamble_Fail定義。

在LS-DYNA中，設置爆轟物氣體的狀態方程和材料參數來建立爆轟物的有限元模型。Jones、Wilkins和Lee(JWL)提出的半經驗狀態方程(EOS)描述了爆轟物氣體的壓力、能量和體積的關系，這個表達式為：

(4)

其中：A、B、R1、R2、ω是常數；E是單位體積能量；V是爆炸沖擊過程氣體產物相對于初始爆炸狀態的相對體積。在LS-DYNA中通過EOS_JWL卡片來設置EOS參數，表2是TNT炸藥的JWL的EOS參數[15]。

表2 JWL參數

通過MAT_HIGH_EXPLOSIVE_BURN設置TNT的物理參數，表3是TNT的物理參數[15]。

表3 TNT物理參數

由于使用對炸藥用網格建模時，難以保證在有限元分析過程中的精度。初始體積法在處理流體模型是有著獨特的優勢，突破物質初始定義網格化的局限，并且在分析過程中能夠保證計算精度[2]。因此在設置好炸藥的狀態方程和物理參數之后，通過初始體積法來定義炸藥的形狀，用 INITIAL_VOLUME_FRACTION_GEOMETRY關鍵字定義炸藥的形狀。

此外，在ALE流固耦合中，主要是用罰函數法來描述耦合關系，用CONSTRAINED_LAGRANGE_IN_SOLID這個關鍵字卡片來設置耦合關系。為了模擬無限空氣和土壤域，在空氣與土壤網格的外表面設置無反射邊界條件。整個有限元模型單元數為664 207，節點數為685 323。

同時使用結構化ALE(S-ALE)算法進行分析，這是LS-DYNA開發的一種新算法。S-ALE算法利用結構化網格的特征來加快求解模型的速度，提升數值分析的效率。S-ALE流固耦合方法很多控制卡片的設置和ALE流固耦合方法是相類似的。在本研究中，針對空氣域和土壤域不再使用實體單元建模，而是重新在LS-DYNA中建立空氣和土壤的結構化網格，分別通過ALE_STRUCTURED_MESH和ALE_STRUCTURED_MESH_CONTROL_POINTS這兩個新的關鍵字卡片來控制，并且通過這個設置在求解模型時調用S-ALE求解器，除了這2個卡片設置外，其他大部分設置與上述的ALE流固耦合方法是一樣的。

3 數值仿真結果與討論

3.1 ALE算法仿真結果

爆炸沖擊下臺架的試樣鋼板的變化如圖4所示。在爆炸沖擊波的作用下，0.7 ms鋼板發生變形向上凸起，并且在4.2 ms時刻鋼板向上凸起的變形最大。圖5是鋼板的最大變形的位移云圖和最大變形點位移曲線，可以看出，在4.2 ms之后，鋼板的向上凸起的變形量略有下降，這是由于鋼板存在彈性變形和塑性變形，彈性變形部分后面逐漸恢復。在8 ms到12.5 ms之間，鋼板的最大變形量逐漸趨于穩定，因此仿真計算時間到12.5 ms。

圖4 ALE算法仿真結果云圖

圖5 ALE算法鋼板仿真結果云圖

3.2 S-ALE算法仿真結果

S-ALE算法仿真結果見圖6所示。從仿真動畫來看，試樣鋼板變形的大致趨勢和ALE算法相類似。在0.7 ms爆炸沖擊波使得鋼板開始發生變形，在4.2 ms鋼板變形量達到最大。鋼板變形云圖與最大變形點位移曲線見圖7所示，8～12.5 ms鋼板變形量趨于穩定。

圖6 S-ALE算法仿真結果云圖

圖7 S-ALE算法試驗鋼板仿真結果

3.3 2種數值模擬方法對比

1) 通過鋼板最大變形量曲線可以看出，ALE算法計算得到的試樣鋼板最終的最大變形量為214.260 mm，S-ALE算法計算得到的試樣鋼板最終的最大變形量為233.171 mm。

2) 圖8顯示2種算法在爆炸沖擊波作用到鋼板時的流固耦合控制，可以看出，ALE算法在仿真過程中出現了爆炸沖擊波流體滲漏，而S-ALE算法則在仿真過程中流體得到有效的控制。

圖8 2種算法流固耦合控制云圖

3) 圖9是2種算法試樣鋼板的內能曲線，從曲線可以看出，2種算法總體趨勢類似，但S-ALE算法的試樣鋼板內能更大，這是因為仿真中S-ALE算法能進行有效的流固耦合控制，使得爆炸沖擊波作用在鋼板上的力更大，鋼板產生更大的變形。

圖9 鋼板內能曲線

4) 2種算法在計算時間上，S-ALE算法需要的時間更短，計算效率更高。

4 臺架試驗

4.1 試驗內容與布置

本次沖擊臺架試驗目的是測試某型號鋼板抗爆炸沖擊性能，并通過試驗結果驗證前期數值分析的精度。

由于爆炸試驗歷時短，可重復性差且危險性高，因此在試驗之前必須做好充足的試驗準備。根據試驗標準的要求，制作的試驗臺架如圖10所示。試驗前還需測量臺架正面與側面的水平度，如圖11所示。試驗臺架除了試樣鋼板外均由Q235材料制造而成。鋼板與臺架用M20的高強度螺栓進行連接固定，兩根支撐橫梁由安裝臺架外部貫穿放置于支撐座上，靶板下表面距地面400 mm。試驗臺架頂部增加配重，配重均勻，臺架重約1.5 t，頂部配重塊約7.5 t。炸藥由TNT制成，當量為8 kg，圓柱狀，埋藏深度為100 mm。炸藥放置在標定的起爆位置之后，起爆體與設置在200 m以外的起爆裝置連接，最后炸藥在距離鋼板500 mm中心處被引爆。

圖10 爆炸沖擊試驗臺架實物圖

圖11 正面與側面水平度測量現場圖

爆炸后臺架的整體狀況如圖12和圖13所示。臺架一側陷入沙坑，橫梁出現折彎，配重塊未脫離；試樣鋼板與臺架的連接螺栓全部斷裂脫離，螺栓孔變形；鋼板背凸形成鼓包，未出現裂紋。

圖12 試驗后臺架實物圖

圖13 試驗后鋼板實物圖

試驗后鋼板發生變形，無裂紋產生，為了測出鋼板的變形量，因此采用3D掃描儀成像系統測量鋼板爆炸后的變形量，測量面為背爆面，掃描成像見圖14所示。經過測量，鋼板的變形量為231.047 mm。

圖14 鋼板3D掃描成像圖

4.2 仿真與試驗對比分析

根據試驗和仿真的分析，對比2種算法與試驗的鋼板變形情況見圖15所示，可以看出，S-ALE算法計算的結果與試驗值更接近。2種算法計算結果與試驗實測結果見表4所示。

圖15 鋼板變形量曲線

表4 仿真結果與試驗結果

ALE算法計算的鋼板變形量與試驗值的誤差為7.27%，而S-ALE算法計算結果的誤差為0.92%；S-ALE算法計算的時間比ALE算法的少4 593 s；對比分析可知，S-ALE算法的計算精度更高，計算時間更短。

5 結論

1) 采用了ALE算法和S-ALE算法模擬炸藥埋爆沖擊臺架試驗，對比了2種算法在仿真過程中的差異。采用ALE算法模擬時，在仿真過程中出現了爆炸沖擊波流體滲漏，而S-ALE算法則很好地避免了流體的滲漏，在模擬爆炸沖擊波作用到鋼板過程中，進行了有效的流固耦合控制，使S-ALE算法計算結果與試驗結果更加接近，計算精度更高。

2) 對比分析2種算法在鋼板變形量和內能的差異，發現在前2 ms內的計算結果吻合度較好，而后面由于ALE算法的流體產生滲漏，導致計算結果產生差異。

3) 與ALE算法相比，S-ALE算法仿真結果與試驗實測結果更吻合，計算精度更高，計算時間更短。因此，S-ALE算法在模擬爆炸沖擊時展現出較大的應用潛力。