PBX9502炸藥爆轟產(chǎn)物的狀態(tài)方程*

趙艷紅，劉海風(fēng)，張廣財

（北京應(yīng)用物理與計算數(shù)學(xué)研究所，北京 100094）

PBX9502炸藥爆轟產(chǎn)物的狀態(tài)方程*

趙艷紅，劉海風(fēng)，張廣財

（北京應(yīng)用物理與計算數(shù)學(xué)研究所，北京 100094）

采用van der Waals等效單組分流體模型和Ross硬球微擾理論軟球修正模型，計算了爆轟氣相產(chǎn)物的狀態(tài)方程；用石墨相、金剛石相、類石墨液相和類金剛石液相4種相態(tài)描述凝聚成分，由Gibbs自由能最小確定了不同狀態(tài)下的凝聚產(chǎn)物相態(tài)。對爆轟產(chǎn)物混合系統(tǒng)采用自由能最小原理，通過化學(xué)平衡方程組求解了炸藥爆轟產(chǎn)物系統(tǒng)的平衡組分。使用該理論計算了高含碳炸藥PBX9502Chapman-Jouguet（CJ）點的爆轟參數(shù)，計算值與實驗值符合很好；同時計算了3條等溫線，并與Sesame庫比較，發(fā)現(xiàn)溫度超過1 000℃時，計算值與Sesame庫的計算結(jié)果比較接近。在計算的5 802K等溫線上發(fā)現(xiàn)了一個拐點，分析發(fā)現(xiàn)是由于在此處游離態(tài)的碳發(fā)生了相變。

爆炸力學(xué)；狀態(tài)方程；化學(xué)平衡方程組；爆轟產(chǎn)物；PBX9502炸藥

1 引 言

爆轟產(chǎn)物的狀態(tài)方程對描述炸藥的爆轟性質(zhì)至關(guān)重要，它是爆轟過程數(shù)值模擬的核心參數(shù)，也是炸藥作功能力的主要表征。但是，在實際的工程應(yīng)用中，爆轟產(chǎn)物的狀態(tài)方程卻一直存在問題，主要是由于描述爆轟產(chǎn)物的物理過程太復(fù)雜了。當(dāng)爆轟產(chǎn)物剛形成時，其狀態(tài)在CJ點附近，溫度高達數(shù)千開，壓力高達數(shù)十吉帕，密度為2～3g／cm3，此時分子之間的相互作用類似于固體或液體的性質(zhì)，隨著產(chǎn)物的膨脹，壓力降到環(huán)境氣壓，分子之間的作用又呈現(xiàn)氣體性質(zhì)。

對爆轟產(chǎn)物的狀態(tài)方程，采用統(tǒng)計物理方法自編類CHEQ程序，對低含碳炸藥PETN進行計算，結(jié)果很好［1］。本文中用類 CHEQ 程序計算高含碳炸藥 PBX9502（C8.27253H6.08519N8O6F0.3662Cl0.12267）CJ點的爆速，與實驗值符合較好。為與Sesame庫的數(shù)據(jù)相比，計算了3條等溫線，發(fā)現(xiàn)在我們關(guān)心的區(qū)域，計算值與Sesame庫的計算結(jié)果比較接近。

2 理論模型簡介

考慮一個有t種可能組分的系統(tǒng)，包括s種氣體，p種固體物質(zhì)，它們由c種不同的化學(xué)元素組成。對一化學(xué)反應(yīng)，反應(yīng)物中各元素的原子個數(shù)與生成物中各元素的原子個數(shù)相等，即滿足化學(xué)反應(yīng)質(zhì)量守恒。由分子式矢量表示得

式中：aij表示第j種元素在第i種組分分子式中的原子個數(shù)，ni為第i種組分的物質(zhì)的量，qj是系統(tǒng)中第j種元素的總的物質(zhì)的量。

化學(xué)平衡要求反應(yīng)物的化學(xué)勢必須等于生成物的化學(xué)勢?；瘜W(xué)平衡時自由能最小，據(jù)此可以導(dǎo)出化學(xué)平衡方程組［2］。

炸藥爆轟產(chǎn)物混合系統(tǒng)的Helmholtz自由能可寫為［3－4］

