胡榮祖 趙鳳起 高紅旭 姚二崗 張 海 王 堯常象宇 趙宏安

(1西安近代化學(xué)研究所燃燒與爆炸技術(shù)重點(diǎn)實(shí)驗(yàn)室,西安710065;2西北大學(xué)數(shù)學(xué)系數(shù)據(jù)分析和計(jì)算化學(xué)研究所,西安710069;3西安交通大學(xué)信息科學(xué)與系統(tǒng)科學(xué)研究所,西安710049;4西北大學(xué)信息科學(xué)與工程學(xué)院,西安710069)

1 引言

放熱分解反應(yīng)體系熱爆炸的臨界溫升速率(dT/dt)Tb,是評(píng)價(jià)含能材料(EMs)安定性和安全性的重要參數(shù).在表達(dá)這個(gè)參數(shù)方面,胡榮祖等從反應(yīng)進(jìn)度和能量變化的關(guān)系,1熱分解過(guò)度到熱爆炸的充分必要條件2和非等溫反應(yīng)的動(dòng)力學(xué)方程3導(dǎo)出了絕熱及近似絕熱條件下,速率方程為dα/dt=A1exp(-Ea1/RT)(1-α)m+A2exp(-Ea2/RT)αn(1-α)p的表觀經(jīng)驗(yàn)級(jí)數(shù)和m=n=p=1的一級(jí)自催化分解反應(yīng)體系(dT/dt)Tb值的估算式;3,4王耘等5導(dǎo)出了絕熱和近似絕熱條件下反應(yīng)體系的dT/dt值的估算式,但對(duì)(dT/dt)Tb值的估算式未作研究.在估算這個(gè)參數(shù)方面,Hu等6-10估算了氮含量分別為11.92%、11.97%、13.54%、13.86%和14.14%的硝化棉(NC)的一級(jí)自催化分解反應(yīng)體系的(dT/dt)Tb值;Eisenreich等11,12估算了NC一級(jí)自催化分解反應(yīng)的動(dòng)力學(xué)參數(shù).在實(shí)測(cè)這個(gè)參數(shù)方面,胡榮祖等13研究了硝仿熱爆炸的爆前加熱、熱分解溫升和時(shí)間的關(guān)系.本工作作為文獻(xiàn)3-12的拓展,報(bào)道了自催化反應(yīng)速率方程分別為:dα/dt=Aexp(-E/RT)α(1-α),dα/dt=Aexp(-E/RT)(1-α)n(1+Kcatα),dα/dt=Aexp(-E/RT)[αa-(1-α)n],dα/dt=A1exp(-Ea1/RT)(1-α)+A2exp(-Ea2/RT)α(1-α),dα/dt=A1exp(-Ea1/RT)(1-α)m+A2exp(-Ea2/RT)αn(1-α)p,dα/dt=Aexp(-E/RT)·(1-α),dα/dt=Aexp(-E/RT)(1-α)n,dα/dt=A1exp(-Ea1/RT)+A2exp(-Ea2/RT)(1-α)和 dα/dt=A1exp(-Ea1/RT)+A2exp(-Ea2/RT)α(1-α)的(dT/dt)Tb估算式的導(dǎo)出途徑和NC(13.54%N)自催化分解反應(yīng)動(dòng)力學(xué)參數(shù)和自催化分解轉(zhuǎn)向熱爆炸時(shí)(dT/dt)Tb值的估算結(jié)果.

2 理論和方法

單位時(shí)間內(nèi)由于EMs熱分解而放出的熱量q1為

式中,Q為熱分解反應(yīng)的焓(J·mol-1),V為EMs的裝填體積(cm3),d為裝填密度(g·cm-3),M為EMs的摩爾質(zhì)量(g·mol-1),dα/dt為機(jī)理函數(shù)為f(α)=α(1-α)的自催化反應(yīng)(Au)速率:

聯(lián)立方程(1)和(2),得

與此同時(shí),單位時(shí)間內(nèi)因傳熱由反應(yīng)區(qū)通過(guò)器壁向四周環(huán)境散失的熱量q2為

式中,k′為傳熱系數(shù)(J·cm-2·K-1·s-1),S表示藥柱表面積(cm2);Tc為按照線性關(guān)系Tc=T0+βt確定的反應(yīng)器壁和空間溫度;β表示線性加熱速率(K·min-1),T0表示熱分析曲線離開(kāi)基線的溫度(K).

熱爆炸發(fā)生時(shí),方程(3)變?yōu)?/p>

此處,αb是相應(yīng)于Tb的α值,Tb是EMs的熱爆炸溫度(K).

方程(4)變?yōu)榇颂?Te0是β→0時(shí)的onset溫度.

根據(jù)Semenov的熱爆炸理論,2熱分解過(guò)渡到熱爆炸的充分必要條件式為

