基于LS-DYNA的大空間柱殼結構爆炸波壓力場分析

高軒能，劉 穎，王書鵬

(華僑大學 土木工程學院，福建泉州 362021)

大空間結構在候車(機)樓、影劇院、體育館和博物館等公共建筑中得到廣泛的應用，已成為我國公共建筑的主要形式之一。由于人員活動頻繁密集，這類建筑也極易成為恐怖爆炸襲擊的目標，造成重大的人員財產損失和惡劣社會影響。如2003年12月5日俄羅斯皮亞季戈爾斯克市車站的恐怖爆炸，造成近50人死亡150多人受傷;2004年1月23日我國鄭州火車站售票廳內的炸彈爆炸，導致2人死亡3人受傷;2004年5月9日俄羅斯車臣共和國格羅茲尼市迪那摩體育場的爆炸，造成32 人當場死亡近 100 人受傷［1，2］。因此，對大空間建筑結構進行防爆和抗爆研究，具有迫切的現實意義。

對于大跨空間結構在爆炸作用下的研究，國外迄今僅見Raftoyiannis［3］關于多點懸吊屋蓋在爆炸作用下的動力性能研究的1篇論文，國內的研究也不多，尚處于探索階段。李海旺等［4］對單層球殼結構在球頂集中沖擊荷載作用下的動力反應及穩定性進行了試驗與數值分析。范峰等［5］對K8型單層球殼在球頂集中沖擊荷載作用下的動力特性、失效模式和失效機理進行了數值模擬計算，建立了桿件破壞形式與網殼失效模式之間的對應關系。周岱等［6］對拉索-網殼結構的動力特性和非線性動力反應進行了研究。宋衛東等［7］對機庫內爆炸效應進行了3維數值模擬分析，計算了沖擊波超壓峰值與爆炸點距離的關系。高軒能等［8］對大跨空間柱殼結構在爆炸荷載下的動力效應及防爆措施等進行了研究，探討了爆炸沖擊波在殼體頂部的聚集效應和對結構動力性能的影響以及泄爆方法等。但對大空間結構進行爆炸作用下的動力響應研究，首先需要確定爆炸沖擊波壓力場在結構上的分布，因大空間結構自振頻率分布密集、振型復雜和計算量大等困難，目前的研究尚少。為此，本文基于有限元程序ANSYS/LS-DYNA，對大空間柱殼結構進行爆炸模擬實驗，在確定合適的模型及參數后對其進行爆炸作用下的動力響應計算。引入POD法［9］，利用提取的結構表面有限測點的沖擊波超壓數據，分析大空間柱殼結構的爆炸沖擊波壓力場分布以及結構幾何參數等因素對沖擊波分布的影響。

1 爆炸沖擊波及其經驗計算公式

圖1 爆炸沖擊波示意圖Fig.1 Blast shock wave sketch map

化學炸藥爆炸后瞬間便轉變成高壓(104MPa～3×104MPa)和高溫(3.5 ×103℃ ～4×103℃)狀態的氣態爆轟產物。氣體急劇膨脹并迫使周圍的空氣離開它原來占據的位置，于是在此氣體的前沿便形成一壓縮空氣層，即爆炸沖擊波(如圖1 所示)［10］。爆炸的瞬時能釋放巨大的能量并產生各種破壞效應，但對建筑物破壞力最強、破壞區域最大的還是爆炸沖擊波的破壞效應。從一些典型的恐怖爆炸事件看，建筑物的爆炸破壞情況主要與爆炸沖擊波的大小、持續時間、爆炸源的位置和建筑物本身的建筑結構特性等因素有關，還與恐怖爆炸襲擊的方式等其它諸多因素有關。由于爆炸過程的短暫性和復雜性，爆炸沖擊波壓力很難通過簡單的理論分析得到。理論與實驗研究發現，炸藥在空氣中爆炸形成的空氣沖擊波，陣面上的壓力并不取決于裝藥的絕對尺寸和質量，而完全取決于離爆炸地點的距離與裝藥半徑的比值、該炸藥的比能和周圍空氣的壓力?？諝鉀_擊波的峰值超壓ΔPf是指沖擊波陣面上的峰值壓力Pf減去空氣中的大氣壓力P0(標準大氣壓)，即:

根據我國國防工程設計規范，化爆沖擊波超壓及正壓作用時間按如下公式計算［11］:

