張坤, 智小琦, 肖游, 王帥

(1.中北大學(xué) 機電工程學(xué)院, 山西 太原 030051;2.中國兵器裝備集團自動化研究所有限公司 智能制造事業(yè)部, 四川 綿陽 621000;3.湖北航天化學(xué)技術(shù)研究所 航天化學(xué)動力技術(shù)重點實驗室, 湖北 襄陽 441003)

0 引言

研究彈藥在熱環(huán)境中的響應(yīng)特性對彈藥在壽命期的安全使用具有重要意義。因為彈藥無論是受熱刺激、強電磁輻射或二者的綜合作用,最終導(dǎo)致彈藥響應(yīng)的均是因含能材料受熱后的自持反應(yīng)所致。目前,模擬彈藥受熱刺激作用的研究方法主要是快速烤燃和慢速烤燃?？焖倏救寄M彈藥在大火作用下的響應(yīng)特性,慢速烤燃模擬彈藥在暗火作用下的響應(yīng)特性。研究手段主要采用試驗和仿真技術(shù)。試驗方法主要采用北大西洋公約組織標(biāo)準(zhǔn)化協(xié)定STANAG 4240(快速烤燃試驗)和STANAG 4382(慢速烤燃試驗),試驗研究內(nèi)容包括彈藥尺寸[1-2]、約束條件[3-4]、裝藥密度及間隙[5-6]、含能材料配方[7]、熔鑄炸藥點火機理[8]以及緩釋結(jié)構(gòu)特性[9]等對彈藥響應(yīng)溫度、響應(yīng)烈度的影響研究,已達到既深又廣的程度;仿真研究與試驗研究既相輔相成又互相補存和促進,能獲得點火時刻的溫度及其梯度、點火點位置、點火前壓力[10]、內(nèi)部機理的分析[11-12]及點火后的瞬時響應(yīng)特性[13]等內(nèi)容。其中關(guān)于點火后的含能材料自持反應(yīng)特性的仿真研究軟件,目前只有美國ALE3D軟件和Uintah軟件。

從上述情況可見,研究彈藥的熱刺激,無論從試驗還是仿真都取得了很大的成就。但是關(guān)于純理論研究點火點位置和點火溫度的文章,目前鮮有報道,而理論研究是彈藥熱刺激研究不可或缺的內(nèi)容,是對熱爆炸理論的完善和補存,具有重要的意義。本文擬采用疊加原理和分離變量方法,首次將炸藥非反應(yīng)性熱傳導(dǎo)與自熱反應(yīng)熱傳導(dǎo)拆分,以一維炸藥模型為基礎(chǔ),研究慢速烤燃的溫度分布及點火點位置的解析解,分析點火點位置的一維變化規(guī)律,以期為烤燃研究奠定理論基礎(chǔ)。

1 一維數(shù)學(xué)理論模型

1.1 理論推演的假設(shè)

為簡化理論計算,需對炸藥的烤燃過程作如下假設(shè):

1)炸藥為凝聚炸藥,炸藥在加熱過程中不發(fā)生相變;

2)彈體外壁均勻升溫,升溫速率為固定值;

3)忽略殼體與炸藥間的接觸熱阻;

4)炸藥的熱導(dǎo)率、比熱容、密度均為恒值;

5)炸藥活化能、指前因子、反應(yīng)熱均為常數(shù);

6)炸藥按照零階動力學(xué)定律分解。

1.2 慢速烤燃一維模型

炸藥的熱傳導(dǎo)方程用式(1)[14]描述:

(1)

(2)

當(dāng)則炸藥厚度L=2δ,其中δ為厚度值的一半。則一維無限平板模型的熱傳導(dǎo)可用式(3)來描述:

(3)

式中:α為熱擴散系數(shù),且α=λ/ρc;S(x,t)為熱傳導(dǎo)方程的源項,用式(4)表示:

(4)

式(3)的定解條件如下:

1)初始條件:T(x,0)=Tr,其中Tr為初始溫度;

2)邊界條件:T(±δ,t)=κt+Tr,其中κ為升溫速率。

2 求解

2.1 溫度場的拆分