式中：Af、Ac分別為流體和固體的Helmholtz自由能，nc為固體的物質(zhì)的量，｛ni｝表示流體中各組元的物質(zhì)的量。流體自由能Af又可表示為不考慮分子間相互作用貢獻的相應(yīng)理想氣體的自由能Aid和體系分子間相互作用貢獻的超額自由能 Aint之和［5－6］，即

對于理想氣體部分的Helmholtz自由能Aid，采用基于Rossiwi（單原子）、Pennigton和Kobe（雙原子）、諧振子－剛性轉(zhuǎn)子近似（多原子）等理論計算給出的擬合函數(shù)結(jié)果［1］。

假定流體中各組元分子間相互作用勢為

則根據(jù)各組元濃度（xi＝ni／）和分子間相互作用參數(shù)（εij，αij，），采用van der Waals等效單組分流體模型（vdw1f），可將混合物等效為一元流體，等效分子勢φ（r）中的參數(shù)（ε，α，r＊）可以通過下列推廣的混合法則計算［7］，即

而（εij，αij）采用推廣的Lorent-Berthelor組合規(guī)則求出［8］，即

式中：lij、mij、kij為修正系數(shù)，除kH2OCO2＝0.965，kH2ON2＝1.03外，其余l(xiāng)ij＝mij＝kij＝1［9］。用求出的炸藥的爆轟產(chǎn)物組分ni，結(jié)合式（5）、（6），可以求出vdwlf模型中等效分子勢φ（r）中的參數(shù)（ε，α，r＊）。

采用Ross硬球微擾理論軟球修正模型，計算炸藥爆轟氣相產(chǎn)物體系等效單組分流體分子間相互作用貢獻的超額自由能Aint。

炸藥爆轟氣相產(chǎn)物分子間相互作用貢獻為

式中：AHS、gPY（r／d，η）和F12（η）分別為硬球超額自由能、Percus-Yevick硬球徑向分布函數(shù)和軟球修正項，η＝πρd3／6，ρ＝N／V，d是硬球直徑，d的選取應(yīng)使式（7）右邊最小，并將此最小值作為體系當(dāng)前狀態(tài)下的Helmholtz自由能Aint。

碳在高壓下的狀態(tài)方程在沖擊波物理中非常重要，由于爆轟是一個瞬態(tài)的高溫高壓過程，很難進行直接觀察和測量，目前還沒有一種方法能夠給出爆轟CJ點上的爆轟產(chǎn)物的真實組成。

在早期的研究中幾乎都把爆轟產(chǎn)物的碳作為石墨處理，后來把爆轟產(chǎn)物中的碳當(dāng)作石墨或金剛石的逐漸增多。最近，L.E.Fried等［9］對碳的4種相態(tài)（石墨、金剛石、類石墨液碳、類金剛石液碳）作了詳細的研究，研究結(jié)果與 M.Togaya［10］和F.P.Bundy等［11］的實驗結(jié)果符合得很好。文獻［9］給出了碳的4種相的Gibbs自由能具體計算公式（適用范圍：0≤p≤600GPa，300K≤T≤15kK）

式中：p0＝100kPa。由式（8），計算出給定（p，T）下碳的石墨相、金剛石相、類石墨液相和類金剛石液相4種相態(tài)的Gibbs自由能，由Gibbs自由能最小原理，從4種相態(tài)中選出炸藥爆轟產(chǎn)物中游離態(tài)的碳最可能存在的相態(tài)，根據(jù)熱力學(xué)性質(zhì)計算其他的熱力學(xué)量。相比于Murnaghan狀態(tài)方程，L.E.Fried在碳的4相狀態(tài)方程中將溫度引了進來，因此更符合實際情況。

根據(jù)熱力學(xué)知識，體系的壓強p、能量E、熵S和各產(chǎn)物的化學(xué)勢μi可由總的自由能A求出

具體的表達式可以參考文獻［1］。

3 計算結(jié)果與分析

采用類CHEQ程序計算了初始密度ρ0＝1.894g／cm3的高含碳炸藥PBX9502在CJ點的爆轟參數(shù)，計算結(jié)果見表1。從表1可以看到，本文中和BKW計算的爆速都與實驗值符合很好（小于國際最好標(biāo)準(zhǔn)5%），但BKW狀態(tài)方程中有4個常數(shù)，這些常數(shù)通過擬合產(chǎn)物組分的Hugoniot線和炸藥爆轟實驗的有關(guān)數(shù)據(jù)而得到。而本文中僅用到各產(chǎn)物的分子間勢，采用統(tǒng)計原理計算的爆轟參數(shù)，其中沒有與炸藥爆轟實驗相關(guān)的可調(diào)參數(shù)。