式中:W為TNT當量，單位kg;R為計算點距爆心的距離，單位m;H為炸藥爆炸時離地面的高度，單位m;ΔPf為爆炸沖擊波超壓峰值，單位MPa，t+為沖擊波超壓正壓作用時間，單位s。

當炸藥在地面爆炸(按剛性地面考慮)時，由于地面的阻擋，空氣沖擊波不是向整個空間傳播，而只向一半無限空間傳播，被沖擊波帶動的空氣量也減少一半，可看作兩倍的裝藥在無限空間爆炸，此時，爆炸沖擊波超壓計算公式為:

亨利奇(Henrych)［12］計算公式則為:

為簡化計算，本文采用線性沖擊波函數［13］:

式中:ΔP-為沖擊波負壓，t-為負壓作用時間。

2 大空間柱殼結構的ANSYS/LS-DYNA爆炸模擬實驗

2.1 LS-DYNA簡介

LS-DYNA程序是由 Hallquist于1976年在美國Lawrence Livermore National Laboratory主持開發完成，是功能齊全的通用大型非線性程序。LS-DYNA的分析功能非常強大，有140多種材料模式庫，包括工程材料、炸藥、推進劑、粘性流體及自定義材料等，以及各種單元模型和50多種接觸方式。它以Lagrange算法為主，兼有ALE和Euler算法;以顯示求解為主，兼有隱式求解功能;以結構分析為主，兼有熱分析、流體-結構耦合功能;以非線性動力分析為主，兼有靜力分析功能。LS-DYNA獨特的算法能十分準確地分析工程結構承受高速碰撞、爆炸、金屬成形等高速高度非線性問題，在工程結構抗爆特性研究、抗爆設計以及爆破拆除等方面有著非常重要的用途［14］。

2.2 材料與結構模型

炸藥采用MAT-HIGH-EXPLOSIVE-BURN材料模型和JWL狀態方程。狀態方程定義的爆炸沖擊壓力為［14］:

式中:E0為單位體積炸藥的初始內能，V為相對體積，A，B，R1，R2，Ω 為基于試驗的輸入參數。TNT 炸藥為0.5 m×0.5 m×0.5 m 立方體。炸藥的材料輸入參數如表1所示，其中ρ為密度，D為爆速，PCJ為爆壓，V0為初始相對體積。

表1 炸藥的材料參數Tab 1.Material parameters of dynamite

空氣采用MAT-NULL材料模型和線性多項式狀態方程EOS-LINEAR-POLYNOMIAL。即:

式中:E為單位體積內能，ρ為質量密度，ρ0為參考質量密度。線性多項式狀態方程描述空氣時遵守Gamma定律，空氣的材料輸入參數如表2所示。

表2 空氣的材料參數Tab 2.Material parameters of air

結構模擬計算模型以公共大跨度建筑中常用的大空間柱殼結構為對象，應用LS-DYNA有限元程序建立模擬大空間結構室內爆炸的計算模型。大空間結構設為三向網格單層柱面鋼網殼，如圖2所示，結構尺寸為40 m×20 m×15 m，矢高為3 m，ALE空間尺寸為42 m×21 m×19 m，矢跨比F/B=1/5，周邊固支，支座下圍護墻體高H=10 m。網殼結構的桿件均為Q235鋼無縫鋼管，焊接球節點，其中軸向直桿截面為φ89 mm×4 mm，斜桿截面為φ140 mm×6 mm。為簡化計算，利用對稱性取1/2模型進行分析，除對稱面外，邊界面為透射邊界以模擬在無限區域的爆炸。本文著重考慮人體或車載炸藥等對大空間結構的爆炸襲擊，爆炸點高度取為離地面高1.2 m。爆炸荷載下不考慮屋面圍護材料等質量，僅考慮桿件及節點自重。計算模型中，屋面殼體結構、地面和墻體視為剛體，以便考慮空氣沖擊波壓力場分布。考慮5種物質的相互作用，選用3種單元?？諝夂驼ㄋ幘捎昧骟w實體Solid164單元，地面和四周圍護墻體采用Shell163單元，結構桿件采用Link160單元。

圖2 大空間鋼柱殼結構計算模型Fig.2 Model of large-space steel cylindrical reticulated shell

2.3 實驗結果及其驗證