當(dāng)沒有化學(xué)反應(yīng)時,炸藥熱傳導(dǎo)問題的解析解即非反應(yīng)性熱傳導(dǎo)比較容易求解。但是,具有化學(xué)反應(yīng)的熱傳導(dǎo)方程即反應(yīng)性熱傳導(dǎo)是極難獲得解析解的。為進一步求解熱傳導(dǎo)方程式(3),需要用到疊加原理[15],將原問題T(x,t)分解為非反應(yīng)性熱傳導(dǎo)項Tt(x,t)和自熱反應(yīng)熱傳導(dǎo)項Ts(x,t)兩個子問題。其中非反應(yīng)熱傳導(dǎo)式為

(5)

式(5)相應(yīng)的定解條件如下:

1)初始條件:Tt(x,0)=Tr;

2)邊界條件:Tt(±δ,t)=κt+Tr。

自熱反應(yīng)熱傳導(dǎo)的方程如下:

(6)

式(6)的定解條件如下:

1)初始條件:Ts(x,0)=0;

2)邊界條件:Ts(±δ,t)=0。

則原方程(式(3))的解T(x,t)可通過T(x,t)=Tt(x,t)+Ts(x,t)求得。

2.2 非反應(yīng)性熱傳導(dǎo)問題的求解

式(5)的熱傳導(dǎo)子問題不包含自熱源項,即屬于沒有化學(xué)反應(yīng)的惰性情況下的熱傳導(dǎo)問題。由于變量x、t相互獨立,使用分離變量法[16]將Tt(x,t)寫為

Tt(x,t)=X(x)U(t)

(7)

式中:X(x)為位置x的函數(shù);U(t)為時間t的函數(shù)。

將式(7)代入式(5),有如下等式:

(8)

(9)

式中:k為常數(shù)。

式(8)根據(jù)邊界條件Tt(±δ,t)=κt+Tr,通過Strum-Liouville理論求解有

Xn(x)=Cnsin[nπ(-x+δ)/(2δ)]

(10)

式中:Cn為常數(shù),X0(x)=C0;固有值λn=nπ,n=0,1,2,3,…。

對式(9)求解,有

Un(t)=Bne(-απ2n2t)/4δ2

(11)

式中:Bn為常數(shù)。

則非反應(yīng)性熱傳導(dǎo)方程式(5)的通解為

(12)

式中:An為常數(shù),且對于Tt(x,t)解函數(shù)首項Tt0(x,t)有Tt0(x,t)=A0。

代入初始條件Tt(x,0)=Tr,進一步求得非反應(yīng)性熱傳導(dǎo)方程的精確解Tt(x,t)為

(13)

隨著時間t的不斷增加,無窮級數(shù)項將逐漸減小,非反應(yīng)性熱傳導(dǎo)所產(chǎn)生的溫度分布將逐漸穩(wěn)定。其中初始溫度Tr作為常數(shù)項,不影響非反應(yīng)性熱傳導(dǎo)所產(chǎn)生的溫度分布。

2.3 自熱反應(yīng)熱傳導(dǎo)問題的求解

對于式(6)的自熱反應(yīng)熱傳導(dǎo)問題,非齊次項(即源項)為Arrhenius方程,直接求解是相當(dāng)復(fù)雜的。通過進一步將源項S(x,t)泰勒展開:

(14)

當(dāng)自熱反應(yīng)產(chǎn)生的溫度Ts小于非反應(yīng)性熱傳導(dǎo)溫度Tt的1%時,源項S(x,t)可取泰勒展開式首項近似代替,作近似后的自熱反應(yīng)熱傳導(dǎo)方程如下:

(15)

式中:

(16)

對于式(15)的求解,采用分離變量法[17]。將Ts(x,t)、S(x,t)分別按固有函數(shù)展開為Dn(t)、Sn(t),于是有

(17)

式中:

(18)

(19)

再將其反變換后,則有

(20)

這樣,可解得自熱反應(yīng)熱傳導(dǎo)方程的解Ts(x,t)為

Ts(x,t)=

(21)

觀察解的結(jié)構(gòu),不難發(fā)現(xiàn)式(21)中Ts(x,t)包含Arrhenius溫度積分,直接求出其精確的結(jié)果極其困難[18]。但Ts(x,t)中常數(shù)顯然不影響對其最高溫度所在位置的求解,即指數(shù)前因子A、反應(yīng)熱Q不影響自熱反應(yīng)熱傳導(dǎo)溫度最高值所在位置,進而不影響炸藥在烤燃中的點火點位置。隨著升溫速率κ的增加,炸藥兩側(cè)升溫加快,則點火點位置向邊界移動。

3 慢速烤燃試驗

3.1 試驗設(shè)計

