NEPE推進劑熱安全性的理論計算與數(shù)值模擬

秦沛文, 王夢蕾,4, 魯銳華, 郭春亮, 黃昌龍, 李 穎, 趙孝彬,3, 彭 松,3

(1.湖北航天化學(xué)技術(shù)研究所, 湖北 襄陽 441003；2. 航天化學(xué)動力技術(shù)重點實驗室, 湖北 襄陽 441003；3. 航天工業(yè)固體推進劑安全技術(shù)研究中心, 湖北 襄陽 441003；4.應(yīng)急救生與安全防護湖北省重點實驗室，湖北 襄陽 441003)

引 言

表征含能材料熱安全性的參數(shù)包括熱爆炸臨界環(huán)境溫度、熱爆炸延滯時間、安全度、熱爆炸概率等[1]。傳統(tǒng)的用于測定熱安全性參數(shù)試驗，如爆發(fā)點、火焰感度測試等，試驗量小、便于實驗室開展。但一些含能材料由于熱作用發(fā)生燃燒或爆炸時，尺寸的影響顯著[2-4]，用經(jīng)典的試驗方法獲得的數(shù)據(jù)難以推斷含能材料在實際狀態(tài)下的潛在危險性[5]，這種現(xiàn)象被國外學(xué)者形象地稱為熱安全性的尺寸效應(yīng)[6]。不同種類含能材料熱安全性尺寸效應(yīng)的顯著程度不同。對于含有大量起始分解溫度較低的硝酸酯組分作為增塑劑的NEPE推進劑，對于同一環(huán)境溫度，尺寸的增大可令熱爆炸延滯期由幾百小時縮短至幾十小時，并且還有繼續(xù)減小的趨勢[7]。采用全尺寸試驗的方法評估大尺寸NEPE推進劑的熱安全性耗費高、危險性大，開展難度極大，至今未見報道。因此，國內(nèi)外學(xué)者采用理論與數(shù)值模擬的方法對含能材料及其裝藥的熱安全性做了不同程度的研究與探索。

Roduit等[8]認為含能材料的熱安全性受藥量的影響很大，并認為用準確的動力學(xué)方程描述含能材料的熱分解過程是評價含能材料熱作用結(jié)果的關(guān)鍵。Krause[6]進行了HNS、NTO、HMX、FOX-12和雙基藥共5種不同尺寸含能材料的熱爆炸試驗，并對通過該試驗獲得的動力學(xué)參數(shù)與文獻中相應(yīng)的參數(shù)作了對比，認為熱爆炸臨界尺寸是表觀活化能、阿倫尼烏斯升溫速率(與升溫速率相同量綱，定義為AQ/Cp)、導(dǎo)熱系數(shù)等參數(shù)的函數(shù)，含能材料的熱安全性可能由以上參數(shù)決定，一旦確定了這些參數(shù)，不同藥量、臨界條件下的熱爆炸臨界環(huán)境溫度便可確定。

路桂娥等[9]進行了藥量范圍5～30kg的9/7單基發(fā)射藥、雙芳-3 16/1雙基發(fā)射藥和273火箭彈推進劑雙芳鎂發(fā)射藥熱自燃實驗，并與小型實驗結(jié)果進行了比較。通過外推求得了貯存40年后不同含能材料的自燃溫度，發(fā)現(xiàn)結(jié)果與小型實驗結(jié)果差別較大。認為主要原因是大尺寸裝藥不容易散熱，容易產(chǎn)生熱積累。

本研究在前期試驗[7]的基礎(chǔ)上，將相關(guān)結(jié)論工程化地應(yīng)用于不同尺寸NEPE推進劑熱安全性的理論計算與數(shù)值模擬，計算了不同尺寸NEPE推進劑的熱安全性，可為填裝NEPE推進劑的固體火箭發(fā)動機的熱安全性評估工作提供一種可參考的計算方法。

1 理論計算