為驗證模型和材料參數取值的可靠性和正確性，應用LS-DYNA對炸藥在無限空中爆炸進行了數值模擬，空氣尺寸為12 m×12 m×12 m，炸藥尺寸為0.2 m×0.2 m ×0.2 m，位于空氣中心。算法采用多物質ALE算法，有限元網格大小取為0.1 m ×0.1 m ×0.1 m，計算模型共劃分8節點六面體單元216 000個，其中炸藥為1個單元，空氣為215 999個單元。邊界面為透射邊界。模擬計算結果與經驗公式［11，12］計算結果的比較如圖3所示。圖3表明，本文模型得到的沖擊波超壓與3個經驗公式結果都很接近，而與亨利奇公式［12］的結果最為相符，在超壓持續時間上，則本文模擬計算結果與亨利奇公式幾乎完全吻合，同時還可看出，我國設計規范公式與亨利奇公式計算結果基本符合，略為偏大。

另一方面，炸藥爆炸后，因爆炸點距離地面較近，沖擊波碰撞地面，地面處空氣質點速度驟然變為零，使質點急劇堆積，壓力和密度驟然升高，達到一定程度后，質點就要向相反方向運動，即形成反射沖擊波。反射沖擊波超壓大大增強，所以大空間結構的一部分構件會受到沖擊波的兩次沖擊作用，并且第2次對結構的沖擊作用遠遠大于第1次。圖4(a)為爆炸點上方h=13 m處的沖擊波壓力曲線，可以看出該點受沖擊波2次沖擊作用，第1次為入射沖擊波，根據文獻［12］該處的入射沖擊波超壓為0.182 MPa(亨利奇公式)和0.205 MPa(我國規范公式)，本文有限元結果為0.217 MPa，兩者分別相差19%和6%，表明本文結果與經驗公式符合較好。第2次為地面反射沖擊波，由圖4(b)質點加速度圖可以看出，該處質點在入射沖擊波作用后，加速度由最大變為零再到負值，隨后受到反射沖擊波的作用，加速度再次升到最大，由于反射沖擊波較入射沖擊波大，質點第2次獲得的加速度也會比第1次大(3.25倍)，之后加速度從最大值變小，期間沖擊波壓力升到最大，然后恢復到標準大氣壓。

如果結構是封閉的且在沖擊作用下不破壞，則沖擊波受到阻擋，空氣質點會再次在屋蓋表面區域堆積，屋蓋表面壓力會隨之增大，但速度較小如圖5所示。

圖6為馬赫反射區水平距爆炸點10 m處的沖擊波壓力曲線，可以看出，該處只承受馬赫反射波1次沖擊作用。參考文獻［11］，地爆超壓仍采用文獻［12］的公式計算，TNT當量分別取為1W和2W，計算結果分別作為地爆超壓的上下限值，分別為0.271 MPa和0.447 MPa，本文結果為 0.341 MPa。

圖5 屋蓋表面處空氣質點壓力及速度Fig.5 Air pressure and velocity near the roof

圖6 馬赫反射區距爆心10 m處的空氣壓力Fig.6 Air pressure in Mach reflection region(10 m away from the explosion center)

從以上模型計算結果的分析可知，計算結果與理論分析相符，各測點的沖擊波超壓和空氣質點速度等數據與經驗公式的計算結果也符合較好，表明模型及參數選取可靠和正確，計算結果可信。

3 數值模擬結果與分析

對于圖2所示大空間柱殼結構，由于炸藥距離地面較近，爆炸后沖擊波受到地面的阻擋，反射沖擊波超壓大大增強，因此，結構的一部分構件會受到沖擊波的2次沖擊作用，并且第2次的沖擊作用會遠遠大于第1次。由于作用于結構表面各點的爆炸荷載時間不同步，大小不相同，使結構表面各點上的沖擊波壓力形成了時空差異。沖擊波壓力場的分布是進行大空間結構防爆和抗爆研究的基礎，為克服計算上的困難，本文引入本征正交分解法(POD法)［15］，提取結構表面n個點處的沖擊波超壓數值并將其進行POD分解，即可得到結構的爆炸沖擊波壓力分布圖。

按前述建立的有限元計算模型對不同跨度和矢跨比的大空間柱殼結構進行爆炸模擬計算后，對爆炸沖擊波超壓進行POD法求解，得到的沖擊波壓力場第1階本征模態分布及其超壓幅值-時間曲線如圖7和圖8所示。圖中x軸表示結構跨度方向(0 m～20 m或40 m)，y軸表示結構長度方向(0 m～20m，對稱取半計算)，結構高度H=10 m，長度L=40 m，炸藥0.5 m ×0.5 m×0.5 m。對于跨度B=20 m，矢跨比F/B分別取 1/2、1/4、1/5，對于跨度B=40 m，矢跨比F/B分別取 1/4、1/5、1/8，其它參數如前所述。