此處(dT/dt)Tb是熱分解轉(zhuǎn)向熱爆炸時(shí)EMs中的臨界溫升速率.

方程(4)對(duì)T微分,得

方程(14)稱為Au自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得E、A、Te0、Tb和αb,就可從方程(14)得到(dT/dt)Tb的值.

類似地,由機(jī)理函數(shù)為 f(α)=(1-α)n(1+Kcatα))速率方程

方程(18)稱為CnB自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得E,A,n,Kcat,Te0,Tb和αb,就可從方程(18)得(dT/dt)Tb值.

由機(jī)理函數(shù)為 f(α)=αa-(1-α)n的自催化反應(yīng)(Bna)的速率方程方程(22)為Bna自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得E,A,Te0,n,a,Tb和αb,就可從方程(22)得(dT/dt)Tb值.

由一級(jí)自催化分解反應(yīng)速率方程

式中,k1=A1exp(-Ea1/RT),k2=A2exp(-Ea2/RT).

方程(26)稱為一級(jí)自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得Ea1,Ea2,A1,A2,Te0,Tb和αb,就可從方程(26)得(dT/dt)Tb值.

由經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)速率方程

方程(30)稱為表觀經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得Ea1,Ea2,A1,A2,Te0,Tb,αb,m,n和p,就可從方程(30)得(dT/dt)Tb值.

由簡(jiǎn)單一級(jí)自催化分解反應(yīng)速率方程

方程(34)為簡(jiǎn)單一級(jí)自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得E、A、Te0和Tb,就可從方程(34)得(dT/dt)Tb值.

由簡(jiǎn)單n級(jí)自催化分解反應(yīng)速率方程

方程(38)稱為簡(jiǎn)單n級(jí)自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得E,A,Te0,n,Tb和αb,就可從方程(38)得(dT/dt)Tb值.

由經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1的自催化分解反應(yīng)速率方程

方程(40)稱表觀經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1的自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得Ea1,Ea2,A1,A2,Te0,Tb和αb,就可從方程(40)得(dT/dt)Tb值.

由經(jīng)驗(yàn)級(jí)數(shù)m=0,n=1,p=1的自催化分解反應(yīng)速率方程

方程(42)稱為表觀經(jīng)驗(yàn)級(jí)數(shù)m=0,n=1,p=1的自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率表達(dá)式.

一旦從熱流曲線解得Ea1,Ea2,A1,A2,Te0,Tb和αb,就可從方程(42)得(dT/dt)Tb值.

為求方程(2)和(31)中的2參數(shù)(A,E),方程(35)中的3參數(shù)(A,E,n),方程(15)、(19)、(23)、(39)和(41)中的4參數(shù)(分別為A,E,n,Kcat;A,E,a,n;A1,A2,Ea1,Ea2;A1,A2,Ea1,Ea2;A1,A2,Ea1,Ea2)和方程(27)中的7參數(shù)(A1,A2,Ea1,Ea2,m,n,p),進(jìn)行了線性最小二乘法和信賴域方法14的數(shù)值模擬.

對(duì)n=1的自催化反應(yīng)速率方程,有

由方程(43)和(44)的左端項(xiàng)與1/T作圖,用最小二乘法從斜率得E,截距得A.

經(jīng)非線性優(yōu)化模型:

對(duì)n級(jí)自催化反應(yīng)速率方程,有

對(duì)CnB速率方程,有

對(duì)Bna速率方程,有

對(duì)一級(jí)自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1的自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=1,p=1的自催化分解反應(yīng)速率方程,有

由最小化均方誤差估計(jì)方程(45)中的3參數(shù)(A,E,n),方程(46)、(47)、(48)、(50)和(51)中的 4 參數(shù)(A,E,n,Kcat;A,E,a,n;A1,A2,Ea1,Ea2;A1,A2,Ea1,Ea2;A1,A2,Ea1,Ea2)和方程(49)中的7參數(shù)(A1,A2,Ea1,Ea2,m,n,p).

為準(zhǔn)確估算NC(13.54%N)放熱分解反應(yīng)的動(dòng)力學(xué)參數(shù),進(jìn)行了邊界約束,限制參數(shù)范圍,將方程(45)-(51)分別改寫(xiě)如下:

對(duì)n級(jí)自催化反應(yīng)速率方程,有

對(duì)CnB速率方程,有

對(duì)一級(jí)自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1的自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=1,p=1的自催化分解反應(yīng)速率方程,有

為敘述簡(jiǎn)便,將具有邊界約束的非線性優(yōu)化問(wèn)題重新描述為

其中,對(duì)n級(jí)自催化反應(yīng)速率方程,有

對(duì)CnB速率方程,有

對(duì)Bna速率方程,有

對(duì)一級(jí)自催化分解反應(yīng)速率,有

u=(1018,1018,106,106)T,i=1,2,…,N