傳統(tǒng)的熱爆炸理論計算方法[10]對熱安全性的尺寸效應(yīng)問題描述不夠精細，不能很好地反映尺寸對含能材料熱安全性參數(shù)的影響，導(dǎo)致大尺寸含能材料熱安全性的理論計算結(jié)果偏差較大[6]。本研究用Semenov理論、Frank-Kamenetskii理論和Thomas理論計算了不同尺寸NEPE推進劑的熱爆炸臨界環(huán)境溫度，并應(yīng)用前期研究中建立的NEPE推進劑尺寸與熱分解反應(yīng)表觀動力學(xué)參數(shù)的數(shù)學(xué)關(guān)系[7]對以上理論算法作了修改，對修改前后得到的計算結(jié)果作了對比，從中優(yōu)選出更適用于NEPE推進劑的算法。在此基礎(chǔ)上，通過進一步計算對比了尺寸與溫度對NEPE推進劑熱安全性的影響。

1.1 熱爆炸理論及其工程化修改

(1)Semenov理論算式及其修改算式為：

(2)Frank-Kamenetskii理論算式及其修改算式為：

(3)Thomas理論算式及其修改算式為：

式(1)～(3)中：V為NEPE推進劑體積，m3；ρ為NEPE推進劑密度，kg/m3,取實測值ρ=1846kg/m3；Q為反應(yīng)熱，J/kg,考慮到NEPE推進劑在熱爆炸臨界溫度附近的放熱主要來源于硝酸酯的熱分解，取混合硝酸酯反應(yīng)熱Q=6641.5×103J/kg(該NEPE推進劑實測真空定容爆熱為7140.5×103J/kg)；Ea為表觀活化能，J/mol；A為指前因子，s-1；R為普適氣體常數(shù)，取R=8.314J/(mol·K)；Tcr為熱爆炸臨界溫度，K；S為NEPE推進劑散熱面積，m2；δcr為Frank-Kamenetskii準數(shù)，對于長徑比為1的藥柱，取δcr=2.77[11]；對于長徑比為1的藥柱，r為藥柱半徑，m；λ為NEPE推進劑導(dǎo)熱系數(shù)，取實測值λ=0.422W/(m·K)；κ為空氣對流換熱系數(shù)，取κ=10W/(m-2·K)；直接采用(工程化修改前)計算時，Ea與A取熱分析法獲得的NEPE推進劑的表觀動力學(xué)參數(shù)Ea=147.1kJ/mol，A=1.96×1013s-1。

采用以上修改前后的熱爆炸理論算法計算直徑不同、長徑比為1的NEPE推進劑藥柱的熱爆炸臨界環(huán)境溫度如圖1所示，其中，曲線起點的橫坐標為10mm，終點的橫坐標為3000mm，相鄰描點間距ΔD=10mm。

圖1 不同直徑時長徑比為1的NEPE推進劑藥柱熱爆炸臨界環(huán)境溫度曲線Fig.1 Curves of thermal explosion critical ambient temperature of NEPE propellant grains with different diameters and length-diameter ratio of 1

由圖1所示的各理論計算結(jié)果與文獻[7]中所述NEPE推進劑尺寸的實測結(jié)果對比可知，修改前的Semenov理論、Frank-Kamenetskii理論和Thomas理論計算得到熱爆炸臨界環(huán)境溫度均偏高，修改算法中，工程化修改的Semenov算法的計算結(jié)果與試驗結(jié)果更為接近。這可能源于以下2個原因：

(1)假設(shè)與實際相符。對于放置于溫度與熱爆炸臨界環(huán)境溫度相差不大環(huán)境中的NEPE推進劑藥柱體系，通過藥柱內(nèi)部不同位置溫度的監(jiān)測[7]，相對于Frank-Kamenetskii和Thomas理論，Semenov理論的均溫假設(shè)與該狀態(tài)下NEPE藥柱內(nèi)部實際的溫度分布吻合。

