0、引言
光輻射是強爆炸的重要殺傷作用之一,通過對光輻射的定性分析,對于了解強爆炸的性質、評估其毀傷作用具有重要意義.早期有關光輻射的研究,建立在試驗測量的基礎上,結合相關的規律性分析,總結給出了光輻射功率走時的經驗關系;在數值模擬方面,采用輻射流體力學方法,對光輻射功率進行了計算分析.在計算過程中,由于光輻射能流的強不穩定性\\(非線性、強烈依賴于物理參數等\\),因而通過一定的方法進行處理:一種是利用輻射熱傳導模型,輻射能流采用相應的擴散近似.
由于假設條件的苛刻,導致光輻射功率計算不夠細致;另一種方法,對整個區域采用輻射流體力學方法直接求解,為避免輻射能流的劇烈振蕩,對能流進行了平滑,這在一定程度上降低了輻射能流的計算精度,對于光輻射功率求解造成了困難.本文基于輻射流體力學計算方法,對光輻射輸運過程采用直接求解.
為提高輻射能流的計算精度,并最大程度的提高計算效率,利用算子分裂方法將方程分裂為對流項和源項,二者獨立求解,有效解決了光輻射計算中的輻射能流不穩定性問題,同時計算效率在原有基礎上也得到了明顯提升.給出的光輻射功率走時曲線,與已有結果和經驗公式一致.
1、計算模型
強爆炸光輻射發展過程,是光輻射與流體耦合作用的過程,因而采用輻射流體力學方法,構建輻流耦合的物理模型.為簡化方程求解,通常采用灰體近似方法,輸運過程以矩方程形式給出.其中,輻射壓強由最大熵Eddington因子方法確定:
流體力學方程組:
其中ρ為空氣密度,u為空氣速度,p為空氣壓強,e為空氣總能量,κ′0為靜止坐標系中考慮了受激發射的空氣吸收系數;Er為輻射能密度,Fr為輻射能流,c為光速;φ=aT4為溫度T下的黑體輻射能密度\\(a=5.67×10-16J/\\(m2·K4\\),為斯蒂芬-波爾茲曼常數\\).
最大熵Eddington因子下的輻射壓強Pr=χEr:
空氣的狀態方程采用實際空氣狀態方程,即:
而eI、ρ分別為空氣內能及密度.
通過以上方程組,在充分保證求解精度的基礎上,能夠完整、全面的反映強爆炸光輻射物理過程.
對于光輻射功率計算,取空氣溫度2000K為火球陣面,取計算得到的陣面處輻射能流,根據式\\(R為火球半徑\\):
P=4πR2Fr\\(5\\)計算得到光輻射功率.
2、求解方法
采用算子分裂法的方案,將輻流方程組分為對流項和源項單獨處理,原輻射流體力學方程組可化為如下形式:
對流項\\(式\\(6\\)第一項\\)通過有限體積法,構造五階WENO格式,數值通量的計算采用局部Lax-Friedrichs方法,從而實現高精度求解.剛性源項\\(式\\(6\\)第二項\\)的處理,最終是化成了一個常微分方程組的求解問題.求解初始條件假定強爆炸總能量集中于等壓火球內,輻射能與空氣內能總合等于爆炸總能量,計算邊界條件采用對稱邊界.
3、計算結果與分析
取當量為1kt、高度在海平面\\(h=0km\\),求解相應的強爆炸火球光輻射輸運過程,空氣初始狀態參量如下:
計算得到的強爆炸火球陣面走時\\(如圖1\\)符合火球的物理過程,與Brode計算結果符合良好,證明了方法的可靠性.
為研究強爆炸光輻射功率走時,需要先確定輻射能流隨時空變化關系\\(式\\(5\\)\\).
圖2、3分別給出了強爆炸光輻射不同時刻\\(4ms前后4ms后\\)輻射能流的時空分布.
從圖中可以看出,輻射能流沒有出現不穩定情況,同時時間步長較原有直接求解大為提高.光輻射能流在不同時刻的空間分布,有兩個階段:早期\\(4ms前\\),輻射能流隨時間和空間位置向外單調減小,說明光輻射功率隨爆炸能量向外發展逐漸減弱;在4ms以后,光輻射能流在空間分布上出現了逐漸增大過程,成為“雙脈沖”.
光輻射能流的第二次脈沖,說明在后期階段,強爆炸光輻射逐漸增大,與火球的“復燃”過程相對應.
輻射能流反映的是光輻射輸運能量的大小,與空氣吸收系數有著重要關系.
下圖給出了4ms前后輻射能流與空氣吸收系數、溫度分布對比關系.
在t=2.88ms時刻,火球的輻射能流在空間分布為單脈沖,對應火球的溫度均在10000K以上,而空氣吸收系數也為單脈沖:輻射能流較大的位置,吸收很弱,而當吸收系數較大時,輻射能流開始減小;在t=8.27ms時刻,輻射能流在空間分布為雙脈沖,在第二脈沖附近,火球溫度已經降到幾千K,此階段,空氣吸收系數在空間分布也出現了兩個階段:
第一階段的分布與2.88ms基本一致,此后吸收系數在3000K附近,卻出現了極小值,這也使得輻射能流出現第二次增大.
空氣吸收系數再次極小的出現,與空氣溫度、密度有關,也是導致火球光輻射第二極大的根本原因.
分析了輻射能流隨時空分布及其與空氣吸收系數的關系之后,利用式\\(5\\)給出光輻射功率走時關系\\(火球陣面為2000K,對應位置為火球半徑\\):在計算的時間范圍內,光輻射功率走時從開始逐漸減小,達到極小值后又開始增大,與實際火球發展過程一致.
其極小值對應于t=4.0ms,根據已有文獻給出的結果,在1kt強爆條件下,光輻射功率極小值約為3.88ms,二者符合較好.
4、結論
采用輻射流體力學方法,直接求解了灰體近似下的強爆炸火球光輻射輸運過程.利用算子分裂方法,實現對流項和強源項獨立求解,保證了求解過程中輻射能流的穩定性,同時也極大提高了方程求解效率,在此基礎上給出了1kt強爆條件下光輻射功率走時曲線.計算結果符合強爆火球發展過程,與已有結果一致.
參考文獻:
[1] 屠琴芬.強爆炸火球的輻射流體力學計算中的幾個問題[R].西安:西北核技術研究所,1986.
[2] 陳健華,王心正,謝龍生等.均勻大氣中的強爆炸一維輻射流體力學數值解[J].爆炸與沖擊,1981,1\\(2\\):37-49.
[3]Eugene M D Symbalisty,John Zinn,Rodnq W Whi-taker.Radflo Physics and Algorithms.LA-12988-MS.1995.
[4]Crowley K Barbara,David H Glenn,Marks E Robert.An analysis of marvel-a nuclear shock-tube experi-ment[J].Journal of Geophysical Research,1971,76\\(14\\):3356-3374.
[5]Marrs R E,Moss W C,Whitlock B.Thermal Radia-tion from Nuclear Detonations in Urban Environ-ments[R].UCRL-TR-231593,June 7,2007.
[6] 田宙,喬登江,郭永輝.不同高度強爆炸早期火球數值研究[J].兵工學報,2009,30\\(8\\):1078-1083.
[7] 高銀軍,田宙,劉峰等.基于多群方法的強爆炸早期火球數值模擬[J].應用物理,2011,2\\(2\\):163-167.
[8] 田宙,喬登江,郭永輝.不同當量強爆炸早期火球現象的數值模擬[J].爆炸與沖擊,2009,29\\(4\\):408-412.
[9]Minerbo Gerald N.Maximum entropy eddington fac-tors[J].J.Quant.Spectrosc.Radiat.Transfer,1977,20\\(6\\):33-42.
[10]閆凱.二維輻射流體動力學方程組的數值求解.西北核技術研究所,碩士學位論文,2011.
[11]閆凱.二維輻射流體動力學方程組的數值求解.第六屆全國青年計算物理學術會議,2011.
[12]喬登江.核爆炸物理概論[M].北京:國防工業出版社,2003.
[13]屠琴芬.高空強爆炸火球的一維輻射流體力學數值模擬[R].西安:西北核技術研究所,1987.
[14]屠琴芬.地下強爆炸早期火球的輻射流體力學總結[R].西安:西北核技術研究所,