設(shè)計尺寸為φ75 mm×295 mm烤燃彈,殼體材料選用45號鋼,兩端帶螺紋的端蓋連接,端蓋和殼體壁厚均為5 mm。內(nèi)部裝藥為甘肅銀光化學(xué)工業(yè)集團有限公司生產(chǎn)的RDX炸藥,單節(jié)藥柱直徑為57.4 mm,長度為47.5 mm,共6節(jié),藥柱總長285 mm,長徑比為5,裝藥密度為1 610 kg/m-3。

為獲得烤燃試驗過程中,炸藥內(nèi)部的溫度的變化情況,在每節(jié)藥柱中間刻槽埋入φ1 mm的上海自動化儀表股份有限公司生產(chǎn)的WRNK191型鎧裝電偶(精度為0.1 K,熱響應(yīng)時間≦3 s)。藥柱與藥柱之間以及藥柱與熱電偶之間均采用蟲膠漆粘接固定。熱電偶的溫度測點由上至下編號1～5號,兩相鄰測點及邊緣測點至端面的間距均為47.5 mm,溫度測點位置如圖1所示。由于采用等長壓制藥柱試驗,無法獲得其他位置的溫度情況。

圖1 溫度測點的位置示意圖Fig.1 Location of temperature measuring points

試驗彈試驗前狀態(tài)如圖2所示。試驗采用圓柱形保溫箱進行加熱,試驗裝置如圖3所示。

圖2 試驗彈體Fig.2 Tested ammunition

圖3 慢速烤燃試驗裝置Fig.3 Slow-cook test setup

使用日本島電公司生產(chǎn)的MR13溫控儀對烤燃彈外壁進行升溫控制,并使用美國Fluke公司的1586A多路測溫儀采集各測點的溫度數(shù)據(jù),采樣周期為1 s。從室溫305 K開始,以3.3 ℃/h的固定速率升溫,直至烤燃彈發(fā)生響應(yīng)。

3.2 試驗結(jié)果

試驗藥柱中5支熱電偶和外壁1支熱電偶在點火時刻的測量結(jié)果如表1所示,試驗中各測點溫度變化曲線如圖4所示。由于炸藥內(nèi)部的自熱反應(yīng)隨溫度的升高而逐漸加劇,臨近響應(yīng)時,各測點的溫度迅速升高。從5條測溫曲線可見,藥柱的兩側(cè)1號和5號熱電偶測得的溫度最高,3號點測得的溫度最低,可見點火區(qū)域位于藥柱兩側(cè),點火時外壁溫度為468.9 K。

表1 點火時刻熱電偶的測量結(jié)果

圖4 溫度變化曲線Fig.4 Temperature history curves

烤燃彈響應(yīng)后,試驗點產(chǎn)生較深的彈坑,加熱套筒與烤燃彈殼體完全碎裂。由于試驗在野外溝壑中進行,響應(yīng)后僅回收到少部分碎片,如圖5所示。根據(jù)響聲及以往同種炸藥的烤燃試驗情況,綜合判斷烤燃彈的響應(yīng)等級為爆轟。

圖5 試驗結(jié)果Fig.5 Test results

4 算例及驗證

4.1 檢驗理論與試驗的一致性

為檢驗理論與試驗的一致性,以試驗藥柱為例,利用表2中的物性參數(shù)和動力學(xué)參數(shù),對式(21)的解析解進行計算。

表2 RDX炸藥物性參數(shù)和動力學(xué)參數(shù)

對于任意長圓柱,單從徑向或軸向考慮,傳熱均屬于一維問題,僅考慮徑向傳熱時,因溫度沿縱軸均勻分布,沒有區(qū)別。而單考慮軸向傳熱時,因溫度分布不均勻,能體現(xiàn)出位置的差別。同理,僅考慮自熱兩個方向也均屬于一維問題,故對藥柱軸向的溫度分布進行計算[19]。

已知炸藥厚度 、初始溫度為305 K、升溫速率κ=3.3 ℃/h。使用Maple數(shù)學(xué)軟件[20]和Clenshaw-Curtis正交方法[21]求解,控制無窮級數(shù)的舍入誤差不大于0.1%。理論解的藥柱沿軸線方向的部分溫度值如表3所示。

表3 沿軸線方向的溫度值

