周豐峻,鄭磊,孫云厚,張鏖,夏明
(1. 中國人民解放軍軍事科學院國防工程研究院,北京 100850;2. 中國人民解放軍92611部隊,廣州 510715)
點源爆炸理論是研究球形裝藥爆炸沖擊波最簡單的模型之一,其假定爆轟產物質量無限小,而裝藥爆炸釋放的能量有限,忽略運動中靜止空氣的影響,因此,不存在任何一個可以表征爆炸的初始線尺度參數,任何一個物理量在空間上的分布隨時間的變化是自相似的,即離爆心不遠處沖擊波后的氣體運動滿足自相似條件 [1]。Taylor [2]、Von Neumann [3]、Bach和Lee [4]以及Sedov [5]各自獨立得到了點源強爆炸自模擬的理論解,與實驗比對可以很好地描述爆炸初始階段的流場特征 [6]。本文基于點源爆炸理論,對平面、柱面以及球面的一維變幾何對稱面(面源、線源和點源)裝藥強爆炸解的自模擬解流場進行了分析。
爆炸沖擊波的傳播是一個非常復雜的不定常過程,但滿足質量守恒、動量守恒和能量守恒關系。由于爆炸力學所考慮的運動變化快,以至流線之間不會發生明顯的動量和能量遷移,因此可以忽略粘性和熱傳導 [7]。
1. 歐拉坐標下的控制方程
歐拉坐標下流場的狀態量主要包括壓力P、密度ρ、質點流動速度u、內能E以及溫度T,這些量值均為自變量x, y , z , t的函數。
質量守恒:

動量守恒:

式(2)中,F為質量力。
能量守恒:

物態方程:

對于一維對稱TNT裝藥爆炸可將介質的運動看作是一維絕熱運動,即介質中各量值只依賴時間和某單一幾何坐標而變化,去掉質量力,即F = 0,可將其在平面直角坐標系、柱坐標系和球形坐標系下建立的歐拉坐標下的控制方程組統一寫為:,u , v , w分別為

式(5)中,r為歐拉位置坐標;t為時間;N為常數,等于0,1和2時分別對應一維平面對稱、一維柱面對稱和一維球面對稱情形。
2. 拉格朗日坐標下的控制方程
為了提高沖擊波陣面局部計算精度,可將方程組寫為拉格朗日形式,控制方程為:

再補充一個運動關系式:

式(6)中,R為拉格朗日位置坐標;u為質點速度;P為壓力;E為內能;ρ為密度;ρ0為爆炸波前初始空氣密度;t為時間;N為相應平面、柱面和球面坐標指數。
為了計算方便,引入質量坐標方程式:

則可將拉格朗日坐標方程組轉化為無量綱形式:

式(9)中,N為常數,等于0,1和2時分別為包含一維平面、柱面和球面的一體無量綱方程組。
1. 真實空氣無量綱形式狀態方程
試驗表明,當溫度大于等于2 000 K或沖擊波超壓大于等于4 MPa時,實際空氣不再滿足理想氣體狀態方程 [7]。隨著溫度與壓力的不斷升高,空氣中各種分子之間會發生化學反應,如氮氧分子開始離解,經反應生成NO和NO2。空氣組分在不同溫度下是不同的,當溫度高于8 000 K后,空氣中各種原子會發生電離反應,空氣組分變得更加復雜;當溫度達到300 000 K后,任何密度下氣體的電離過程都可以完成。在溫度升高過程中發生的空氣中各原子電離、分子離解及激發的各種分子振動均需要消耗能量,因此,在高溫條件下,空氣的狀態方程與理想氣體狀態方程有顯著不同 [8,9]。俄羅斯學者Kuznetzov在試驗基礎上,采用統計物理學和沖擊波氣體動力學的理論方法得到了空氣熱力學函數表和沖擊絕熱線的計算表格,Brode [10]對表格數據進行擬合得到了真實空氣狀態方程的擬合公式 :


式(10)~式(16)中,e,θ,π,η分別為無量綱內能、無量綱溫度、無量綱壓力和無量綱密度。
2. TNT爆轟產物狀態方程
Brode [10]從能量角度出發,對TNT裝藥爆轟產物的狀態方程進行了詳細研究,并給出了與真實空氣狀態方程形式相同的描述方程:

式(19)中,α1= 1.76,b1= 52.4892,c1= 7.3。
利用量綱分析方法中的Π定理,可得到由主定參量初始密度ρ0、大氣聲速C0、爆炸所釋放的能量E、球形坐標r和時間t構成的面源、線源和點源TNT裝藥強爆炸的獨立自模擬解變量分別為ξ =若沖擊波陣面的位置為rs,則ξ取值為ξ0,則可得到面源、線源和點源TNT對稱裝藥強爆炸自模擬解的控制方程組分別如下:
面源裝藥:

線源裝藥:

點源裝藥:

邊界條件:在波陣面上,f(1) = g(1) = h(1) = 1;在中心處,f(0) = 0。其中,f,g和h分別為無量綱速度、無量綱密度和無量綱壓力。
式(21)~式(23)為封閉方程組,可利用其進行求解。
通過計算可得到面源、線源和點源裝藥強爆炸自模擬數值解,結果如圖1~圖3所示。圖中橫、縱坐標的數值區間(0~1)均為無量綱形式,其中橫坐標的相對位置,即r/rs,為無量綱形式的距爆炸中心的距離。
從圖1~圖3可以看出,面源、線源和點源裝藥強爆炸解參數隨波陣面相對位置分布的特征基本相同,均呈現類似的變化趨勢。質點速度隨著r/rs下降基本上呈線性下降趨勢;壓力在r/rs從1降低到0的過程中先下降比較快,當r/rs降到0.5之后直到爆炸中心,壓力基本穩定在0.36~0.37;密度隨r/rs降低很快接近至0。
試驗表明,自模擬解 [7]可以很好地描述爆炸的初始階段。因此,利用真實空氣條件下面源、線源和點源裝藥強爆炸的自模擬解,作為TNT裝藥爆炸近區流場特性數值計算的初始條件是精確合理的。
選取爆轟波傳至裝藥表面的時刻為初始計算時刻,爆轟波的波陣面將整個流場分為兩個部分(見圖4)。
根據守恒定律,在激波陣面上必須滿足如下關系:
質量守恒:

動量守恒:

圖1 面源裝藥強爆炸解參數分布圖

圖2 線源裝藥強爆炸解參數分布圖

圖3 點源裝藥強爆炸解參數分布圖

能量守恒:

式(24)~式(26)中,P0、u0、ρ0、E0分別表示激波前未擾動空氣的壓力、質點速度、密度和內能;用P、u、ρ、E表示激波陣面后擾動空氣的壓力、質點速度、密度和內能;D為沖擊波速度。

圖4 初始計算時刻流場分布示意圖
對式(24)~式(26)無量綱化后,可以解出單值對應的激波關系。對偏微分方程組(9)進行差分近似,為了保持接觸間斷面處和激波跳躍間斷面處計算的貫通,采用人工粘性法,在出現壓力項處附加人工粘性項ζ,并將其無量綱化:ρ = ζ/η,于是將方程組(9)(N = 2時為球形裝藥爆炸)變為無量綱差分格式:

式(27)中,l為差分網格整格點,n為時間整步長,e為空間半格l–1/2和時間整格點n + 1處。
其中粘性項為:

關于能量方程的差分,由于參數比較復雜,因此,采用新的方程式,從無量綱狀態方程出發,將e進行代換可得到新的形式:

式(28)中,ξ = lnη。
從而得到能量的差分格式:

因上式是一個隱式的差分格式,因此采用二次迭代進行求解。建立激波后裝配方法:首先進行激波探查,用三次樣條插值求ξ的最大值點,將此位置定為激波陣面位置λs,求得激波陣面上的壓力,然后調用波陣面參數的求解子程序求得波陣面上的其他參數。數值計算的時間步長由Courant條件和激波區域的差分方程穩定性條件共同確定。
1. 空氣和TNT裝藥參數
靜止空氣的參數見表1。

表1 靜止空氣參數
TNT裝藥參數見表2。

