PBX9502炸药爆轰产物的状态方程*

赵艳红，刘海风，张广财

（北京应用物理与计算数学研究所，北京 100094）

PBX9502炸药爆轰产物的状态方程*

赵艳红，刘海风，张广财

（北京应用物理与计算数学研究所，北京 100094）

采用van der Waals等效单组分流体模型和Ross硬球微扰理论软球修正模型，计算了爆轰气相产物的状态方程；用石墨相、金刚石相、类石墨液相和类金刚石液相4种相态描述凝聚成分，由Gibbs自由能最小确定了不同状态下的凝聚产物相态。对爆轰产物混合系统采用自由能最小原理，通过化学平衡方程组求解了炸药爆轰产物系统的平衡组分。使用该理论计算了高含碳炸药PBX9502Chapman-Jouguet（CJ）点的爆轰参数，计算值与实验值符合很好；同时计算了3条等温线，并与Sesame库比较，发现温度超过1 000℃时，计算值与Sesame库的计算结果比较接近。在计算的5 802K等温线上发现了一个拐点，分析发现是由于在此处游离态的碳发生了相变。

爆炸力学；状态方程；化学平衡方程组；爆轰产物；PBX9502炸药

1 引 言

爆轰产物的状态方程对描述炸药的爆轰性质至关重要，它是爆轰过程数值模拟的核心参数，也是炸药作功能力的主要表征。但是，在实际的工程应用中，爆轰产物的状态方程却一直存在问题，主要是由于描述爆轰产物的物理过程太复杂了。当爆轰产物刚形成时，其状态在CJ点附近，温度高达数千开，压力高达数十吉帕，密度为2～3g／cm3，此时分子之间的相互作用类似于固体或液体的性质，随着产物的膨胀，压力降到环境气压，分子之间的作用又呈现气体性质。

对爆轰产物的状态方程，采用统计物理方法自编类CHEQ程序，对低含碳炸药PETN进行计算，结果很好［1］。本文中用类 CHEQ 程序计算高含碳炸药 PBX9502（C8.27253H6.08519N8O6F0.3662Cl0.12267）CJ点的爆速，与实验值符合较好。为与Sesame库的数据相比，计算了3条等温线，发现在我们关心的区域，计算值与Sesame库的计算结果比较接近。

2 理论模型简介

考虑一个有t种可能组分的系统，包括s种气体，p种固体物质，它们由c种不同的化学元素组成。对一化学反应，反应物中各元素的原子个数与生成物中各元素的原子个数相等，即满足化学反应质量守恒。由分子式矢量表示得

式中：aij表示第j种元素在第i种组分分子式中的原子个数，ni为第i种组分的物质的量，qj是系统中第j种元素的总的物质的量。

化学平衡要求反应物的化学势必须等于生成物的化学势。化学平衡时自由能最小，据此可以导出化学平衡方程组［2］。

炸药爆轰产物混合系统的Helmholtz自由能可写为［3－4］

式中：Af、Ac分别为流体和固体的Helmholtz自由能，nc为固体的物质的量，｛ni｝表示流体中各组元的物质的量。流体自由能Af又可表示为不考虑分子间相互作用贡献的相应理想气体的自由能Aid和体系分子间相互作用贡献的超额自由能 Aint之和［5－6］，即

对于理想气体部分的Helmholtz自由能Aid，采用基于Rossiwi（单原子）、Pennigton和Kobe（双原子）、谐振子－刚性转子近似（多原子）等理论计算给出的拟合函数结果［1］。

假定流体中各组元分子间相互作用势为

则根据各组元浓度（xi＝ni／）和分子间相互作用参数（εij，αij，），采用van der Waals等效单组分流体模型（vdw1f），可将混合物等效为一元流体，等效分子势φ（r）中的参数（ε，α，r＊）可以通过下列推广的混合法则计算［7］，即

而（εij，αij）采用推广的Lorent-Berthelor组合规则求出［8］，即

式中：lij、mij、kij为修正系数，除kH2OCO2＝0.965，kH2ON2＝1.03外，其余lij＝mij＝kij＝1［9］。用求出的炸药的爆轰产物组分ni，结合式（5）、（6），可以求出vdwlf模型中等效分子势φ（r）中的参数（ε，α，r＊）。

采用Ross硬球微扰理论软球修正模型，计算炸药爆轰气相产物体系等效单组分流体分子间相互作用贡献的超额自由能Aint。

炸药爆轰气相产物分子间相互作用贡献为

式中：AHS、gPY（r／d，η）和F12（η）分别为硬球超额自由能、Percus-Yevick硬球径向分布函数和软球修正项，η＝πρd3／6，ρ＝N／V，d是硬球直径，d的选取应使式（7）右边最小，并将此最小值作为体系当前状态下的Helmholtz自由能Aint。

碳在高压下的状态方程在冲击波物理中非常重要，由于爆轰是一个瞬态的高温高压过程，很难进行直接观察和测量，目前还没有一种方法能够给出爆轰CJ点上的爆轰产物的真实组成。