外壁溫度為484.6 K時,理論解與測點測量值的溫度分布如圖6所示。各測點測量得到的溫度分布趨勢與理論解的溫度分布趨勢相符,但是由于試驗值是三維結(jié)果,而理論僅是一維結(jié)果,因此兩曲線的彎曲程度有所不同(即對于試驗藥柱,沿軸向的傳熱屬于一維問題[19])。

圖6 溫度分布Fig.6 Temperature distribution

為直觀地比較理論解的相對溫度分布與試驗測得的相對溫度分布,將圖6中試驗測得的溫度值、理論解中截取與試驗測量范圍(-0.095～0.095 m)對應(yīng)區(qū)間的溫度值,再根據(jù)各自的溫度幅值分別作如下歸一化處理:

(22)

式中:xN為歸一化后的溫度值;xi為歸一化前的溫度值;xmax為歸一化前數(shù)據(jù)中的最高溫度值;xmin為歸一化前數(shù)據(jù)中的最低溫度值。通過歸一化將試驗測得的溫度值、理論解與試驗對應(yīng)位置的溫度值映射到0～1之間,得到各個位置的相對溫度分布,如圖7所示。理論解各個位置的溫度之間的比例與試驗結(jié)果一致,即理論解的相對溫度分布與試驗測得的相對溫度分布相符,進一步驗證了理論解的溫度分布的正確性。

圖7 相對溫度分布Fig.7 Relative temperature distribution

4.2 數(shù)值驗證

為進一步驗證理論解的正確性,使用數(shù)值方法對解析解式(21)進行驗證。Suceka[22]對熱傳導(dǎo)方程式(1)進行了有限差分近似。對于無限大平板而言,關(guān)于位置x的偏導(dǎo)數(shù)可以通過以下差分來近似:

(23)

式中:Δx為將厚度為L的炸藥被均勻分為n份后的寬度;上下標(biāo)用于區(qū)分網(wǎng)格對應(yīng)的時間和位置,上標(biāo)j為時間序列,下標(biāo)i為位置序列。

對時間t的偏導(dǎo)數(shù)進行同樣的有限差分近似:

(24)

式中:Δt為時間增量。

使用有限差分替代熱傳導(dǎo)方程式(1)中的偏導(dǎo)數(shù),即得到炸藥熱傳導(dǎo)離散方程:

(25)

同理,建立非反應(yīng)性熱傳導(dǎo)離散方程如下:

(26)

對于式(25)與式(26)具有相同定解的條件如下:

通過式(25)可以獲得任意時刻炸藥內(nèi)部的溫度分布;將式(25)與式(26)同時刻的溫度結(jié)果作差,即可獲得該時刻炸藥自熱反應(yīng)熱傳導(dǎo)的溫度分布。

以上兩種顯式方法的穩(wěn)定性由空間步長和時間步長的比值決定。對于熱傳導(dǎo)方程,如果滿足以下穩(wěn)定條件[23],則該方法是穩(wěn)定的:

(27)

穩(wěn)定性條件保證了有限差分解的穩(wěn)定性,該條件僅在炸藥分解緩慢的情況下,能保證解的精確性。對于慢速烤燃,點火時自熱反應(yīng)產(chǎn)生的熱量很低,故適用于該判據(jù)。

計算選取RDX炸藥為研究對象,計算中使用的物性參數(shù)和動力學(xué)參數(shù)的值如表2所示,炸藥厚度、初始溫度Tr為305 K、升溫速率κ=3.3 ℃/h。代入以上參數(shù),對熱傳導(dǎo)離散方程式(25)、熱傳導(dǎo)離散方程式(26)賦予初始條件與邊界條件,以時間增量Δt=1 s循環(huán)迭代離散方程(即j+1時刻的溫度分布由j時刻的溫度分布計算得出),并將各個時刻的溫度分布進行輸出。

點火點的溫度變化如圖8所示。在D點(加熱時間約為177 000 s時)發(fā)生點火,自熱反應(yīng)逐漸加劇,溫度不再可控,進而發(fā)生響應(yīng)。

圖8 點火點溫度的理論解與數(shù)值解Fig.8 Theoretical and numerical results of the temperature of ignition point

炸藥點火溫度的理論計算結(jié)果與采用文獻[22]方法的數(shù)值解具有良好的匹配度,理論計算值為468.8 K;數(shù)值計算值為468.6 K,誤差不超過1%。

