基于LS-DYNA的大空間柱殼結(jié)構(gòu)爆炸波壓力場(chǎng)分析

高軒能，劉 穎，王書(shū)鵬

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

大空間結(jié)構(gòu)在候車(機(jī))樓、影劇院、體育館和博物館等公共建筑中得到廣泛的應(yīng)用，已成為我國(guó)公共建筑的主要形式之一。由于人員活動(dòng)頻繁密集，這類建筑也極易成為恐怖爆炸襲擊的目標(biāo)，造成重大的人員財(cái)產(chǎn)損失和惡劣社會(huì)影響。如2003年12月5日俄羅斯皮亞季戈?duì)査箍耸熊囌镜目植辣?，造成?0人死亡150多人受傷;2004年1月23日我國(guó)鄭州火車站售票廳內(nèi)的炸彈爆炸，導(dǎo)致2人死亡3人受傷;2004年5月9日俄羅斯車臣共和國(guó)格羅茲尼市迪那摩體育場(chǎng)的爆炸，造成32 人當(dāng)場(chǎng)死亡近 100 人受傷［1，2］。因此，對(duì)大空間建筑結(jié)構(gòu)進(jìn)行防爆和抗爆研究，具有迫切的現(xiàn)實(shí)意義。

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

1 爆炸沖擊波及其經(jīng)驗(yàn)計(jì)算公式

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

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

根據(jù)我國(guó)國(guó)防工程設(shè)計(jì)規(guī)范，化爆沖擊波超壓及正壓作用時(shí)間按如下公式計(jì)算［11］:

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

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

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

為簡(jiǎn)化計(jì)算，本文采用線性沖擊波函數(shù)［13］:

式中:ΔP-為沖擊波負(fù)壓，t-為負(fù)壓作用時(shí)間。

2 大空間柱殼結(jié)構(gòu)的ANSYS/LS-DYNA爆炸模擬實(shí)驗(yàn)

2.1 LS-DYNA簡(jiǎn)介

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

2.2 材料與結(jié)構(gòu)模型

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

式中:E0為單位體積炸藥的初始內(nèi)能，V為相對(duì)體積，A，B，R1，R2，Ω 為基于試驗(yàn)的輸入?yún)?shù)。TNT 炸藥為0.5 m×0.5 m×0.5 m 立方體。炸藥的材料輸入?yún)?shù)如表1所示，其中ρ為密度，D為爆速，PCJ為爆壓，V0為初始相對(duì)體積。

表1 炸藥的材料參數(shù)Tab 1.Material parameters of dynamite

空氣采用MAT-NULL材料模型和線性多項(xiàng)式狀態(tài)方程EOS-LINEAR-POLYNOMIAL。即:

式中:E為單位體積內(nèi)能，ρ為質(zhì)量密度，ρ0為參考質(zhì)量密度。線性多項(xiàng)式狀態(tài)方程描述空氣時(shí)遵守Gamma定律，空氣的材料輸入?yún)?shù)如表2所示。

表2 空氣的材料參數(shù)Tab 2.Material parameters of air

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

圖2 大空間鋼柱殼結(jié)構(gòu)計(jì)算模型Fig.2 Model of large-space steel cylindrical reticulated shell

2.3 實(shí)驗(yàn)結(jié)果及其驗(yàn)證

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

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

如果結(jié)構(gòu)是封閉的且在沖擊作用下不破壞，則沖擊波受到阻擋，空氣質(zhì)點(diǎn)會(huì)再次在屋蓋表面區(qū)域堆積，屋蓋表面壓力會(huì)隨之增大，但速度較小如圖5所示。

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

圖5 屋蓋表面處空氣質(zhì)點(diǎn)壓力及速度Fig.5 Air pressure and velocity near the roof

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

從以上模型計(jì)算結(jié)果的分析可知，計(jì)算結(jié)果與理論分析相符，各測(cè)點(diǎn)的沖擊波超壓和空氣質(zhì)點(diǎn)速度等數(shù)據(jù)與經(jīng)驗(yàn)公式的計(jì)算結(jié)果也符合較好，表明模型及參數(shù)選取可靠和正確，計(jì)算結(jié)果可信。

3 數(shù)值模擬結(jié)果與分析

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

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

從圖7、圖8中可以看出，沖擊波本征模態(tài)分布隨結(jié)構(gòu)的矢跨比減小而變得分散，時(shí)間坐標(biāo)的峰值隨矢跨比減小而減小?？缍葹?0 m時(shí)，矢跨比F/B=1/2結(jié)構(gòu)的第1階本征模態(tài)相應(yīng)的超壓幅值時(shí)間曲線與另兩種矢跨比的曲線差別較大，曲線有兩個(gè)波峰，其第2個(gè)波峰較第1個(gè)波峰小很多，曲線的負(fù)壓區(qū)不明顯，而矢跨比F/B=1/4和1/5的超壓幅值時(shí)間曲線出現(xiàn)了3次波峰，且負(fù)壓區(qū)明顯?？缍葹?0 m時(shí)，第1階本征模態(tài)相應(yīng)的超壓幅值時(shí)間曲線均出現(xiàn)兩個(gè)波峰，第1個(gè)波峰隨矢跨比減小而減小，而第2個(gè)波峰隨矢跨比減小而增大。

4 結(jié)論

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

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

［2］高軒能，王書(shū)鵬.大空間柱面網(wǎng)殼結(jié)構(gòu)在爆炸荷載下的動(dòng)力響應(yīng)［J］.振動(dòng)與沖擊，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］李海旺，郭 可，魏劍偉，等.撞擊載荷作用下單層球面網(wǎng)殼動(dòng)力響應(yīng)模型實(shí)驗(yàn)研究［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］周 岱，劉紅玉，李春祥.拉索-網(wǎng)殼結(jié)構(gòu)的動(dòng)力特性和非線性動(dòng)力反應(yīng)［J］.振動(dòng)與沖擊，2002，21(1):7-11.

［7］宋衛(wèi)東，吳開(kāi)騰，寧建國(guó).機(jī)庫(kù)內(nèi)爆炸效應(yīng)的數(shù)值模擬研究［J］.兵工學(xué)報(bào)，2007，28(11):1329-1334.

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

［9］張建勝，武 岳，沈世釗.單層網(wǎng)殼結(jié)構(gòu)風(fēng)振響應(yīng)的模態(tài)耦合效應(yīng)分析［J］.振動(dòng)與沖擊，2006，25(6):39-42.

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

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

［12］Henrych J等，熊建國(guó)譯.爆炸動(dòng)力學(xué)及其應(yīng)用［M］.北京:科學(xué)出版社，1987.

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

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

［15］陳 波.大跨屋蓋結(jié)構(gòu)等效靜風(fēng)荷載精細(xì)化理論研究［D］.哈爾濱:哈爾濱工業(yè)大學(xué)，2006.