(2)改進空間更大。相較于其余理論，Semenov理論更為原始、基礎(chǔ)，并未體現(xiàn)對尺寸效應(yīng)的過多考慮，因而改進空間更大。表觀動力學(xué)參數(shù)是影響熱安全性計算結(jié)果的最敏感變量，傳統(tǒng)的熱安全性計算采用基于熱分析試驗獲得的該參數(shù)作為參照，且該參數(shù)一般基于其熱分解峰值溫度求得。但對于NEPE推進劑，此種方法存在諸多不足。首先，實際尺寸的NEPE推進劑的熱爆炸臨界溫度與熱分析測得的峰值溫度差距較大，峰值溫度附近主要為HMX和部分AP的熱分解[11-12]，而熱爆炸多發(fā)生在硝酸酯遷移、分解階段對應(yīng)的溫度區(qū)間內(nèi)，此時推進劑內(nèi)部主要發(fā)生硝酸酯的遷移、分解，幾乎沒有HMX和AP的分解放熱；其次，熱分析試樣量為毫克級，大部分硝酸酯還來不及熱分解便迅速遷移至周圍氮氣氛圍中，進而被氣流帶走，因此不能反映硝酸酯遷移對動力學(xué)參數(shù)的影響。

結(jié)合試驗中得出的NEPE推進劑尺寸與熱分解反應(yīng)表觀動力學(xué)參數(shù)的數(shù)學(xué)關(guān)系對Semenov理論進行修改，能夠?qū)⒗碚撚嬎闩c實際試驗結(jié)果結(jié)合起來，使計算結(jié)果更接近試驗結(jié)果。該曲線還反映了隨著尺寸增大，熱爆炸臨界環(huán)境溫度由130℃逐漸向60℃靠近，且尺寸越小，增大尺寸導(dǎo)致熱爆炸臨界環(huán)境溫度降低的效果越明顯。

1.2 尺寸與溫度對推進劑放熱速率的影響

考慮到NEPE推進劑在熱爆炸臨界溫度附近的放熱主要來源于硝酸酯的熱分解[13-14]及NEPE推進劑尺寸與熱分解反應(yīng)表觀動力學(xué)參數(shù)的數(shù)學(xué)關(guān)系[7]，將單位體積NEPE推進劑的放熱速率用公式(4)表示：

(4)

式中：Qv為藥柱整體放熱速率;α為反應(yīng)深度(已反應(yīng)的百分數(shù));n為反應(yīng)級數(shù);ΔrH為NEPE推進劑在恒定環(huán)境溫度下熱分解反應(yīng)焓變，NE代表混合硝酸酯。

根據(jù)式(4)求出不同尺寸NEPE推進劑藥柱在不同溫度環(huán)境中單位體積的放熱量，見圖2，曲線中相鄰描點的實際間距為ΔD=10mm，ΔT=1℃。

圖2 不同尺寸時長徑比為1的NEPE推進劑在不同溫度下的單位體積放熱量Fig.2 Heat release per unit volume of NEPE propellant grains with different sizes and length-diameter ratio of 1 at different temperature ranges

由圖2可知，對于NEPE推進劑，溫度對放熱速率的影響比尺寸對放熱速率的影響更為顯著，NEPE推進劑的熱安全性對溫度的變化更為敏感。

2 數(shù)值模擬

采用有限元法進行NEPE推進劑的熱安全性分析，可以通過調(diào)整邊界條件、熱載荷并采用試探法確定NEPE推進劑的熱爆炸臨界環(huán)境溫度，與理論計算相比，效率較高，假設(shè)條件較少，更適用于結(jié)構(gòu)復(fù)雜推進劑的熱安全性評估。除此之外，還可通過采取不同的建模方式，方便地模擬不同熱載荷條件下NEPE推進劑內(nèi)部的溫度分布，確定裝藥的熱爆炸臨界環(huán)境溫度、環(huán)境溫度變化時裝藥內(nèi)部溫度響應(yīng)及高溫環(huán)境下裝藥的熱爆炸延滯時間等問題[16]。

本研究將NEPE推進劑尺寸與熱分解反應(yīng)表觀動力學(xué)參數(shù)的數(shù)學(xué)關(guān)系[7]應(yīng)用于數(shù)值模擬，采用有限元方法計算了NEPE推進劑的熱安全性。

2.1 數(shù)值模擬的模型與方法

NEPE推進劑在恒溫環(huán)境中的貯存過程可以看作是一個具有內(nèi)部化學(xué)反應(yīng)放熱、邊界對流換熱的瞬態(tài)熱傳導(dǎo)過程，假設(shè)推進劑內(nèi)部各處放熱速率相同且表面各處熱交換均勻一致，根據(jù)傅里葉熱傳導(dǎo)定律，NEPE推進劑的瞬態(tài)導(dǎo)熱方程見式(5)：