表2 TNT裝藥參數
2. 結果
采用本文建立的計算方法可得到下述結果:
(1)超壓
圖5為無量綱峰值超壓ΔPs、峰值動壓qs隨無量綱激波陣面距離位置的變化曲線,從圖中可看出,本文的計算結果與Baker [11]的計算結果吻合度較好。
(2)接觸間斷面和主激波陣面的時–空曲線
圖6為主激波(MS)、二次激波(2S)以及接觸間斷面(CS)無量綱位置λ隨無量綱時間τ的變化規律。從圖中可看到Friedman [12]發現的二次激波現象,而出現的多個后續激波是由于波的反射產生的。從圖6可得到沖擊波的影響范圍遠遠大于爆轟產物接觸間斷面的擴散范圍。
(3)峰值超壓、峰值動壓隨激波陣面距離變化曲線
圖7和圖8分別為不同裝藥量球形TNT裝藥在真實空氣中爆炸產生沖擊波的主激波峰值超壓及峰值動壓隨激波陣面距離的變化曲線圖(采用的是指數坐標)。由圖7、圖8可知,對于同一裝藥量,球形TNT裝藥隨著距離的不斷增大,超壓峰值和動壓峰值不斷下降。不同裝藥量的炸藥爆炸超壓、動壓隨距離變化符合相似律關系。隨著裝藥量的增加,相同位置處的超壓峰值和動壓峰值相應增加。

圖5 峰值超壓、峰值動壓隨激波半徑的變化曲線(實心點為本文結果,空心點為Baker的計算結果)

圖6 沖擊波與接觸間斷面的時–空曲線
(4)沖擊波到達時間隨激波陣面距離變化的曲線
圖9為不同裝藥量球形TNT裝藥在真實空氣中爆炸產生的沖擊波到達時間隨激波陣面距離的變化規律。可見隨著裝藥量的增加,球形TNT裝藥爆炸沖擊波到達相同位置處的時間不斷降低。不同裝藥量的沖擊波到達時間變化規律大致相同,可見同樣滿足相似律關系。
(5)沖擊波超壓沖量隨激波陣面距離變化的曲線

圖8 不同裝藥量球形TNT裝藥爆炸峰值動壓隨激波陣面距離變化的曲線

圖9 不同裝藥量球形TNT裝藥爆炸沖擊波到達時間隨激波陣面距離變化的曲線
圖10為不同裝藥量球形TNT裝藥爆炸沖擊波超壓沖量隨激波陣面距離的變化曲線。從圖中可知,對于同一裝藥量,球形TNT裝藥隨著激波陣面距離的不斷增大,沖量先逐漸降低,當到達一定距離后,在此區域會出現一個平緩甚至稍微增加的階段,之后再開始衰減,這是由于二次激波追趕上了主激波使得正壓作用沖量增加造成的;另一方面,隨著裝藥量的增加,相同位置處的沖量也逐漸增加。
(6)沖擊波的正壓作用時間隨激波陣面距離變化的曲線
圖11是不同裝藥量的球形TNT裝藥空中爆炸沖擊波的正壓作用時間隨激波陣面距離變化的曲線。從圖中可知,對于同一裝藥量,爆炸沖擊波隨激波陣面距離的變化規律與圖10類似,正壓作用時間也存在一個轉折點,它先隨距離的增加而降低,但在某個位置處又開始隨距離的增大而增大。這主要是因為存在二次激波的原因,隨著沖擊波不斷向外傳播,主激波的峰值超壓不斷降低,且主激波越晚到達的點的正壓作用時間就會越少,但是由于二次激波的不斷追趕,當到達某個位置后,主激波的正壓還未衰減到負壓,而二次激波就追趕上來,從而使得正壓作用持續時間增加。對于1~50 kg的TNT裝藥,轉折點的位置彼此相差不大,基本在4 m左右;而對于100~1 000 kg的大裝藥量藥球,轉折點離爆心的距離隨裝藥量的增大而增大。另外,對于裝藥量為1~50 kg的藥球,在固定的位置點,轉折點前后均有正壓作用時間隨裝藥量的增大而增大的趨勢;但對于裝藥量為100~1 000 kg的藥球,在轉折點前,正壓作用時間同樣隨裝藥量的增大而增大,但是在轉折點后則相反,正壓作用時間隨裝藥量的增加反而降低(見圖11中的箭頭顯示)。