在早期的研究中几乎都把爆轰产物的碳作为石墨处理，后来把爆轰产物中的碳当作石墨或金刚石的逐渐增多。最近，L.E.Fried等［9］对碳的4种相态（石墨、金刚石、类石墨液碳、类金刚石液碳）作了详细的研究，研究结果与 M.Togaya［10］和F.P.Bundy等［11］的实验结果符合得很好。文献［9］给出了碳的4种相的Gibbs自由能具体计算公式（适用范围：0≤p≤600GPa，300K≤T≤15kK）

式中：p0＝100kPa。由式（8），计算出给定（p，T）下碳的石墨相、金刚石相、类石墨液相和类金刚石液相4种相态的Gibbs自由能，由Gibbs自由能最小原理，从4种相态中选出炸药爆轰产物中游离态的碳最可能存在的相态，根据热力学性质计算其他的热力学量。相比于Murnaghan状态方程，L.E.Fried在碳的4相状态方程中将温度引了进来，因此更符合实际情况。

根据热力学知识，体系的压强p、能量E、熵S和各产物的化学势μi可由总的自由能A求出

具体的表达式可以参考文献［1］。

3 计算结果与分析

采用类CHEQ程序计算了初始密度ρ0＝1.894g／cm3的高含碳炸药PBX9502在CJ点的爆轰参数，计算结果见表1。从表1可以看到，本文中和BKW计算的爆速都与实验值符合很好（小于国际最好标准5%），但BKW状态方程中有4个常数，这些常数通过拟合产物组分的Hugoniot线和炸药爆轰实验的有关数据而得到。而本文中仅用到各产物的分子间势，采用统计原理计算的爆轰参数，其中没有与炸药爆轰实验相关的可调参数。

表1 PBX9502炸药爆轰CJ点的爆轰参数Table 1 Detonation parameters at the CJ point for PBX9502explosive

近代，由于工程技术的需要，科学家们建立了许多状态方程数据库，著名的有美国的Sesame库和俄罗斯的DBMS库。为进一步验证爆轰产物物态方程的正确性，在图1中比较了1 160、2 320和5 802K等3条等温线上本文中和Sesame库计算的压力随密度的变化。从图1中可以看到，本文计算结果与Sesame库的计算结果比较接近。

在图1中，本文中计算的5 802K等温线上发现，压强有一突变。由于气相产物状态方程是连续变化的，不会引起突变，所以可能是固态产物碳的状态方程发生了改变。在图2中，将游离态的碳设为单一的石墨相（Graphite）、金刚石相（Diamond）、类石墨液碳相（Liquid 1）和类金刚石液碳相（Liquid 2），并与碳的4相态状态方程（CHEQ）相比，可以看出，在密度小于2.9g／cm3时，4相态状态方程与类金刚石液碳相计算出的压强相等；密度大于2.9g／cm3时，与金刚石相计算出的压强相等。因此，在密度为2.9g／cm3附近，游离态的碳由类金刚石液碳转化为固态金刚石，从而引起热力学量的突变。

图1 等温线上压力随密度的变化Fig.1Isothermal pressure-density curves calculated by the two different methods

图2 碳取不同的相态时等温线上压力随密度的变化Fig.2 Isothermal pressure-density curves in the case of caebon at different states

碳的相变，不仅引起压强的突变，还引起各爆轰产物组分的变化。在图3中，列举了5 802K等温线上主要产物随密度的变化。从图3可以看到，在5 802K等温线上，随着产物密度的增加，炸药中的N元素主要以N2的形式存在，H元素主要以H2O的形式存在，H2O含量的增加，致使炸药中与C元素结合的O元素减少，CO逐渐减少，游离态的碳析出。这些曲线随密度的变化应是光滑的，但在密度为2.9g／cm3附近，都有一突变，这是由于碳的相变引起的。

图3 等温线上各产物物质的量随密度的变化Fig.3 Variation of the amounts of substance with density along isothermal curves for detonation products

4 结 论

（1）本文中计算的CJ点的爆轰参数与实验值符合很好。

（2）本文中计算的等温线与Sesame库的计算结果比较接近。

（3）本文中计算的5 802K等温线上发现了一个拐点，分析发现是因为在此处游离态的碳发生了相变，由类金刚石液碳转化为固态金刚石，从而引起热力学量及产物组分的突变。

［1］赵艳红，刘海风，张弓木.基于统计物理的爆轰产物物态方程研究［J］.物理学报，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］刘海风，陈栋泉，张世泽.爆轰产物物态方程及CHBr3相变的理论研究［J］.高压物理学报，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］.高压物理学报，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］刘福生，陈先猛，陈攀森，等.液态CO2高温高密度状态方程研究［J］.高压物理学报，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］.高压物理学报，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］.四川师范大学学报（自然科学版），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

（责任编辑 曾月蓉）

O381 国标学科代码：130·35

A

1001-1455（2010）06-0647-05

2009-08-20；

2009-11-06

中国工程物理研究院科学技术发展基金项目（2008B0201019）

赵艳红（1977— ），女，硕士，助理研究员。