李科斌, 李曉杰,2, 閆鴻浩, 王小紅, 楊晨琛
(1.大連理工大學(xué) 工程力學(xué)系,遼寧 大連 116024;2.工業(yè)裝備結(jié)構(gòu)分析國家重點(diǎn)實(shí)驗(yàn)室,遼寧 大連 116024)
空中爆炸是爆炸力學(xué)研究中的一個基本現(xiàn)象,深入分析爆炸空氣沖擊波的傳播規(guī)律,特別是凝聚炸藥空中爆炸時的近場特性,對于研究武器試驗(yàn)、爆破工程、炸藥工業(yè)等領(lǐng)域的安全防護(hù)問題有重要意義,如露天爆破、爆炸焊接時沖擊波等有害效應(yīng)的控制,爆炸容器等結(jié)構(gòu)的抗爆防護(hù)設(shè)計等。
眾所周知,空中爆炸時,爆轟產(chǎn)物高速膨脹強(qiáng)烈壓縮周圍的空氣形成空氣沖擊波,而表征空中爆炸近場特性的特征參量包括峰值超壓Δpm、比沖量i、正壓作用時間t+、沖擊波到達(dá)時間t′等(如圖1所示),為計算這些參量,學(xué)者們通過大量試驗(yàn)提出了各自的計算公式[1-8],這些經(jīng)驗(yàn)公式具有不同的適用范圍,相互之間也存在一定的差異。
另一方面,當(dāng)外界空氣環(huán)境(如壓力、溫度、濕度、風(fēng)力等)改變時,上述特征參量的規(guī)律將會隨之發(fā)生變化[9],即滿足薩克斯比例定律[10],Dewey等[11]通過調(diào)節(jié)高度模擬器中的環(huán)境壓力和溫度,測定了多種尺寸裸彭托利特炸藥的爆炸波參數(shù),得出的結(jié)論證實(shí)了薩克斯定律,該定律同樣由Ericsson等[12-13]的實(shí)驗(yàn)所證實(shí);Veldman等[14]研究了環(huán)境壓力對反射沖擊波超壓和沖量的影響;Silnikov等[15]對比分析了高空中和常壓環(huán)境下爆炸沖擊波的參數(shù)變化及對飛行器的破壞影響。另一方面,Jack等[16]在炸藥近旁和模擬高空情況時的試驗(yàn)指出,上述條件下爆炸波的性質(zhì)將會發(fā)生改變,薩克斯定律此時不能用于特征參量的換算,如何準(zhǔn)確而又完整地描述爆炸近場特征參量的變化規(guī)律是本文將要解決的主要問題之一。
圖1 典型的空中爆炸p-t曲線Fig.1 A typical pressure-time curve for air explosion
AUTODYN有限元程序可很好地用于模擬爆炸沖擊波的傳播過程[17-20],本文將基于該平臺建立空中爆炸的球?qū)ΨQ楔形模型,在與標(biāo)準(zhǔn)狀態(tài)經(jīng)驗(yàn)公式對比的基礎(chǔ)上,引入特征比例距離分析不同真空度下空中爆炸近場的特征參量變化規(guī)律,結(jié)論可為高空爆炸時目標(biāo)的毀傷效應(yīng)以及爆破工程的沖擊波防護(hù)提供參考。
本文采用AUTODYN程序中的球?qū)ΨQ楔形模型來模擬自由場中球形裝藥爆炸沖擊波的傳播,建立如圖2所示的計算模型。模型采用多物質(zhì)Euler算法,選用與實(shí)際對應(yīng)的ANFO炸藥,炸藥半徑與網(wǎng)格尺寸之比λ=100,根據(jù)文獻(xiàn)[18]和[19]的分析結(jié)果,所用計算模型和網(wǎng)格尺寸足以保證計算精度,從炸藥邊界開始設(shè)置一系列監(jiān)測點(diǎn)(編號1#~96#)用于特征參量的結(jié)果輸出。
圖2 自由場空中爆炸楔形計算模型Fig.2 Wedge model of free-space air explosion
1.2.1 炸藥
采用AUTODYN標(biāo)準(zhǔn)材料模型庫里的ANFO炸藥,其狀態(tài)方程為經(jīng)典的JWL狀態(tài)方程,具體為
(1)
式中:P為爆壓,V為爆轟產(chǎn)物相對比容(爆轟產(chǎn)物體積與炸藥初始體積之比),E0為炸藥單位體積的初始內(nèi)能,A、B、R1、R2、ω為炸藥的特征參數(shù)(常數(shù)),本文ANFO炸藥各參數(shù)取值如表1所示。
1.2.2 空氣
空氣的狀態(tài)方程近似采用理想氣體狀態(tài)方程,可描述為
表1 ANFO炸藥參數(shù)Tab.1 Parameters of ANFO
(2)
式中:P為空氣壓力,γ為理想氣體絕熱指數(shù),取1.4;ρ為空氣密度;ρ0為標(biāo)準(zhǔn)狀態(tài)下的空氣密度,取1.225×10-3g/cm3;E0為空氣的初始比內(nèi)能,取2.068×105J/m3。本文通過改變空氣的初始密度近似實(shí)現(xiàn)不同的初始壓力,即實(shí)現(xiàn)不同真空度的空氣環(huán)境。
為檢驗(yàn)?zāi)M結(jié)果的可靠性,通過對比經(jīng)驗(yàn)公式計算的結(jié)果分析數(shù)值計算的精度。根據(jù)薩克斯比例定律可知,包含超壓Δpm、比沖量i、正壓時間t+和某些空氣參數(shù)的無量綱群
特定地的是(Rp01/3/E1/3)的單一函數(shù)。其中,R為爆心距離,p0為外界空氣初壓,c0是環(huán)境空氣的聲速,E是炸藥的總能量(通常用裝藥質(zhì)量W代替)。當(dāng)炸藥種類、裝藥密度和空氣狀態(tài)不變時,即為一般形式的爆炸相似定律,也稱霍普金森比例定律,即
基于以上規(guī)律,學(xué)者們通過大量試驗(yàn)各自得到了計算沖擊波峰值超壓Δpm、正壓作用時間t+和比沖量i的經(jīng)驗(yàn)公式,但這些經(jīng)驗(yàn)公式的適用范圍和計算結(jié)果互有差異,以爆炸波峰值超壓為例,圖3給出了不同沖擊波峰壓經(jīng)驗(yàn)公式計算得到的結(jié)果與數(shù)值模擬結(jié)果之間的關(guān)系(R0為裝藥半徑)。
圖3 不同沖擊波峰壓經(jīng)驗(yàn)公式與模擬結(jié)果的關(guān)系
Fig.3 Comparison of simulated peak pressure with calculated results by empirical formulas
從圖3中可以看出:
① 各經(jīng)驗(yàn)公式在R>10~12R0后的計算結(jié)果趨于一致,因此經(jīng)驗(yàn)計算不適用于近場,其適用范圍一般為R>10~12R0;
② 數(shù)值計算結(jié)果與J. Henrych公式最為吻合,在R>6R0時,誤差在10%以內(nèi),在R>12R0時,誤差小于5%;
③ 在裝藥近旁(R<1.6R0),峰壓與距離之間存在其他對應(yīng)關(guān)系。
圖4和圖5則分別給出了利用不同經(jīng)驗(yàn)公式計算得到的正壓作用時間和比沖量。
圖4 不同經(jīng)驗(yàn)公式的正壓作用時間計算結(jié)果Fig.4 Calculated results of positive pressure duration by empirical formulas
可以看出,不同經(jīng)驗(yàn)公式關(guān)于正壓作用時間和比沖量的計算結(jié)果存在很大差別,無法進(jìn)行參照對比,因此這兩個參量的模擬結(jié)果不適合與經(jīng)驗(yàn)公式作比較,文獻(xiàn)[21]中通過數(shù)值計算得到了如圖6和圖7所示的不同密度PETN炸藥的比沖量和正壓作用時間變化曲線,實(shí)線1代表密度1.6 g/cm3的PETN,虛線2密度為0.4 g/cm3,此處Re=(E/pM)1/3,tM=(E/pM)1/3/(pM/ρM)1/2,iM=pMtM。圖中同時給出了本文常壓下ANFO炸藥的模擬結(jié)果,雖然與文獻(xiàn)選用的炸藥種類不同,但兩者在變化規(guī)律上基本一致。
圖5 不同經(jīng)驗(yàn)公式的比沖量計算結(jié)果Fig.5 Calculated results of specific impulse by empirical formulas
圖6 比沖量的數(shù)值計算結(jié)果對比Fig.6 Comparison of simulated specific impulse curves
圖7 正壓作用時間的數(shù)值計算結(jié)果對比Fig.7 Comparison of simulated positive pressure duration curves
通過以上分析和對比可知,本文采用的數(shù)值方法是可行的,利用AUTODYN有限元程序計算空中爆炸的結(jié)果具有較高置信度。
圖8 不同真空度下凝聚炸藥爆炸的峰值超壓曲線
Fig.8 Peak overpressure curves of condensed explosive at different vacuum degrees
通過觀察圖8各曲線發(fā)現(xiàn),爆炸峰壓的變化呈現(xiàn)明顯的區(qū)間性,裝藥近旁和遠(yuǎn)處具有不同的變化規(guī)律,峰壓首先沿p0=0 atm(即完全真空)曲線變化,在某一距離后,爆炸峰壓的衰減即為空氣沖擊波的衰減,這是因?yàn)檠b藥近區(qū)的特征參量受爆轟產(chǎn)物的影響。圖9顯示了1標(biāo)準(zhǔn)大氣壓下距爆炸中心不同距離處壓力隨時間變化的曲線群,在沖擊波傳播的初始階段(4#Gauss點(diǎn)),壓力出現(xiàn)了兩個峰值,分別對應(yīng)爆轟產(chǎn)物-空氣界面前端的沖擊波峰值和界面后端的產(chǎn)物壓力峰值,在最初階段,空氣沖擊波峰壓小于爆轟產(chǎn)物峰壓(圖8中臨界點(diǎn)前的初始段),隨著爆轟產(chǎn)物的急劇膨脹,沖擊波峰壓迅速超過產(chǎn)物峰壓,因此圖8中各真空度峰壓曲線的大部分體現(xiàn)的是空氣沖擊波峰壓的衰減過程。這樣,通過對爆轟產(chǎn)物自由膨脹曲線(p0=0 atm)和空氣沖擊波峰壓衰減曲線分別擬合,就可以得到爆炸場全范圍內(nèi)峰值超壓的表達(dá)式。
圖9 不同距離處壓力-時間變化曲線Fig.9 Pressure-time curve at different distances
從圖9中還可以看出,爆轟產(chǎn)物對壓力的影響逐漸減小,圖中56#監(jiān)測點(diǎn)距離爆炸中心12R0,該處的壓力衰減基本趨于平滑,此處也正是爆轟產(chǎn)物極限膨脹體積的位置,而大部分經(jīng)驗(yàn)公式的適用范圍也正是在極限膨脹體積之后。
為分析外界初始空氣狀態(tài)對爆炸特征參量的影響,我們定義特征比例距離
(3)
式中:fd=(ρ0/ρM)1/3=[(p0/pM)/(T0/TM)]1/3,稱為初始狀態(tài)比例系數(shù),下標(biāo)0表示初始空氣狀態(tài),M表示參考狀態(tài)(pM=0.101 325 MPa,TM=273 K);R為爆心距離,m;W為炸藥TNT當(dāng)量,kg。
圖10給出了以特征比例距離Z表示的不同空氣環(huán)境的峰壓變化規(guī)律,計算結(jié)果與相關(guān)實(shí)驗(yàn)結(jié)果[22]符合良好。需要說明的是,由于本文旨在研究初始壓力的影響,因此計算結(jié)果中空氣的初始溫度都取定為288 K。
圖10 以特征比例距離表示的沖擊波峰壓曲線
Fig.10 Peak overpressure of blast wave with the scaled distanceZas variables
考察圖10中本文計算結(jié)果可以發(fā)現(xiàn),在該坐標(biāo)系中,對于同一炸藥,不同初始空氣壓力的沖擊波峰壓衰減趨于同一變化規(guī)律,擬合后的方程形式為
(4)
式中:Δpm為峰值超壓,MPa;p0為空氣域初始壓力,MPa。對于本文所選炸藥,α=5.565,β=1.078。
分析圖8中p0=0 atm的曲線,爆轟產(chǎn)物自由擴(kuò)散時峰壓的變化規(guī)律可擬合為
(5)
于是,我們就可以得到不同初始空氣壓力下全范圍峰壓的分段函數(shù)
(6)
其中,α、β、λ和ξ由炸藥的性質(zhì)所決定,對于AUTODYN程序提供的ANFO炸藥,上述系數(shù)分別為5.565、1.078、3.902和5.970。Rs為不同初始壓力對應(yīng)的臨界距離,它與真空度的關(guān)系如圖11所示。
圖11 不同真空度的峰壓臨界距離Fig.11 Critical distance vs vacuum degree
擬合后得
(7)
圖12 不同真空度下的比沖量變化Fig.12 Specific impulse curves of various vacuum degrees
為具體分析比沖量極值點(diǎn)的位置,我們以1標(biāo)準(zhǔn)大氣壓時的情況為例,給出了八個相繼時刻下20R0范圍內(nèi)的相對壓力分布(如圖13所示)。
圖13 1 atm下不同時刻的壓力分布曲線Fig.13 Pressure distributions at different times
首先,我們對t=0.085 ms時刻(分布4)進(jìn)行單獨(dú)考察。圖14給出了該時刻下的壓力分布曲線及云圖,圖中產(chǎn)物-空氣界面前是空氣沖擊波陣面,而在爆轟產(chǎn)物內(nèi)部,從界面反射回的稀疏波又從爆心產(chǎn)生了又一次反射,稀疏波尾部的壓力達(dá)到最小值,而主沖擊波后的壓力還沒能正常衰減到稀疏波尾部的狀態(tài),因此在稀疏波尾部就形成了如圖所示的間斷面,該間斷面最初是在運(yùn)用特征線法進(jìn)行爆炸數(shù)值計算中發(fā)現(xiàn)的[23],一些文獻(xiàn)中闡述此間斷面時稱其為“二次沖擊”(secondary shock)[24-25]。
圖14 不同真空度下的比沖量變化Fig.14 Specific impulse curves of various vacuum degrees
繼續(xù)觀察圖13,間斷面在初始階段由于爆轟產(chǎn)物的推動作用由內(nèi)向外移動,直到分布7(間斷面到達(dá)約10R0處)之后才開始向爆心運(yùn)動。圖12中,曲線波谷對應(yīng)圖13中的分布4,該時刻稀疏波尾部的壓力等于初壓p0;而曲線的波峰對應(yīng)圖13中的分布7,該時刻則是間斷面?zhèn)鞑シ较虻霓D(zhuǎn)折點(diǎn),即由向外移動開始向爆心運(yùn)動。其他真空度下的比沖量規(guī)律類似,這種非單調(diào)性變化都反映了爆轟產(chǎn)物中稀疏波尾部間斷面形成和傳播的波動特性。
不同真空度下正壓作用時間的變化規(guī)律如圖15所示,可以發(fā)現(xiàn),正壓作用時間隨比例距離Z呈分段變化規(guī)律,各折點(diǎn)位置同樣與稀疏波尾部的間斷面相關(guān),其中A點(diǎn)對應(yīng)圖13中的分布4,即此時稀疏波尾部壓力等于p0;在此之前,真空度小的曲線始終位于真空度大的曲線上方,此后,各曲線重疊在一起,即可用特征比例距離Z為變量的同一方程來描述不同真空度下的正壓作用時間變化規(guī)律。折點(diǎn)B則對應(yīng)圖13中的分布7,即此時間斷面開始向爆心運(yùn)動;折點(diǎn)C對應(yīng)爆轟產(chǎn)物膨脹的極限體積處,而在BC之間的范圍內(nèi),正壓作用時間與特征比例距離呈線性關(guān)系。
通過程序給定的監(jiān)測點(diǎn)的超壓數(shù)據(jù)可以得到峰壓的時程曲線,即為空氣沖擊波波陣面的運(yùn)動軌跡,如圖16所示;而通過追蹤爆轟產(chǎn)物界面的運(yùn)動軌跡,則可以得到產(chǎn)物界面運(yùn)動的時程曲線,如圖17所示,兩者經(jīng)特征化后的曲線也相應(yīng)地給出,圖中tM含義與前述相同(能量E用質(zhì)量W代替)。
圖16 空氣沖擊波波陣面時程曲線Fig.16 Time-history curves of shock wave front
分析圖16和圖17可知,在常規(guī)坐標(biāo)系中,真空度越大,沖擊波陣面和產(chǎn)物界面運(yùn)動速度越快,且在一定距離后,沖擊波速度各自趨于一恒定值,產(chǎn)物界面速度則趨于0,即爆轟產(chǎn)物膨脹至極限體積,而完全真空狀態(tài)時,爆轟產(chǎn)物膨脹速度為一恒定值(U∞≈6 600 m/s)。不同真空度下產(chǎn)物界面和沖擊波陣面時程曲線在相應(yīng)的特征坐標(biāo)系中都各自歸結(jié)為一條運(yùn)動軌跡,即在以Z表示的系統(tǒng)中,盡管外界的真空環(huán)境不同,但產(chǎn)物界面和沖擊波陣面的運(yùn)動規(guī)律都是相同的。
圖17 爆轟產(chǎn)物界面時程曲線Fig.17 Time-history curves of detonation product interface
不同真空度下爆轟產(chǎn)物的極限膨脹體積與初始空氣壓力的關(guān)系如圖18所示,可以推算,它的垂直漸近線為Y軸,水平漸近線為X軸;而極限體積與臨界距離之間呈線性關(guān)系(圖19),即
(8)
圖18 極限膨脹體積與初始壓力的關(guān)系Fig.18 Ultimate expansion volume vs. initial pressure
圖19 極限膨脹體積與臨界距離的關(guān)系Fig.19 Ultimate expansion volume as a function of critical distance
本文基于有限元程序AUTODYN-2D,建立了空中爆炸的一維球?qū)ΨQ楔形模型,通過改變空氣域的初始參數(shù),得到了不同真空度下爆炸近場特征參量的變化規(guī)律,主要結(jié)論如下:
(1)常壓下的模擬結(jié)果與峰壓經(jīng)驗(yàn)公式符合良好,尤其與J. Henrych最為吻合,在R>12R0時,誤差小于5%,因此數(shù)值計算具有較高的置信度。
(3)不同真空度下比沖量i和正壓作用時間t+的變化具有非單調(diào)性,它們的極值點(diǎn)和折點(diǎn)與稀疏波尾部間斷面密切相關(guān),當(dāng)稀疏波尾部壓力等于p0時,i和t+達(dá)到極小值,而當(dāng)間斷面的傳播方向由外向爆心轉(zhuǎn)變時,i達(dá)到極大值,t+則出現(xiàn)折點(diǎn)。不同真空度在稀疏波尾部壓力等于p0之后具有相同的比沖量和正壓作用時間;此外,不同真空度下的爆轟產(chǎn)物界面和沖擊波陣面時程曲線各自重合,極限體積與臨界距離之間呈線性關(guān)系。