從圖7、圖8中可以看出，沖擊波本征模態分布隨結構的矢跨比減小而變得分散，時間坐標的峰值隨矢跨比減小而減小。跨度為20 m時，矢跨比F/B=1/2結構的第1階本征模態相應的超壓幅值時間曲線與另兩種矢跨比的曲線差別較大，曲線有兩個波峰，其第2個波峰較第1個波峰小很多，曲線的負壓區不明顯，而矢跨比F/B=1/4和1/5的超壓幅值時間曲線出現了3次波峰，且負壓區明顯?？缍葹?0 m時，第1階本征模態相應的超壓幅值時間曲線均出現兩個波峰，第1個波峰隨矢跨比減小而減小，而第2個波峰隨矢跨比減小而增大。

4 結論

綜上所述，基于有限元軟件ANSYS/LS-DYNA建立的大空間柱殼結構爆炸作用數值模擬計算模型，與Henrych經驗公式計算結果符合較好，表明計算模型和材料參數取值正確而適用，可用于大空間柱殼結構的爆炸動力響應分析。引入POD法，利用提取的結構表面有限測點的沖擊波超壓數據，解決了沖擊波荷載的時空差異性及結構表面壓力場分布問題。通過對不同跨度和不同矢跨比的大空間柱殼結構在爆炸荷載作用下的沖擊波壓力場的數值模擬計算及分布規律的研究，可知結構的矢跨比等幾何參數變化對大空間柱殼結構的沖擊波壓力場分布影響很大，封閉大空間結構的頂蓋處可產生沖擊波壓力場積聚現象，不利于結構的防爆和抗爆，所得結果給大空間柱殼結構在爆炸荷載作用下的動力特性和結構的抗爆設計提供了基礎。需要說明的是，由于所探討問題的復雜性，本文研究僅僅是初步的，對于大空間柱殼結構在爆炸荷載作用下壓力場分布的精確分析，尚需考慮材料的非彈性和結構的大變形影響，做更深入的研究。

［1］錢七虎.反爆炸恐怖安全對策［M］.北京:科學出版社，2005.

［2］高軒能，王書鵬.大空間柱面網殼結構在爆炸荷載下的動力響應［J］.振動與沖擊，2009，28(10):68-73.

［3］Raftoyiannis I G，Spyrakos C C，Michaltsos G T.Behavior of suspended roofsunderblastloading［J］. Engineering Structures，2007，29:88-100.

［4］李海旺，郭 可，魏劍偉，等.撞擊載荷作用下單層球面網殼動力響應模型實驗研究［J］.爆炸與沖擊，2006，26(1):39-45.

［5］ Fan F，Wang D Z，Zhi X D，et al.Failure modes for singlelayer reticulated domes under impact loads［J］.Transactions of Tianjin University，2008，14:545-550.

［6］周 岱，劉紅玉，李春祥.拉索-網殼結構的動力特性和非線性動力反應［J］.振動與沖擊，2002，21(1):7-11.

［7］宋衛東，吳開騰，寧建國.機庫內爆炸效應的數值模擬研究［J］.兵工學報，2007，28(11):1329-1334.

［8］高軒能，王書鵬，江 媛.爆炸荷載下大空間結構的沖擊波壓力場分布及泄爆措施研究［J］.工程力學，2010，27(4):226-233.

［9］張建勝，武 岳，沈世釗.單層網殼結構風振響應的模態耦合效應分析［J］.振動與沖擊，2006，25(6):39-42.

［10］張守中.爆炸基本原理［M］.北京:國防工業出版社，1988.

［11］ 李翼祺，馬素貞.爆炸力學［M］.北京:科學出版社，1992.

［12］Henrych J等，熊建國譯.爆炸動力學及其應用［M］.北京:科學出版社，1987.

［13］林大超，白春華，張 奇.空氣中爆炸時爆炸波的超壓函數［J］.爆炸與沖擊，2001，21(1):41-46.

［14］ LS-DYNA keyword user's manual［Z］.Version 971，LSTC，March，2005.

［15］陳 波.大跨屋蓋結構等效靜風荷載精細化理論研究［D］.哈爾濱:哈爾濱工業大學，2006.