表1 PBX9502炸藥爆轟CJ點的爆轟參數(shù)Table 1 Detonation parameters at the CJ point for PBX9502explosive

近代，由于工程技術(shù)的需要，科學(xué)家們建立了許多狀態(tài)方程數(shù)據(jù)庫，著名的有美國的Sesame庫和俄羅斯的DBMS庫。為進一步驗證爆轟產(chǎn)物物態(tài)方程的正確性，在圖1中比較了1 160、2 320和5 802K等3條等溫線上本文中和Sesame庫計算的壓力隨密度的變化。從圖1中可以看到，本文計算結(jié)果與Sesame庫的計算結(jié)果比較接近。

在圖1中，本文中計算的5 802K等溫線上發(fā)現(xiàn)，壓強有一突變。由于氣相產(chǎn)物狀態(tài)方程是連續(xù)變化的，不會引起突變，所以可能是固態(tài)產(chǎn)物碳的狀態(tài)方程發(fā)生了改變。在圖2中，將游離態(tài)的碳設(shè)為單一的石墨相（Graphite）、金剛石相（Diamond）、類石墨液碳相（Liquid 1）和類金剛石液碳相（Liquid 2），并與碳的4相態(tài)狀態(tài)方程（CHEQ）相比，可以看出，在密度小于2.9g／cm3時，4相態(tài)狀態(tài)方程與類金剛石液碳相計算出的壓強相等；密度大于2.9g／cm3時，與金剛石相計算出的壓強相等。因此，在密度為2.9g／cm3附近，游離態(tài)的碳由類金剛石液碳轉(zhuǎn)化為固態(tài)金剛石，從而引起熱力學(xué)量的突變。

圖1 等溫線上壓力隨密度的變化Fig.1Isothermal pressure-density curves calculated by the two different methods

圖2 碳取不同的相態(tài)時等溫線上壓力隨密度的變化Fig.2 Isothermal pressure-density curves in the case of caebon at different states

碳的相變，不僅引起壓強的突變，還引起各爆轟產(chǎn)物組分的變化。在圖3中，列舉了5 802K等溫線上主要產(chǎn)物隨密度的變化。從圖3可以看到，在5 802K等溫線上，隨著產(chǎn)物密度的增加，炸藥中的N元素主要以N2的形式存在，H元素主要以H2O的形式存在，H2O含量的增加，致使炸藥中與C元素結(jié)合的O元素減少，CO逐漸減少，游離態(tài)的碳析出。這些曲線隨密度的變化應(yīng)是光滑的，但在密度為2.9g／cm3附近，都有一突變，這是由于碳的相變引起的。

圖3 等溫線上各產(chǎn)物物質(zhì)的量隨密度的變化Fig.3 Variation of the amounts of substance with density along isothermal curves for detonation products

4 結(jié) 論

（1）本文中計算的CJ點的爆轟參數(shù)與實驗值符合很好。

（2）本文中計算的等溫線與Sesame庫的計算結(jié)果比較接近。

（3）本文中計算的5 802K等溫線上發(fā)現(xiàn)了一個拐點，分析發(fā)現(xiàn)是因為在此處游離態(tài)的碳發(fā)生了相變，由類金剛石液碳轉(zhuǎn)化為固態(tài)金剛石，從而引起熱力學(xué)量及產(chǎn)物組分的突變。

［1］趙艷紅，劉海風(fēng)，張弓木.基于統(tǒng)計物理的爆轟產(chǎn)物物態(tài)方程研究［J］.物理學(xué)報，2007，56（8）：4791-4797.

ZHAO Yan-hong，LIU Hai-feng，ZHANG Gong-mu.Equation of state of detonation products based on statistical mechanical theory［J］.Acta Physica Sinica，2007，56（8）：4791-4797.

［2］劉海風(fēng)，陳棟泉，張世澤.爆轟產(chǎn)物物態(tài)方程及CHBr3相變的理論研究［J］.高壓物理學(xué)報，1996，10（4）：284-290.

LIU Hai-feng，CHEN Dong－quan，ZHANG Shi-ze.Equation of state of detonation products and the possible phase transition for CHBr3［J］.Chinese Journal of High Pressure Physics，1996，10（4）：284-290.