(5)

式中：C為比熱容，取實測值C=1733J/(kg·K)；qv為單位體積NEPE推進劑的放熱速率，作如下處理：

(1)對于未作工程化修改的算法，qv可采用式(6)替代：

(2)對于已作工程化修改的算法，qv采用式(4)替代。

當(dāng)環(huán)境溫度為熱爆炸臨界環(huán)境溫度時，推進劑外表面處換熱達到熱平衡，根據(jù)牛頓冷卻定律(忽略空氣熱傳導(dǎo)與熱輻射)，邊界條件見式(7)：

(7)

式中：Tw和Tf分別為固體推進劑外表面溫度和環(huán)境溫度。

采用有限元法對不同裝藥尺寸的NEPE推進劑的熱安全性進行數(shù)值模擬。首先，根據(jù)推進劑裝藥尺寸建立1∶1的計算模型。然后，根據(jù)式(6)或式(7)編寫用戶自定義函數(shù)(qv)。由于理論上熱爆炸臨界環(huán)境溫度處，熱爆炸延滯期為無限長。因此采用二分法逐步逼近原理，由高到低逐步調(diào)整環(huán)境溫度，并在各預(yù)設(shè)環(huán)境溫度下，逐步延長貯存時間，觀察裝藥內(nèi)部的溫度分布。并以此為依據(jù)，逐步縮小熱爆炸臨界環(huán)境溫度范圍，從而最終確定熱爆炸臨界環(huán)境溫度范圍，同時可粗略獲得熱爆炸延滯期范圍。

以裝藥1/4模型截面設(shè)置溫度對稱邊界條件，初始溫度為15℃，并采用8節(jié)點三維熱實體單元(SOLID70)進行網(wǎng)格劃分，所涉及數(shù)值模擬單元數(shù)目均不少于6×104，物性參數(shù)與本研究理論計算取值相同。對于瞬態(tài)熱分析，載荷子步越多，結(jié)果反而越差。因此，當(dāng)所模擬的貯存時間較長時，為減少載荷步數(shù)，提高計算精度和速度，將單位體積放熱速率(qv)、空氣對流換熱系數(shù)(α)等所有與時間有關(guān)的輸入?yún)?shù)，均轉(zhuǎn)化為采用天(d)計量的單位，如λ=0.422W/(m·K)=0.422×3600×24J/(d·m·K)，其余物理量均采用國際單位制。載荷步長根據(jù)需要模擬的貯存時間長短靈活調(diào)整，兼顧計算效率與精度，如當(dāng)貯存時間不超過1000d時，載荷步長取1d，超過1000d時，載荷步長根據(jù)貯存時間可取10d和100d等。

2.2 某NEPE推進劑裝藥熱安全性的數(shù)值模擬

對圖3所示NEPE推進劑裝藥的熱安全性進行數(shù)值模擬。

圖3 用于計算的NEPE推進劑裝藥尺寸簡化圖Fig.3 Simplified diagram of NEPE propellant grain size for calculation

該NEPE推進劑藥柱外包覆EPDM絕熱層和鋼殼體，絕熱層厚度為2mm，殼體厚度為3mm。由于絕熱層和鋼殼體厚度較薄，自身不產(chǎn)熱，相對于推進劑本體，其熱阻在熱安全性計算中可忽略不計，因此通過計算NEPE推進劑本體的安全性代表整體的安全性。

該NEPE推進劑藥柱在臨近熱爆炸臨界環(huán)境溫度的環(huán)境中貯存不同時間后的溫度分布如圖4所示。

圖4 基于數(shù)值模擬獲得的某NEPE推進劑藥柱在臨近熱爆炸臨界溫度環(huán)境中貯存不同時間后的內(nèi)部溫度分布Fig.4 The interior temperature distribution of one NEPE propellant grain stored with different time near thermal explosion critical ambient temperature based on numerical simulation

