炸藥爆轟熱力學(xué)數(shù)據(jù)非線性擬合方法研究?

杜明燃 張金龍 梁琳娜 黃亮亮 韓體飛 劉 鋒

安徽理工大學(xué)化學(xué)工程學(xué)院(安徽淮南，232001)

引言

研究炸藥性能時，爆轟參數(shù)一直都是備受關(guān)注的課題，它是預(yù)測炸藥性能和理論設(shè)計炸藥的基礎(chǔ)[1]。為計算炸藥爆轟參數(shù)，國內(nèi)外學(xué)者提出了BKW、JWL、LJD、JCZ3和 WCA 等物理模型[2-7]，并發(fā)展了相關(guān)計算程序，對理論研究爆轟參數(shù)作出了重大貢獻(xiàn)。研究發(fā)現(xiàn)，利用BKW[3]、WCA和LJD等狀態(tài)方程計算爆轟參數(shù)時，爆轟參數(shù)的求解過程是一個反復(fù)迭代的過程，求解爆溫的熱力學(xué)參數(shù)若滿足一種函數(shù)關(guān)系，對設(shè)計程序自循環(huán)求解爆轟參數(shù)具有重要意義。

以往研究者在研究爆溫時多采用Kast平均熱容法[8-9]，這種計算方法基于平均熱容建立，在精確計算中存在不足，導(dǎo)致爆轟參數(shù)計算得不準(zhǔn)確。從理論角度分析，能量守恒法[9]計算爆溫具有較高的精度，然而，由于理論數(shù)據(jù)的缺失，能量守恒最終確定的爆轟溫度一般根據(jù)理論數(shù)據(jù)估算得出；同時，在編程計算中發(fā)現(xiàn)，爆轟產(chǎn)物一般十余種，若不對熱力學(xué)參數(shù)進(jìn)行處理，能量守恒法需要反復(fù)迭代大量數(shù)據(jù)，程序繁瑣[10]且容易出錯。因此，利用擬合的方法確定一種函數(shù)關(guān)系變得非常必要，在此基礎(chǔ)上再建立求解爆轟參數(shù)的數(shù)學(xué)模型，對求解爆轟參數(shù)及爆炸物品的應(yīng)用具有一定參考意義。

本文中，采用非線性擬合的方法處理爆轟產(chǎn)物熱力學(xué)數(shù)據(jù)，比較非線性擬合法與Kast平均熱容法計算物質(zhì)的焓變(HT-H298)與理論值的差別。在保證非線性擬合方程精度的基礎(chǔ)上，依據(jù)此非線性擬合方程建立求解爆轟參數(shù)的數(shù)學(xué)模型，模型選取BKW狀態(tài)方程作為爆轟產(chǎn)物的狀態(tài)方程，利用最小自由能原理和Hugoniot關(guān)系，最終編程實(shí)現(xiàn)此過程。利用此程序，求解RDX和硝基甲烷(CH3NO2)的爆溫、爆壓、爆速和爆轟產(chǎn)物組成，驗證此方法的可行性。

1 炸藥爆轟溫度的確定方法

1.1 Kast平均熱容法

計算炸藥爆轟溫度常用的方法是爆轟產(chǎn)物平均熱容法[8]。爆轟產(chǎn)物平均熱容法是利用在一定溫度范圍內(nèi)物質(zhì)的平均熱容計算爆轟溫度，物質(zhì)的平均熱容由經(jīng)驗公式確定。爆轟產(chǎn)物平均熱容法的計算公式如下：

式中：QV為炸藥定容爆熱，kJ/mol；為爆轟產(chǎn)物的平均熱容；T為爆溫，K；T0=298.15 K；t為溫差；A、B為平均熱容的線性系數(shù)。

1.2 能量守恒法

能量守恒法[9]利用物質(zhì)能量守恒確定爆溫，能量守恒法滿足如下計算公式：

1.3 非線性擬合法