對(duì)經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1的自催化分解反應(yīng)速率方程,有

對(duì)經(jīng)驗(yàn)級(jí)數(shù)m=0,n=1,p=1的自催化分解反應(yīng)速率方程,有

定義矩陣D(x)為D(x)? diag[|vi(x)|-1/2],即,D-2是一對(duì)角矩陣且其第i個(gè)對(duì)角元素為|vi(x)|.考慮非線性方程組

設(shè)xk∈int(Φ),式(59)的牛頓迭代滿足

其中,M(x)?B(x)+C(x),C(x)?D(x)diag[g(x)]Jv(x)D(x),

業(yè)已證明,xk是式(59)的局部最小值點(diǎn),當(dāng)且僅當(dāng)sk是下式的解

這里Δk為正數(shù),表示信賴域的大小.因此,通過(guò)求解子問(wèn)題(70)來(lái)確定xk是否為目標(biāo)函數(shù)的最小值點(diǎn)是合理的,進(jìn)而,當(dāng)信賴域大小Δk足夠大時(shí),局部最小值的鄰域中式(67)的牛頓迭代正好是信賴域子問(wèn)題(70)的解.

解決具有邊界約束的非線性優(yōu)化問(wèn)題的信賴域方法模型為:

取x0∈int(Φ),循環(huán)k=0,1,…

(1)計(jì)算f(xk),gk,Hk及Ck;

(2)根據(jù)式(70)計(jì)算sk且滿足xk+sk∈int(Φ).

(3)計(jì)算

(5)修正模型yk,矩陣Dk及信賴域大小Δk.

信賴域大小Δk的修正為:給定0＜μ＜η＜1,γ1＜1＜γ2及Λl＞0,

3 計(jì)算實(shí)例

依據(jù)表1和表2中NC(13.54%N)放熱分解反應(yīng)的T-轉(zhuǎn)化率(α)和β-Te關(guān)系數(shù)據(jù),用方程(14),(18),(22),(26),(30),(34),(38),(40)和(42)分析表1和表2中數(shù)據(jù),得表3中Au,CnB,Bna,一級(jí)/表觀經(jīng)驗(yàn)級(jí)數(shù)、簡(jiǎn)單一級(jí)和n級(jí)及表觀經(jīng)驗(yàn)級(jí)數(shù)m=0,n=0,p=1和m=0,n=1,p=1的自催化分解反應(yīng)動(dòng)力學(xué)參數(shù)和自催化分解反應(yīng)轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率(dT/dt)Tb值.據(jù)此,發(fā)現(xiàn)9個(gè)微分方程中Au速率方程[方程(2)]算得的E和A值在EMs分解反應(yīng)動(dòng)力學(xué)參數(shù)正常范圍[E=80-250 kJ·mol-1,lg(A/s-1)=7-30]17外,簡(jiǎn)單一級(jí)自催化分解反應(yīng)方程[方程(31)]求E和A的線性相關(guān)系數(shù)為0.9157,遠(yuǎn)離0.98,所得(dT/dt)Tb為負(fù)值,屬不合理,被排除.同時(shí)滿足E和A值在正常范圍內(nèi),(dT/dt)Tb為正值的方程有7個(gè),在這7個(gè)方程中,以計(jì)算值與實(shí)驗(yàn)值的相對(duì)誤差為判據(jù),宜做描述NC(13.54%N)自催化分解反應(yīng)過(guò)程的方程依次為:方程(27)、(23)、(41)、(39)、(19)、(35)、(15),對(duì)應(yīng)誤差依 次 為:0.0206、0.0255、0.0299、0.109、0.1698、0.1932、0.7620,相應(yīng)熱爆炸臨界溫升速率依次為:0.103、0.169、0.189、0.00002、0.615、0.020、0.266 K·s-1.由此認(rèn)為,用相對(duì)誤差最小的經(jīng)驗(yàn)級(jí)數(shù)自催化反應(yīng)速率方程描述NC(13.54%N)熱分解過(guò)程是可取的,熱分解轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率為0.103 K·s-1是可接受的.

表1 用DSC測(cè)得的NC(13.54%N)的熱分解數(shù)據(jù)8Table 1 Thermal decomposition data8of NC(13.54%N)determined by DSC

表2 NC(13.54%N)熱爆炸臨界溫度(Tb)的計(jì)算值Table 2 Calculated values of the critical temperature(Tb)of thermal explosion for NC(13.54%N)

表3 用方程(14),(18),(22),(26),(30),(34),(38),(40)和(42)分析表1和表2中數(shù)據(jù)的結(jié)果Table 3 Results for analyzing the data in Tables 1 and 2 by Eqs.(14),(18),(22),(26),(30),34),(38),(40),and(42)

4 結(jié)論