數(shù)值模擬結(jié)果表明，該裝藥在78℃環(huán)境中分別貯存200d和400d后內(nèi)部溫度基本無變化，可認為已經(jīng)達到熱平衡；在79℃環(huán)境中貯存150d后顯然已發(fā)生熱失控，因此，根據(jù)數(shù)值模擬結(jié)果推測該NEPE推進劑裝藥的熱爆炸臨界環(huán)境溫度可能為78～79℃，在79℃環(huán)境中熱爆炸延滯期可能為100～150d。

為驗證該數(shù)值模擬結(jié)果的合理性，將圖3所示裝藥放入60℃烘房中貯存1752h(73d)，現(xiàn)場狀態(tài)如圖5所示。

圖5 某NEPE推進劑裝藥60℃貯存試驗現(xiàn)場照片F(xiàn)ig.5 Test photo of one NEPE grains tored at 60℃

然后，采用工業(yè)CT對對該裝藥進行X射線無損檢測分析，結(jié)果見圖6。探傷完成后對該裝藥進行“挖藥”實驗，探傷檢測和現(xiàn)場“挖藥”結(jié)果表明，該裝藥未出現(xiàn)脫粘，內(nèi)部未見氣孔、裂紋，裝藥完整性未受到破壞。結(jié)合Ф100mm×100mm裝藥90℃貯存試驗結(jié)果[7]，可以推測該裝藥的熱爆炸臨界環(huán)境溫度介于60～90℃，能在一定程度上反映以上數(shù)值模擬結(jié)果的合理性，但并不充分。

圖6 某NEPE推進劑裝藥貯存試驗后X射線照片F(xiàn)ig.6 X-ray photo of one NEPE propellant grain after storage test

3 不同尺寸NEPE推進劑熱安全性的理論計算與數(shù)值模擬

采用上述計算方法對不同尺寸NEPE推進劑裝藥的熱安全性進行理論計算與數(shù)值模擬，結(jié)果見表1。

由表1可知，對于同一尺寸NEPE推進劑，原始的數(shù)值模擬和Semenov理論計算結(jié)果基本吻合，說明二者采用的基本原理并無差別。對于工程化修改的數(shù)值模擬和Semenov理論計算結(jié)果，二者有一定的差距，并且對于較小尺寸NEPE推進劑，差距較大；對于較大尺寸NEPE推進劑，差距較小。

表1 不同尺寸NEPE推進劑熱爆炸臨界環(huán)境溫度的數(shù)值模擬與理論計算結(jié)果

注：Tn為數(shù)值模擬結(jié)果；Tt為理論計算結(jié)果。

由文獻[7]中的試驗結(jié)果，Ф100mm×100mm和Ф150mm×150mm NEPE推進劑藥柱在90℃環(huán)境中發(fā)生了熱爆炸，可知原始的熱爆炸理論算法和有限元算法得出的NEPE推進劑的熱爆炸臨界溫度高于試驗結(jié)果。本研究所采取的工程化修改措施對NEPE推進劑熱安全性的理論計算和數(shù)值模擬結(jié)果均起到了一定程度的修正作用。

4 結(jié) 論

(1)對于NEPE推進劑，原始的Semenov、Frank-Kamenetskii、Thomas理論算法和有限元算法得出的熱爆炸臨界環(huán)境溫度較實際偏高。其中，工程化修改的Semenov算法得出的計算結(jié)果與試驗結(jié)果更為接近。長徑比為1的藥柱在直徑由10mm增至3000mm過程中，熱爆炸臨界環(huán)境溫度由130℃逐漸向60℃靠近，且尺寸越小，增大尺寸導(dǎo)致熱爆炸臨界環(huán)境溫度降低的效果越明顯。

(2)通過理論計算或數(shù)值模擬評估裝藥熱安全性，若計算結(jié)果與實際差別大，可先從熱爆炸試驗獲得熱分解反應(yīng)表觀動力學(xué)參數(shù)與尺寸之間的數(shù)學(xué)關(guān)系、再將該關(guān)系應(yīng)用于熱爆炸理論計算或數(shù)值模擬仿真修正誤差。對于尺寸效應(yīng)明顯的含能材料，應(yīng)增大試驗樣本尺寸變化范圍并增加樣本數(shù)。