能量守恒法具有嚴(yán)格的理論依據(jù)，但能量守恒法與Kast平均熱容法相比也具有一個重要缺陷：理論數(shù)據(jù)不足，導(dǎo)致部分計算只能停留在估算階段，給爆溫的確定帶來誤差，導(dǎo)致計算其他爆轟參數(shù)也出現(xiàn)誤差。因此，依據(jù)能量守恒法建立模擬物質(zhì)的熱力學(xué)關(guān)系的方程變得十分必要。依據(jù)式(4)、式(5)和式(6)，可以得出以下假設(shè)公式：

式中：A1、A2、A3、A4為待擬合系數(shù)，可以根據(jù)最小二乘法編程計算得到，編程依據(jù)熱化學(xué)參數(shù)可由文獻(xiàn)[11]獲得。

式(7)、式(8)確定后，可以依據(jù)式(3)建立確定爆溫的能量守恒公式。表1給出了爆轟產(chǎn)物對應(yīng)系數(shù)的擬合結(jié)果，爆轟產(chǎn)物中的碳單質(zhì)視為金剛石[12-18]。表1適應(yīng)范圍為298～5 000 K。為比較線性擬合法與平均熱容法計算方法的準(zhǔn)確度，分別利用Kast平均熱容法和非線性擬合法計算CO2、H2O和N2的焓變(HT-H298)，共7組數(shù)據(jù)。表2列出兩種方法計算所得的3種主要產(chǎn)物的焓變(HT-H298)與理論值[11]。

表2顯示，對于CO2、H2O和N2，Kast平均熱容法的計算誤差遠(yuǎn)遠(yuǎn)大于本工作所采用的非線性擬合法。表2中，非線性擬合法的最大誤差為0.5%；Kast平均熱容法的計算最小誤差為3.1%，最大誤差為12.5%。

2 求解爆轟參數(shù)過程

2.1 BKW狀態(tài)方程

BKW狀態(tài)方程[10]最初是由Becker提出的，后來經(jīng)過Kistiakowsky和Wilson等的修正，最后確定為以下形式：

表1 式(7)和式(8)系數(shù)的擬合值Tab.1 Fitting values of coefficient for Eq-7 and Eq-8

表2 焓變(HT-H298)的計算與理論比較Tab.2 Calculation results and theoretical values of(HT-H298) kJ/mol

式中：p為爆轟壓力，Pa；β為常數(shù)；X=κ∑xiki/V(T＋θ)α；V為爆轟氣體摩爾體積，κ、α、θ為經(jīng)驗常數(shù)，ki為第i種物質(zhì)的余容，xi為這種物質(zhì)的摩爾分?jǐn)?shù)。

2.2 求解爆轟參數(shù)的數(shù)學(xué)模型及模擬算法

利用最小自由能原理建立非線性方程組。根據(jù)質(zhì)量守恒和物質(zhì)自由能公式，再運(yùn)用拉格朗日乘數(shù)法把條件極值轉(zhuǎn)化為非條件極值，可以得出與爆轟產(chǎn)物系統(tǒng)具有相同極值條件的公式：

式中：G、g、g0分別為爆轟體系自由能、物質(zhì)自由能和標(biāo)準(zhǔn)狀況下自由能，kJ/mol；Nk為1 kg爆炸產(chǎn)物含有第k種元素的物質(zhì)的量；λ為拉格朗日因子；k為第k種元素序號；λk為第k種元素的拉格朗日因子；i為第i種氣體產(chǎn)物序號；j為第j種固體產(chǎn)物序號；aik為每摩爾第i種氣相組分含有k元素的物質(zhì)的量，mol；ajk為每摩爾第j種凝聚相組分含有k元素的物質(zhì)的量，mol；ng為氣體產(chǎn)物總物質(zhì)的量；nig為第i種氣體的物質(zhì)的量；njs為第j種固體的物質(zhì)的量；l為炸藥所含元素的種類數(shù)；m和n分別為氣體和固體的種類。

式(10)對應(yīng)極值條件的方程組為非線性方程組，為求解此方程組，根據(jù)本小組研究的牛頓迭代法求解非線性方程組的思路[19]，利用以下模型求解爆轟參數(shù)和爆轟產(chǎn)物組成。