(1)提出了從不同升溫速率條件下的DSC曲線數(shù)據(jù)計(jì)算/確定EMs自催化反應(yīng)的動(dòng)力學(xué)參數(shù)和自催化分解轉(zhuǎn)向熱爆炸時(shí)的臨界溫升速率(dT/dt)Tb的方法.

(2)非等溫條件下NC(13.54%N)熱分解過(guò)程可用表觀經(jīng)驗(yàn)級(jí)數(shù)自催化分解反應(yīng)動(dòng)力學(xué)方程描述:

(3)NC(13.54%N)自催化分解反應(yīng)過(guò)渡到熱爆炸時(shí)的臨界溫升速率值為0.103 K·s-1.

(1) Deng,Y.Chem.J.Chin.Univ.1985,6(7),621.[鄧 郁.高等學(xué)?；瘜W(xué)學(xué)報(bào),1985,6(7),621.]

(2) Semenov,N.N.On Some Problem of Kineties and Reactivity;Ind.,AN SSSR:Moscow,1958;p 421.

(3) Hu,R.Z.;Gao,S.L.;Zhao,F.Q.;Shi,Q.Z.;Zhang,T.L.;Zhang,J.J.Thermal Analysis Kinetics;Science Press:Beijing,2008;pp 322-334.[胡榮祖,高勝利,趙鳳起,史啟楨,張同來(lái),張建軍.熱分析動(dòng)力學(xué).北京:科學(xué)出版社,2008:322-334.]

(4) Hu,R.Z.;Zhang,H.;Xia,Z.M.;Guo,P.J.;Gao,S.L.;Shi,Q.Z.;Lu,G.E.;Jiang,J.Y.Chin.J.Energ.Mater.2003,11(3),130.[胡榮祖,張 海,夏志明,郭鵬江,高勝利,史啟楨,路桂娥,江勁勇.含能材料,2003,11(3),130.]

(5)Wang,Y.;Feng,C.G.;Zheng,X.Chin.J.Energ.Mater.2000,8(3),119.[王 耘,馮長(zhǎng)根,鄭 曉.含能材料,2000,8(3),119.]

(6) Hu,R.Z.;Guo,P.J.;Song,J.R.;Zhang,H.;Xia,Z.M.;Ning,B.K.;Fang,Y.;Shi,Q.Z.;Liu,R.;Luo,G.E.;Jiang,J.Y.Chin.J.Explo.Prop.2003,26(2),53.

(7) Ning,B.K.;Hu,R.Z.;Zhang,H.,Xia,Z.M.;Guo,P.J.;Liu,R.;Luo,G.E.;Jiang,J.Y.Thermochim.Acta 2004,416(1-2),47.doi:10.1016/j.tca.2003.11.029

(8)Hu,R.Z.;Guo,P.J.;Gao,S.L.;Zhang,H.;Xia,Z.M.;Ning,B.K.;Fang,Y.;Shi,Q.Z.;Liu,R.Chinese Journal of Polymer Science 2003,21(3),285.

(9) Guo,P.J.;Hu,R.Z.;Zhang,H.,Xia,Z.M.;Song,J.R.;Gao,S.L.;Ning,B.K;Liu,R.;Lu,G.E.;Jiang,J.Y.Chem.Res.Chin.Univ.2004,20(2),163.

(10)Zhang,H.;Xia,Z.M.;Guo,P.J.;Hu,R.Z.;Gao,S.L.;Ning,B.K.;Fang,Y.;Shi,Q.Z.;Liu,R.Journal of Hazardous Materials 2002,94(3),205.doi:10.1016/S0304-3894(02)00118-8

(11) Eisenreich,N.Beitrag Zur Kinetic Thermischer Zersetzungseaktionen(ThermoanalytischeAuswertung Der Zersetzung Von Nitrocellulose).Ph.D.Dissertation,Technical University of Munich,Munich,1978.

(12) Eisenreich,N.;Pfeil,A.Thermochim.Acta 1983,61(1-2),13.doi:10.1016/0040-6031(83)80300-1

(13) Hu,R.Z.;Song,Q.C.;Xie,J.J.;Liang,Y.J.Explosion and Shock Waves 1987,7(3),217.[胡榮祖,松全才,謝俊杰,梁燕軍.爆炸與沖擊,1987,7(3),217.]

(14)Thomas,F.C.;Li,Y.Y.SIAM Journal on Optimization 1996,6(2),418.doi:10.1137/0806023

(15) Ozawa,T.Bulletin of the Chemical Society of Japan 1965,38(1),1881.

(16) Zhang,T.L.;Hu,R.Z.;Xie,Y.;Li F.P.Thermochim.Acta 1994,244(3),171.

(17)Hu,R.Z.;Yang,Z.Q.;Ling,Y.J.Thermochim.Acta 1988,123,135.doi:10.1016/0040-6031(88)80017-0