點火時刻(外壁溫度為480.7 K時),炸藥總體溫度(即非反應(yīng)性熱傳導(dǎo)溫度與自熱反應(yīng)熱傳導(dǎo)溫度之和)分布與僅自熱反應(yīng)所產(chǎn)生的溫度分布如 圖9 所示,點火點位置位于炸藥兩側(cè),與試驗結(jié)果相符?？傮w溫度分布主要受升溫邊界的影響,由于點火前達到分解能壘的分子數(shù)很少,自熱反應(yīng)所產(chǎn)生的溫度較低。自熱反應(yīng)所產(chǎn)生的溫度的理論計算結(jié)果與數(shù)值解相比,最大誤差為8.8%(誤差主要來自于式(14)對源項的近似處理);自熱反應(yīng)溫度最高值所在位置的計算結(jié)果與數(shù)值解相比,最大誤差不超過1%。

圖9 炸藥厚度L=0.285 m時的理論計算結(jié)果與數(shù)值計算結(jié)果Fig.9 Theoretical and numerical results with explosive thickness L=0.285 m

5 結(jié)果與討論

自熱反應(yīng)溫度最高值的位置是由炸藥自熱反應(yīng)與邊界的熱散失共同作用的結(jié)果。圖10為不同厚度條件下,炸藥自熱反應(yīng)溫度最高值隨時間的變化情況。由圖10可見,隨著時間的增加,自熱反應(yīng)溫度最高值所在的相對位置(即|x/δ|)將趨于穩(wěn)定,直至發(fā)生點火。

圖10 不同厚度炸藥自熱反應(yīng)溫度最高值所在的相對位置Fig.10 Relative location of the maximum value point of self-heating reaction temperature of explosiveswith different thicknesses

在炸藥厚度較薄的情況下,自熱反應(yīng)溫度最高位置將迅速移動至中心;隨著炸藥厚度L的增加,位置移動的相對量逐漸減小,當(dāng)炸藥厚度L>0.3 m時,從50 000 s直至發(fā)生點火,自熱反應(yīng)溫度最高位置的相對移動量不大于2%。顯然可以認(rèn)為在烤燃試驗進行到一定時間之后,溫度相對分布趨于穩(wěn)定,自熱反應(yīng)溫度最高值位置隨時間不再發(fā)生改變,直至點火。這樣,可以通過理論計算,預(yù)則最高溫度及其所在位置,為降低響應(yīng)等級的結(jié)構(gòu)設(shè)計提供便利。

點火時,自熱反應(yīng)溫度最高值所在位置至邊界的距離變化如圖11所示。當(dāng)炸藥厚度L>0.3 m時,隨著炸藥厚度的進一步增加,自熱反應(yīng)溫度最高值所在位置至邊界的距離變化不大于2%,且最高值所在位置至邊界的距離逐漸趨于穩(wěn)定。即隨著炸藥厚度的增加,邊界的熱散失條件逐漸成為影響點火點位置的重要因素。

圖11 自熱反應(yīng)溫度最高值所在位置至邊界的距離Fig.11 Distance from the maximum value point of self-heating reaction temperature to boundary

當(dāng)炸藥厚度L>0.3 m時,從邊界至中心的自熱反應(yīng)溫度梯度如圖12所示,隨著炸藥厚度的增加,點火點附近的自熱反應(yīng)溫度梯度趨于恒定。隨著炸藥厚度的增加,自熱反應(yīng)溫度分布的規(guī)律相似。

圖12 炸藥邊界到中心的自熱反應(yīng)溫度梯度Fig.12 Temperature gradient of self-heating reaction from boundary to center of explosives

6 結(jié)論

本文通過理論解析的方法研究一維凝聚炸藥慢速烤燃的點火點位置及點火溫度。得出以下主要結(jié)論:

1)將炸藥非反應(yīng)性熱傳導(dǎo)與反應(yīng)性熱傳導(dǎo)拆分,通過理論方程及計算可得出慢烤條件下一維凝聚炸藥內(nèi)部的溫度分布情況及點火點位置,并與試驗及數(shù)值計算結(jié)果進行對比,驗證了理論方程及計算的正確性。

2)通過理論計算可得到一維炸藥點火時刻僅自熱反應(yīng)所產(chǎn)生的溫度分布和自熱反應(yīng)溫度最高值所在位置;獲得點火時刻一維凝聚炸藥溫度最高值所在位置及溫度梯度沿厚度的變化規(guī)律。為烤燃試驗點火溫度及位置的預(yù)測提供了可靠的理論依據(jù)。