首先，需要把 ?G(n，λk)/?nig、?G(n，λk)/?njs和?G(n，λk)/?λk看作是關(guān)于nig、njs和λk的方程組，令

把式(11)、式(12)、式(13)和式(14)代入式(16)，可求得F′(x)，進(jìn)一步求得F′(x)-1，并記解為x(k＋1)，對應(yīng)牛頓迭代法公式：x(k＋1)=xk-F′(xk)-1F(xk)，求解上述公式要先假設(shè)一組爆溫T′和爆壓p′，此溫度和壓力可以根據(jù)一般經(jīng)驗計算值確定，再取經(jīng)驗計算所得爆炸產(chǎn)物組分量并作一定變換后，作為初始近似根x1，設(shè)置計算精度ε=0.000 001，即xk和x(k＋1)對應(yīng)元素的誤差絕對值。

炸藥爆轟氣相產(chǎn)物分子體系必須滿足Hugoniot關(guān)系[10，20]：

爆速表達(dá)式為

根據(jù)BKW狀態(tài)方程、以上求解模型和本文中設(shè)計的爆轟產(chǎn)物非線性擬合方程，計算爆溫T和爆壓p，與假設(shè)爆溫T′和爆壓p′比較，反復(fù)修正，直到滿足假設(shè)值與計算所得值相同。

以上求解爆轟參數(shù)的數(shù)學(xué)模型利用自編程序?qū)崿F(xiàn)，計算過程詳見本課題組前期研究工作[19]，其程序算法框圖見圖1。

圖1 計算模擬程序圖Fig.1 Computer simulation program chart

3 爆轟參數(shù)計算結(jié)果與分析

利用非線性擬合方法求得密度為1.80 g/cm3的RDX和密度為1.14 g/cm3的CH3NO2的爆轟產(chǎn)物組成，如表3所示，爆轟參數(shù)如表4所示。

計算過程BKW狀態(tài)方程參數(shù)[10]α取0.5，β取0.16，θ取400，κ取10.91。碳的單質(zhì)生成物為金剛石[12-18]。RDX和CH3NO2的生成熱分別是-70.66 kJ/mol和 113.09 kJ/mol。

為對比本文中計算方法的精確性，表4還分別給出爆轟參數(shù)的實(shí)驗結(jié)果和依據(jù)Kast計算爆溫所得結(jié)果。

表4顯示，采用非線性擬合法作為計算爆溫的方法計算爆轟參數(shù)，比使用Kast平均熱容法計算的數(shù)值更加接近實(shí)驗值。

對于RDX，非線性擬合法計算爆溫、爆壓和爆速與實(shí)驗值的誤差分別為3.3%、0.4%和0.3%；Kast平均熱容法計算爆溫、爆壓和爆速與實(shí)驗值的誤差分別為6.4%、3.4%和2.0%。

表3 每摩爾炸藥CJ點(diǎn)爆轟產(chǎn)物組成Tab.3 Detonation products of explosives at CJ point mol

表4 RDX和CH3NO2的爆溫、爆壓和爆速Tab.4 Detonation temperature，pressure and velocity of RDX and CH3NO2

對于CH3NO2，非線性擬合計算的爆溫值也比Kast平均熱容法更接近實(shí)驗值。

4 結(jié)論

1)采用最小二乘法法得到爆轟產(chǎn)物熱力學(xué)參數(shù)的非線性擬合方程，與Kast平均熱容法比較，非線性擬合法在計算物質(zhì)焓變(HT-H298)時顯示出更好的計算精確性。

2)在此基礎(chǔ)上，利用自編程序?qū)崿F(xiàn)爆轟產(chǎn)物組成、爆壓、爆溫和爆速的計算。密度為1.80 g/cm3的RDX和密度為1.14 g/cm3的CH3NO2的爆轟參數(shù)計算結(jié)果顯示，采取非線性擬合法計算的爆轟參數(shù)比采用Kast平均熱容法計算的結(jié)果更加接近實(shí)驗值，顯示出更好的應(yīng)用前景，對炸藥爆轟性能預(yù)測具有一定參考價值。