圖10 不同裝藥量球形TNT裝藥爆炸沖擊波超壓沖量隨激波陣面距離變化的曲線
由于空氣中爆炸的近區空氣與理想氣體之間存在較大差異,因此,本文采用真實空氣狀態方程來描述爆炸近區流場的特性。基于本文所建立的理論,結合人工粘性技術和激波裝配技術,對TNT裝藥真實空氣中爆炸流場參數進行了計算,結果與已有文獻研究結果相吻合,并得到以下結論:

圖11 不同裝藥量球形TNT裝藥爆炸沖擊波正壓作用時間隨激波陣面距離變化的曲線
(1)TNT裝藥爆炸之后會產生一個初始超壓為60 MPa的強沖擊波(圖5中無量綱波陣面位置λs= 1對應的峰值超壓),強度隨沖擊波的不斷傳播而逐漸衰減,并且在主激波之后存在一個不斷向外傳播的接觸間斷面。
(2)爆炸近區最典型的物理特征是存在一個爆轟產物區,在半徑小于3倍裝藥半徑的區域內,壓力強度遠遠大于主激波的壓力值。
(3)當爆轟波到達裝藥與空氣的接觸面時,在形成向外傳播的空氣沖擊波的同時,還會形成一個向內傳播的激波,向內傳播的激波在爆心處發生碰撞后會得到加強,然后繼續向外傳播,在到達接觸面后又形成一個新的向內傳播和向外傳播的激波,從而形成定點超壓波形負壓區的二次、三次激波。
由此可見,本文提出的計算分析方法可以很好地模擬和描述TNT裝藥爆炸近區在真實空氣中的傳播規律。本文的研究成果可以對爆炸力學和防護工程研究起到重要的參考作用。
[1] 盧芳云. 一維不定常流體動力學教程 [M]. 北京∶ 科學出版社,2006.
Lu F Y. Tutorial of one-dimensional unsteady fl ow dynamics [M].Beijing∶ China Science Publishing &Media Ltd., 2006.
[2] Taylor G. The formation of a blast wave by a very intense explosion. I. Theoretical discussion [J]. Proceeding of the Royal Society of London. 1950 (201)∶ 159–174.
[3] Von Neumann J. The point source solution [R]. Washington, D.C.∶National Defense Research Committee, 1947.
[4] Bach G G, Lee J H S. An analytical solution for blast waves [J].AIAA Journal, 1970 (8)∶ 271–275.
[5] Sedov L I. Similarity and dimensional methods in mechanics [M].New York∶ CRC Press, 1965.
[6] 湯文輝. 沖擊波物理 [M]. 北京∶ 科學出版社, 2011.
Tang W H. Shock wave physics [M]. Beijing∶ China Science Publishing &Media Ltd., 2011.
[7] Л.П.奧爾連科. 爆炸物理學(原書第三版) [M]. 孫承緯譯. 北京∶ 科學出版社, 2011.
Л.П.ОРЛЕНКО. Explosion physics [M]. Translated by Sun C W. Beijing∶ China Science Publishing &Media Ltd., 2011.
[8] 喬登江. 核爆炸物理概論 [M]. 北京∶ 原子能出版社, 1988.Qiao D J. Introduction to nuclear explosion physics [M]. Beijing∶Atomic energy Press, 1988.
[9] 周豐峻, 陳葉青, 任輝啟. 爆炸近區空氣沖擊波規則反射和非規則反射 [C].應用力學進展論文集. 北京∶ 中國力學學會, 1999∶1–12.
Zhou F J, Chen Y Q, Ren H Q. The regular reflection and nonregular re fl ection of shockwave in the close-in fi eld of TNT charge explosion [C]. Proceedings of applied mechanics progress. Beijing∶The Chinese Society of Theoretical and Applied Mechanics, 1999∶1–12.
[10] Brode H L. Blast wave from a spherical charge [J]. Physics of Fluids, 1959, 2 (2)∶ 217–229.
[11] Baker W E. Explosives in air [M]. Austin∶ University of Texas Press, 1973.
[12] Friedman M P. A simplified analysis of spherical or cylindrical blast waves [J]. Journal of Fluid Mechanics, 1961, 11 (1)∶ 1–15.