［3］楊向東，謝文，武保劍，等.液氮的沖擊壓縮理論計算［J］.高壓物理學(xué)報，1998，12（1）：1-7.

YANG Xiang-dong，XIE Wen，WU Bao-jian，et al.Theoretical calculation for the hugoniot curves of liquid nitrogen［J］.Chinese Journal of High Pressure Physics，1998，12（1）：1-7.

［4］劉福生，陳先猛，陳攀森，等.液態(tài)CO2高溫高密度狀態(tài)方程研究［J］.高壓物理學(xué)報，1998，12（1）：28-33.

LIU Fu-sheng，CHEN Xian-meng，CHEN Pan-sen，et al.Equation of sate of liquid CO2at high temperatures and high densities［J］.Chinese Journal of High Pressure Physics，1998，12（1）：28-33.

［5］楊向東，胡棟，經(jīng)福謙.炸藥爆轟產(chǎn)物液氮、液氦和水狀態(tài)方程研究［J］.高壓物理學(xué)報，1999，13（2）：93-102.

YANG Xiang-dong，HU Dong，JING Fu-qian.Studies of EOS for detonation products：Liquid nitrogen，liquid helium and water［J］.Chinese Journal of High Pressure Physics，1999，13（2）：93-102.

［6］李德華，楊繽維，程新路，等.液 H2O沖擊壓縮特性的理論計算［J］.四川師范大學(xué)學(xué)報（自然科學(xué)版），2005，42（1）：108-111.

LI De-hua，YANG Bin-wei，CHENG Xin-lu，et al.Theoretical calculated of shock-compression properties for liquid water［J］.Journal of Sichuan University（Natural Science Edition），2005，42（1）：108-111.

［7］Ree F H.Simple mixing rule for mixtures with exp-6interactions［J］.Journal of Chemical Physics，1983，78（1）：409-415.

［8］Ree F H.A statistical mechanical theory of chemically reacting multiphase mixtures：Application to the detonation properties of PETN［J］.Journal of Chemical Physics，1984，81（3）：1251-1263.

［9］Fried L E，Howard W M.Explicit Gibbs free energy equation of state applied to the carbon phase diagram［J］.Physical Review B，2000，61（13）：8734-8743.

［10］Togaya M.Pressure dependences of the melting temperature of graphite and the electrical resistivity of liquid carbon［J］.Physical Review Letters，1997，79（13）：2474-2477.

［11］Bundy F P，Bovenkerk H P，Strong H M.Diamond-graphite equilibrium line from growth and graphitization of diamond［J］.Journal of Chemical Physics，1961，35（2）：383-391.

［12］Mader C L.Numerical modeling of detonation［M］.London：University of California Press，1979.

Equation of state of detonation products for PBX9502explosive＊

ZHAO Yan-hong，LIU Hai-feng，ZHANG Guang-cai
（Beijing Institute of Applied Physics and Computational Mathematics，Beijing100094，China）

The equation of state of gas detonation products was described by Ross’s modification of hard-sphere variation theory and the improved one-fluid van der Waals mixture model.The Gibbs free energy of dissociated carbon was calculated for the most probable state，which was determined by distinguishing the following four states of carbon：graphite，diamond，graphitelike and diamondlike.The equilibrium compositions of detonation products are calculated by solving chemical equilibrium equations based on minimizing free energy.The detonation properties at the CJ point of PBX9502explosive were calculated with this theory.The results show satisfactory agreement with the experimental data.Comparison of the isothermal pressure-density curves displays that the results calculated with the present theory are in good agreement with those based on Sesame database at high temperature.A slope is found at 5 802Kisotherm.The reason is that the carbon phase is changed at this point.

mechanics of explosion；equation of state；chemical equilibrium equations；detonation products；PBX9502explosive

20August 2009；Revised 6November 2009

ZHAO Yan-hong，zhao＿yanhong＠iapcm.ac.cn

（責(zé)任編輯 曾月蓉）

O381 國標(biāo)學(xué)科代碼：130·35

A

1001-1455（2010）06-0647-05

2009-08-20；

2009-11-06

中國工程物理研究院科學(xué)技術(shù)發(fā)展基金項目（2008B0201019）

趙艷紅（1977— ），女，碩士，助理研究員。