劉海，何遠(yuǎn)航

（北京理工大學(xué)爆炸科學(xué)與技術(shù)國家重點(diǎn)實(shí)驗(yàn)室，北京100081）

梯恩梯/萘共晶初始熱分解的反應(yīng)分子動(dòng)力學(xué)模擬

劉海，何遠(yuǎn)航

（北京理工大學(xué)爆炸科學(xué)與技術(shù)國家重點(diǎn)實(shí)驗(yàn)室，北京100081）

ReaxFF反應(yīng)力場在沖擊起爆、爆轟等問題中的應(yīng)用多圍繞常規(guī)含能材料。利用增加了長程修正項(xiàng)Elg的ReaxFF/lg反應(yīng)力場對(duì)梯恩梯/萘（TNT/C10H8）共晶初始高溫?zé)岱纸膺M(jìn)行模擬，并通過TNT單晶比較了C10H8對(duì)共晶整體反應(yīng)特征的影響，利用指數(shù)函數(shù)及反應(yīng)速率方程擬合得到共晶初級(jí)吸熱反應(yīng)和次級(jí)放熱反應(yīng)的活化能分別為35.7 kcal/mol和56.1 kcal/mol.初級(jí)吸熱反應(yīng)的活化能與TNT單晶基本相同，而次級(jí)放熱反應(yīng)的活化能則遠(yuǎn)遠(yuǎn)高于TNT單晶，并且同溫度條件下共晶次級(jí)反應(yīng)放熱量小于TNT單晶。擬合得到的反應(yīng)物的衰減速率表明，C10H8的加入將抑制共晶內(nèi)TNT的分解。產(chǎn)物識(shí)別分析顯示初始產(chǎn)物為NO2、NO和HONO，并且通過NO2和TNT—NO2，NO和TNT—NO，HONO和TNT—HONO分布數(shù)量的比較論證了固相共晶和TNT單晶的初始反應(yīng)路徑為雙分子反應(yīng)機(jī)制。共晶熱分解最終主要產(chǎn)物為N2，H2O，CO2和CO.并且由于共晶內(nèi)C10H8分子C—C鍵斷裂所需的能量高于C—H鍵斷裂所需的能量，C10H8分子的初始熱分解路徑為C—H鍵斷裂，并且形成的H原子將促進(jìn)共晶內(nèi)H2O的產(chǎn)率高于TNT單晶內(nèi)H2O的產(chǎn)率。

兵器科學(xué)與技術(shù)；共晶；高溫?zé)岱纸?；反?yīng)動(dòng)力學(xué)；ReaxFF反應(yīng)力場

0 引言

極端條件下固相含能材料的物理化學(xué)響應(yīng)問題是爆轟學(xué)的核心內(nèi)容之一，并且，近幾十年的實(shí)驗(yàn)測試和理論分析多圍繞常規(guī)含能材料沖擊起爆和爆轟波的宏觀問題，如爆壓、爆速、爆轟波結(jié)構(gòu)以及穩(wěn)定性等［1］。極端溫度/壓力條件下固相含能材料的初始反應(yīng)路徑是研究沖擊起爆問題的基礎(chǔ)，同樣也與含能材料爆轟波能量輸出結(jié)構(gòu)直接關(guān)聯(lián)。由于固相含能材料的沖擊起爆及爆轟問題涉及極短的時(shí)間和極小的空間尺度，這給現(xiàn)有的實(shí)驗(yàn)技術(shù)帶來了很大的挑戰(zhàn)。近期，各國學(xué)者利用美國加州理工學(xué)院Goddard小組提出的ReaxFF反應(yīng)力場分子動(dòng)力學(xué)方法［2］已經(jīng)就常規(guī)含能材料高溫?zé)岱纸獾某跏挤磻?yīng)路徑［3-5］以及沖擊起爆機(jī)制［6-16］進(jìn)行了廣泛而深入地研究。

包含氧化劑和燃料組分的含能材料主要通過激發(fā)后釋放的化學(xué)能產(chǎn)生高溫/高壓燃燒產(chǎn)物以達(dá)到毀傷、推進(jìn)等效果。但限于沖擊感度、密度、熔點(diǎn)、穩(wěn)定性等因素，當(dāng)前僅部分含能材料可用于戰(zhàn)斗部和推進(jìn)劑。共晶，作為一種改性技術(shù)，主要通過分子間非共價(jià)鍵將不同種類的含能材料分子結(jié)合在同一晶格中以實(shí)現(xiàn)含能材料的高能低感特性［17］。國內(nèi)外在共晶含能材料領(lǐng)域的研究尚處于初步的實(shí)驗(yàn)室合成階段，但已經(jīng)取得了一定的成果（六硝基六氮雜異伍茲烷/梯恩梯（CL20/ TNT）［18-19］、六硝基六氮雜異伍茲烷/奧克托今（CL20/HMX）［20］、六硝基六氮雜異伍茲烷/三氟甲基苯（CL20/BTF）［21］、奧克托今/三氨基三硝基苯（HMX/TATB）［22］、梯恩梯/萘（TNT/C 10H8）［23］），實(shí)驗(yàn)結(jié)果業(yè)已表明上述共晶含能材料在爆速及安定性方面具有優(yōu)異表現(xiàn)。當(dāng)前，分子動(dòng)力學(xué)技術(shù)在熱分解以及沖擊起爆等問題中的應(yīng)用多面向常規(guī)含能材料，針對(duì)共晶含能材料的相關(guān)研究卻鮮有發(fā)表［24-26］。本文利用ReaxFF反應(yīng)力場分子動(dòng)力學(xué)方法對(duì)TNT/C10H8共晶的初始熱分解問題進(jìn)行了研究，此結(jié)果可為優(yōu)化共晶含能材料內(nèi)分子配比、合成高能低感新型共晶含能材料等問題提供參考。

1 模型及計(jì)算細(xì)節(jié)

TNT/C10H8共晶初始單胞結(jié)構(gòu)來自X射線衍射測試結(jié)果［23］。初始單胞內(nèi)含有2個(gè) TNT和2個(gè)C10H8分子。隨后構(gòu)建4×2×2共晶超胞（見圖1），晶胞內(nèi)TNT和C10H8的分子數(shù)量均為32個(gè)，共計(jì)1 248個(gè)原子。共晶初始密度為1.495 g/cm3，其中C10H8的密度為0.539 g/cm3.首先利用共軛梯度法通過系統(tǒng)勢能收斂對(duì)超胞進(jìn)行能量最小化。隨后利用Berendsen熱浴方法對(duì)超胞進(jìn)行溫度和壓力弛豫，以獲得常溫（300 K）下的晶體結(jié)構(gòu)，此過程持續(xù)5 ps.為了比較共晶內(nèi)C10H8對(duì)TNT熱分解的影響，相應(yīng)地構(gòu)建2×2×2正交TNT超胞結(jié)構(gòu)（見圖1），TNT初始單胞結(jié)構(gòu)來自劍橋晶體結(jié)構(gòu)數(shù)據(jù)庫（CCDC）［27］，初始密度為1.704 g/cm3，初始單胞內(nèi)含有8個(gè)TNT分子。為了消除密度差異對(duì)熱分解結(jié)果造成的影響，在構(gòu)建TNT超胞前，對(duì)TNT單胞進(jìn)行體積放大以獲得近似于TNT/C10H8共晶的密度（1.510 g/cm3），超胞內(nèi)共計(jì)64個(gè)TNT，1 344個(gè)原子。能量最小化以及溫度和壓力弛豫過程同上。最后對(duì)兩含能材料晶體施加2 000 K～4 500 K高溫以模擬含能材料爆轟Chapmam-Jouguet溫度，此過程持續(xù)25 ps，整個(gè)模擬過程的時(shí)間步長為0.1 fs，高溫?zé)岱纸膺^程采用集成于lammps程序［28-29］內(nèi)的ReaxFF/lg［30］分子動(dòng)力學(xué)進(jìn)行模擬。

ReaxFF/lg是在ReaxFF勢函數(shù)的基礎(chǔ)上增加了長程修正項(xiàng)Elg.ReaxFF/lg勢函數(shù)在含能材料晶體結(jié)構(gòu)以及狀態(tài)方程的預(yù)測相比于ReaxFF勢函數(shù)更加接近實(shí)驗(yàn)值［30］。當(dāng)前，ReaxFF/lg勢函數(shù)已經(jīng)在含能材料感度分析［31］、熱分解［24］、熱沖擊［32］、沖擊壓縮［33］、沖擊起爆［6-7］等問題得到應(yīng)用。

圖1　TNT/C10H8共晶和TNT超胞模型Fig.1 Models of TNT/C10H8cocrystal and TNT supercell

2 結(jié)果及分析

2.1勢能演化及曲線擬合

含能材料高溫?zé)岱纸膺^程中，系統(tǒng)將首先吸收熱量使得系統(tǒng)勢能（PE）升高，至峰值后，系統(tǒng)勢能將隨著放熱化學(xué)反應(yīng)的進(jìn)行而開始衰減［3-5］。圖2 為TNT/C10H8共晶和TNT高溫?zé)岱纸膺^程中系統(tǒng)勢能隨時(shí)間的演化曲線。當(dāng)T=2 000 K時(shí)，TNT/ C10H8共晶和TNT分子吸收熱量后，在有限的時(shí)間內(nèi)并沒有發(fā)生放熱化學(xué)反應(yīng)。當(dāng)外界溫度升至3 000 K以上，兩種含能材料在經(jīng)過短暫的“平衡-誘導(dǎo)期”后，系統(tǒng)勢能開始衰減，且隨著溫度的升高，系統(tǒng)勢能開始衰減的時(shí)刻逐漸提前，并且衰減速率逐漸增加，這表明溫度的增加將促進(jìn)放熱化學(xué)反應(yīng)的進(jìn)行。比較同溫度條件下兩種材料勢能衰減曲線，TNT勢能衰減量大于TNT/C10H8共晶，這表明同種溫度條件下，TNT次級(jí)化學(xué)反應(yīng)的放熱量大于TNT/C10H8共晶的放熱量。

含能材料高溫?zé)岱纸膺^程中，勢能的變化和系統(tǒng)內(nèi)化學(xué)反應(yīng)進(jìn)程緊密關(guān)聯(lián)。各溫度條件下，兩種含能材料高溫?zé)岱纸膺^程中系統(tǒng)勢能隨時(shí)間的演化分布表明，系統(tǒng)首先吸收熱量，勢能快速達(dá)到峰值。此階段對(duì)應(yīng)分子由平動(dòng)模態(tài)轉(zhuǎn)換為振動(dòng)模態(tài)，同時(shí)伴隨著初始化學(xué)反應(yīng)［34］，此階段可視為“平衡-誘導(dǎo)期（tE-I）”［3］，且為吸熱反應(yīng)階段。隨后，伴隨著次級(jí)放熱化學(xué)反應(yīng)的進(jìn)行，勢能開始衰減并逐漸達(dá)到平衡，并且系統(tǒng)勢能衰減曲線可用下述指數(shù)函數(shù)［3-5］進(jìn)行擬合：

圖2　TNT/C10H8共晶及TNT熱分解勢能演化曲線Fig.2 Evolution curves of potential energy of thermal decomposition for TNT/C10H8cocrystal and TNT supercell

式中：U（t）為t時(shí)刻系統(tǒng)勢能（kcal/mol）；U∞為勢能漸進(jìn)平衡值（kcal/mol）；ΔUexo為次級(jí)反應(yīng)放熱量（kcal/mol）；tE-I為平衡-誘導(dǎo)期（ps）；τ為特征反應(yīng)時(shí)間（ps，反應(yīng)速率的倒數(shù)）。部分?jǐn)M合曲線如圖3黑線所示，擬合得到的各動(dòng)力學(xué)以及熱力學(xué)參數(shù)如表1.

擬合結(jié)果顯示，隨著溫度的升高，含能材料高溫?zé)岱纸獾摹捌胶?誘導(dǎo)期”和特征反應(yīng)時(shí)間越短，也就是初級(jí)吸熱反應(yīng)，勢能增加至峰值的時(shí)間以及次級(jí)放熱反應(yīng)達(dá)到平衡所需的時(shí)間越少。并且，同溫度條件下，TNT高溫?zé)峤獾姆艧崃烤黠@高于TNT/ C10H8共晶。

系統(tǒng)勢能峰值為高溫?zé)岱纸膺^程中含能材料吸熱和放熱反應(yīng)的分界點(diǎn)。為了比較TNT/C10H8共晶和TNT高溫?zé)岱纸馕鼰岷头艧岱磻?yīng)的活化能，這里采用經(jīng)典 Arrhenius方程，同時(shí)利用表1中 TNT/ C10H8共晶和TNT高溫?zé)峤鈩菽芩p曲線（TNT/C10H8共晶：3000 K、3 500 K、4 000 K、4 500 K；TNT：3 000 K、3 500 K、4 000 K）擬合得到“平衡-誘導(dǎo)期”、特征反應(yīng)時(shí)間進(jìn)行擬合：

式中：Ea為活化能；A為指前因子；R為理想氣體常數(shù)。變換（2）式得到：lnτ=（Ea/RT）-lnA，可見lnτ 與1/T呈線性關(guān)系。擬合結(jié)果為TNT/C10H8共晶和TNT高溫?zé)峤獾奈鼰岷头艧岱磻?yīng)階段的活化能分別為

圖3　兩種含能材料高溫?zé)岱纸膺^程吸熱和放熱反應(yīng)階段的活化能Fig.3 Activation energy of endothermic and exothermic stages of thermal decomposition for cocrystal and TNT single crystal

表1　利用指數(shù)函數(shù)擬合得到的各動(dòng)力學(xué)以及熱力學(xué)參數(shù)Tab.1 Kinetic and thermodynamic parameters by exponential function fitting

吸熱反應(yīng)階段系統(tǒng)勢能快速增加，并且發(fā)生初始化學(xué)反應(yīng)。對(duì)于TNT/C10H8共晶，晶胞內(nèi)TNT分子上—NO2鍵的斷裂，表明已經(jīng)越過共晶熱分解初級(jí)反應(yīng)的能量勢壘，并且此過程需要較小的時(shí)間尺度，因此共晶和TNT熱分解初級(jí)吸熱反應(yīng)的活化能非常接近，且與Furman的模擬結(jié)果［34］基本一致。放熱反應(yīng)階段主要發(fā)生次級(jí)反應(yīng)并逐漸形成最終產(chǎn)物，TNT/C10H8共晶放熱反應(yīng)的活化能為TNT的2倍［34］，此結(jié)果表明C10H8的加入抑制了共晶高溫?zé)岱纸獯渭?jí)放熱反應(yīng)的進(jìn)行。TNT放熱反應(yīng)階段活化能較Fruman模擬結(jié)果［34］高出35.6%.誤差原因可能是由于當(dāng)溫度較高時(shí)，在25 ps的模擬時(shí)間內(nèi)可以較好地描述初級(jí)和次級(jí)反應(yīng)，而當(dāng)溫度相對(duì)較低時(shí)，熱分解平衡則需要較長的時(shí)間，這將導(dǎo)致溫度較低時(shí)的特征反應(yīng)時(shí)間存在誤差。另外，從結(jié)構(gòu)上來說，共晶含能材料多為類三明治結(jié)構(gòu)，相對(duì)惰性的C10H8分子的加入，降低了TNT的濃度，這將會(huì)使得相對(duì)TNT單晶來說，共晶發(fā)生次級(jí)化學(xué)反應(yīng)的活化能增加，這與黑索今（RDX）高溫?zé)岱纸膺^程中發(fā)生次級(jí)化學(xué)反應(yīng)的活化能隨著RDX初始密度的降低而增加的結(jié)論基本一致［3］。

2.2反應(yīng)物衰減速率比較

圖4為TNT/C10H8共晶熱分解過程中，超胞內(nèi)TNT及C10H8分子數(shù)量衰減曲線比較。整體來說，TNT和C10H8分子均隨著溫度的升高，完全分解所需的時(shí)間越少。當(dāng)溫度為2 000 K時(shí)，共晶內(nèi)TNT以及C10H8分子均部分“分解”，對(duì)比初始反應(yīng)路徑（見圖5），共晶以及TNT吸收熱量后可能形成大的分子團(tuán)簇，至模擬后期，TNT分子脫落NO2使得TNT分子的數(shù)量進(jìn)一步減少。當(dāng)溫度為2 500 K、3 000 K 和3 500 K時(shí)，共晶內(nèi)TNT完全分解所需的時(shí)間分別為18.1 ps、10.2 ps和3.8 ps，而C10H8分子完全分解的時(shí)間均大于25 ps，并且當(dāng)溫度為4 000 K和4 500 K時(shí)，共晶內(nèi)TNT分子較C10H8分子均提前完全分解。為了定量比較高溫?zé)岱纸膺^程中共晶內(nèi)TNT和C10H8衰減速率以及C10H8對(duì)共晶內(nèi)TNT衰減速率的影響。這里采用遞減的指數(shù)函數(shù)對(duì)共晶內(nèi)TNT和C10H8以及TNT單晶衰減曲線進(jìn)行擬合［4］：

式中：C0為晶體中TNT和C10H8分子的初始數(shù)量（共晶內(nèi)TNT和C10H8：32，TNT單晶：64）；C（t）為t時(shí)刻TNT或C10H8分子的數(shù)量；k1為衰減速率（ps-1）。通過（3）式分別對(duì)共晶中TNT和C10H8分子以及TNT單晶的分解演化曲線進(jìn)行擬合（如圖4中黑色線所示）得到衰減速率k1如表2所示。

共晶內(nèi)TNT和C10H8分子衰減速率均隨著溫度的升高而加快，并且TNT分子的衰減速率均高于C10H8.比較共晶內(nèi)TNT分子和TNT單晶內(nèi)分子數(shù)量的衰減速率，同溫度條件下，TNT單晶的衰減速率高于共晶內(nèi)TNT，因此，C10H8對(duì)共晶內(nèi)TNT分子的熱分解有抑制作用。

表2　反應(yīng)物衰減速率比較Tab.2 Attenuation rates of reactants

圖4　反應(yīng)物分子數(shù)量衰減曲線及擬合Fig.4 Decay curves of reactant molecules and curve fitting

2.3初始反應(yīng)路徑

1）NO2的出現(xiàn)預(yù)示著硝基型含能材料高溫?zé)岱纸膺^程的能量勢壘被克服［3］。氣相硝基型含能材料熱分解的活化能主要由NO2的脫落決定［34］。而固相硝基型含能材料在脫去NO2的同時(shí)還將發(fā)生一系列雙分子鏈?zhǔn)椒磻?yīng)。Cohen［38］利用可測性設(shè)計(jì)（DFT）方法分析了氣相TNT高溫?zé)岱纸膺^程，并得出可能的3種初始反應(yīng)路徑：①C—NO2鍵斷裂形成NO2（g）以及C7H5O4N2；② 甲基上C—H鍵斷裂和鄰位NO2上的O結(jié)合形成H2O和C7H3O5N3（DNAn）；③TNT分子上—NO2異構(gòu)形成O—NO，最后形成NO（g）.并且當(dāng)溫度高于1 700 K時(shí)，氣相TNT分解的主要初始反應(yīng)路徑為—NO2鍵斷裂，其他的路徑可以忽略。

圖5　TNT/C10H8共晶和TNT單晶主要初級(jí)產(chǎn)物演化分布Fig.5 Distributions of main primary products for TNT/C10H8cocrystal and TNT single crystal

產(chǎn)物識(shí)別分析表明，NO2、NO和HONO為兩種含能材料的主要初級(jí)產(chǎn)物。圖5為TNT/C10H8共晶和TNT單晶高溫?zé)岱纸獬跏挤磻?yīng)路徑的比較。對(duì)于共晶，NO2、NO、HONO等初級(jí)產(chǎn)物主要來自于共晶內(nèi)TNT分子。為了直觀地比較共晶和單晶熱分解的初始反應(yīng)路徑以及共晶內(nèi)C10H8對(duì)TNT熱分解的影響，縱坐標(biāo)為各主要初級(jí)產(chǎn)物數(shù)量與初始TNT分子數(shù)量的比值（共晶：初級(jí)產(chǎn)物量/TNT分子數(shù)（32），TNT單晶：初級(jí)產(chǎn)物量/TNT分子數(shù)（64））。初級(jí)產(chǎn)物的演化分布顯示兩種含能材料在吸熱階段發(fā)生初始熱分解，且初級(jí)產(chǎn)物不斷累積，至峰值后，各產(chǎn)物逐漸地加入到次級(jí)反應(yīng)中，進(jìn)而形成最終的N2、H2O、CO2和 CO.NO2和 TNT—NO2，NO 和TNT—NO，HONO和TNT—HONO分布數(shù)量上的較大差異表明，并非TNT分子直接脫落形成各主要初級(jí)產(chǎn)物，而是 TNT-H/TNT+H與 C10H8分子（共晶）和TNT-H與TNT+H分子（TNT單晶）首先結(jié)合隨后脫落形成初級(jí)產(chǎn)物。反應(yīng)分子動(dòng)力學(xué)模擬已經(jīng)論證固相TNT高溫?zé)岱纸獾姆磻?yīng)路徑為相鄰的TNT分子間轉(zhuǎn)移H原子，形成TNT-H和TNT+H，此過程需要克服的能量勢壘為43 kcal/mol［34］，小于單分子TNT直接脫落—NO2所需的活化能（61.7 kcal/mol）［38］。并且TNT+H將進(jìn)一步分解形成HONO.對(duì)比共晶和TNT單晶內(nèi)各主要初級(jí)產(chǎn)物的演化分布，同溫度條件下，共晶熱分解形成的初級(jí)產(chǎn)物量整體上小于TNT熱分解形成的產(chǎn)物數(shù)量，主要原因是由于固相TNT熱分解的初始路徑為雙分子反應(yīng)機(jī)制，C10H8的加入，將限制TNT分子間的H原子轉(zhuǎn)移，進(jìn)而對(duì)初級(jí)產(chǎn)物的分布造成影響。熱分解初期C10H8不會(huì)發(fā)生開環(huán)斷裂，而是與共晶內(nèi)TNT形成的各產(chǎn)物結(jié)合，使得初級(jí)產(chǎn)物的量進(jìn)一步減少，這主要是由于C—C鍵離解能 （BDE）（100～150 kcal/mol），遠(yuǎn)遠(yuǎn)大于 C—N鍵（約 70 kcal/mol）、O—N（約50 kcal/mol）［38］。

2.4最終氣相產(chǎn)物分析

液態(tài)TNT高溫?zé)岱纸庾罱K氣相產(chǎn)物為H2O、CO2、N2、H2、NH3、CO、OH及大的團(tuán)簇［39］。本文模擬結(jié)果顯示共晶以及固相TNT高溫?zé)岱纸獾闹饕罱K產(chǎn)物為CO、CO2、N2和H2O，其中TNT高溫?zé)岱纸獾闹饕罱K產(chǎn)物分布與Furman的模擬結(jié)果［34］一致。比較TNT爆轟溫度條件下（T=3 500 K）熱分解最終產(chǎn)物與密封真空環(huán)境下 TNT（密度為1.53 g/cm3）爆轟實(shí)驗(yàn)［40］測得的產(chǎn)物進(jìn)行比較，產(chǎn)物類型一致，但前者最終產(chǎn)物的量整體上小于后者。究其原因，實(shí)驗(yàn)值為等熵膨脹后得到的氣體產(chǎn)物量，而本文采用NVT系綜，固定體積，且有限的模擬時(shí)間將共同造成最終產(chǎn)物的數(shù)量產(chǎn)生明顯差異。

比較至模擬結(jié)束時(shí)刻，兩種含能材料高溫?zé)岱纸庑纬傻闹饕罱K產(chǎn)物的數(shù)量，TNT單晶內(nèi)，平均每個(gè)TNT分子產(chǎn)生的CO、CO2和N2的數(shù)量高于共晶內(nèi)每個(gè)TNT分子產(chǎn)生的數(shù)量，而H2O的產(chǎn)量則相反。這主要是由于共晶內(nèi)C10H8分子上C—H鍵斷裂所需的能量為 C—H（約100 kcal/mol），小于C—C鍵斷裂需要的能量（約100～150 kcal/mol）［38］，因此，共晶高溫?zé)岱纸膺^程中，C10H8分子開環(huán)斷裂前，C—H鍵斷裂，脫落形成的H原子，將促進(jìn)共晶內(nèi)H2O的形成。另外，兩種含能材料高溫?zé)岱纸庾罱K產(chǎn)物（N2、H2O、CO和CO2）的演化分布同樣滿足指數(shù)函數(shù)特征［4］：

式中：C∞為至模擬結(jié)束時(shí)各主要產(chǎn)物分子數(shù)量的漸進(jìn)分布值；k2為產(chǎn)物形成的速率常數(shù)（ps-1）。利用上述指數(shù)函數(shù)對(duì)兩種含能材料熱分解最終產(chǎn)物分布曲線進(jìn)行擬合得到N2，H2O以及CO2的產(chǎn)生速率如表3所示（CO數(shù)量較少，未對(duì)其分布曲線擬合，擬合曲線如圖6所示）。

表3　主要產(chǎn)物產(chǎn)生速率比較Tab.3 Comparison of production rates for main products

圖6　主要最終產(chǎn)物演化分布及曲線擬合Fig.6 Time evolution curves of main final products and curve fitting

擬合結(jié)果顯示，TNT/C10H8共晶及TNT單晶高溫?zé)岱纸庾罱K產(chǎn)物的產(chǎn)生速率均隨著溫度的升高而加快。當(dāng)溫度為2000 K、2500 K時(shí)，至模擬結(jié)束時(shí)，TNT/C10H8共晶高溫?zé)峤猱a(chǎn)物未出現(xiàn)N2和CO2.當(dāng)溫度為3 000 K～4 500 K，H2O分子的產(chǎn)生速率＞N2＞CO2.對(duì)于TNT單晶，H2O分子的產(chǎn)生速率＞N2≈CO2.并且，各溫度條件下共晶內(nèi)H2O的產(chǎn)生速率均高于TNT單晶。此現(xiàn)象可用上述共晶內(nèi)C10H8分子上C—H鍵斷裂形成的H原子促進(jìn)了H2O的產(chǎn)生速率加以解釋。

3 結(jié)論

ReaxFF分子動(dòng)力學(xué)方法已經(jīng)在含能材料熱分解、沖擊起爆、爆轟穩(wěn)定性等問題得到廣泛應(yīng)用，且較實(shí)驗(yàn)技術(shù)得到了更多的微觀細(xì)節(jié)，但當(dāng)前的研究對(duì)象多為常規(guī)含能材料。本文利用增加了長程修正項(xiàng)Elg的反應(yīng)力場分析了TNT/C10H8共晶初始高溫?zé)岱纸馇闆r，并通過TNT單晶高溫?zé)岱纸膺M(jìn)行了對(duì)比分析。研究結(jié)果表明，同溫度條件下，TNT單晶放熱量高于TNT/C10H8共晶，并且利用經(jīng)典反應(yīng)速率方程擬合得到共晶高溫?zé)岱纸獬跫?jí)吸熱和次級(jí)放熱反應(yīng)的活化能分別為35.7 kcal/mol和56.1 kcal/mol.共晶熱分解初級(jí)吸熱階段的活化能和TNT單晶基本相同，這主要是由于共晶內(nèi)TNT分子上—NO2鍵斷裂，共晶的熱分解能量勢壘已經(jīng)越過。次級(jí)放熱反應(yīng)的活化能明顯高于TNT單晶，這主要是由于固相含能材料高溫?zé)岱纸獾某跏挤磻?yīng)路徑為雙分子反應(yīng)機(jī)制，C10H8的加入抑制了共晶熱分解的放熱反應(yīng)。此結(jié)果可以與擬合得到的共晶和TNT單晶的反應(yīng)放熱以及共晶內(nèi)TNT和TNT單晶內(nèi)分子數(shù)量的衰減速率相互驗(yàn)證。共晶含能材料的熱分解特征主要受晶體內(nèi)材料自身物性以及二者相互作用機(jī)制的影響，并且具有類三明治結(jié)構(gòu)的共晶與純物理混合結(jié)構(gòu)的含能材料熱分解特征具有較明顯的差別［41］。

產(chǎn)物識(shí)別分析顯示，TNT/C10H8共晶和TNT共晶的初級(jí)反應(yīng)產(chǎn)物為NO2、NO和HONO，并且通過NO2和TNT—NO2，NO和TNT—NO，HONO和TNT—HONO分布數(shù)量的比較論證了固相共晶和TNT單晶的初始反應(yīng)路徑為雙分子反應(yīng)機(jī)制。共晶和TNT單晶熱分解的主要最終產(chǎn)物均為N2、H2O、CO2和CO.這主要是由于C10H8中C—C鍵斷裂需要較高的能量，因此，C10H8熱分解的初始反應(yīng)路徑為C—H鍵斷裂，并且此過程將會(huì)促進(jìn)共晶內(nèi)H2O的形成。

當(dāng)前，利用ReaxFF分子動(dòng)力學(xué)通過施加Chapman-Jouget點(diǎn)溫度模擬含能材料爆轟 Chapman-Jouget點(diǎn)后方的初始反應(yīng)路徑和最終產(chǎn)物，但相對(duì)爆轟來說，此過程并沒有涉及壓力條件，另外，還需要含能材料爆轟Chapman-Jouget點(diǎn)密度也是一個(gè)重要因素。因此，高溫?zé)岱纸饽M僅能對(duì)含能材料爆轟波內(nèi)的反應(yīng)情況提供參考。除此之外，利用ReaxFF分子動(dòng)力學(xué)模擬缺乏自由邊界條件，并且對(duì)系統(tǒng)體積加以固定，這將與實(shí)驗(yàn)結(jié)果造成差異，因此，未來需要建立真實(shí)的初始和邊界條件，同時(shí)不斷完善勢函數(shù)以模擬含能材料爆轟波結(jié)構(gòu)內(nèi)的反應(yīng)情況。

（References）

［1］ Fickett W，Davis W C.Detonation：theory and experiment［M］. New York：Dover Publications，2000：44-51.

［2］ van Duin A C T，Dasgupta S，Lorant F，et al.ReaxFF：a reactive force field for hydrocarbons［J］.Journal of Physical Chemistry A，2001，105（41）：9396-9409.

［3］ Strachan A，Kober E M，van Duin A C T，et al.Thermal decomposition of RDX from reactive molecular dynamics［J］.Journal of Chemical Physics，2005，122（5）：054502.

［4］ Rom N，Zybin S V，van Duin A C T，et al.Density-dependent liquid nitromethane decomposition：molecular dynamics simulations based on ReaxFF［J］.Journal of Physical Chemistry A，2011，115（36）：10181-10202.

［5］ Zhang L Z，Zybin S V，van Duin A C T，et al.Carbon cluster formation during thermal decomposition of octahydro-1，3，5，7 tetranitro-1，3，5，7-tetrazocine and 1，3，5-triamino-2，4，6-trinitrobenzene high explosives from ReaxFF reactive molecular dynamics simulations［J］.JournalofPhysicalChemistryA，2009， 113（40）：10619-10640.

［6］ Shan T R，Wixom R R，Mattsson A E，et al.Atomistic simulation of orientation dependence in shock-induced initiation of pentaerythritol tetranitrate［J］.Journal of Physical Chemistry B，2013，117（3）：928-936.

［7］ Wen Y S，Xue X G，Zhou X Q，et al.Twin induced sensitivity enhancement of HMX versus shock：a molecular reactive force field simulation［J］.Journal of Physical Chemistry C，2013，117（46）：24368-24374.

［8］ Ge N N，Wei Y K，Song Z F，et al.Anisotropic responses and initial decomposition of condensed-phase β-HMX under shock loadings via molecular dynamics simulations in conjunction with multiscale shock technique［J］.Journal of Physical Chemistry B，2014，118（29）：8691-8699.

［9］ Budzien J，Thompson A P，Zybin S V.Reactive molecular dynamics simulations of shock through a single crystal of pentaerythritol tetranitrate［J］.Journal of Physical Chemistry B 2009，113（40）：13142-13151.

［10］ Cai Y，Zhao F P，An Q，et al.Shock response of single crystal and nanocrystalline pentaerythritol tetranitrate：Implications to hotspot formation in energetic materials［J］.Journal of Chemical Physics，2013，139（16）：164704.

［11］ Zhang L Z，Zybin S V，Van Duin A C T，et al.Modeling high rate impact sensitivity of perfect RDX and HMX crystals by ReaxFF reactive dynamics［J］.Journal of Energetic Materials，2010，28（1）：92-127.

［12］ He L，Sewell T D，Thompson D L.Molecular dynamics simulations of shock waves in oriented nitromethane single crystals［J］. Journal of Chemical Physics，2011，134（12）：124506.

［13］ An Q，Zybin S V，Goddard III W A，et al.Elucidation of the dynamics for hot-spot initiation at nonuniform interfaces of highly shocked materials［J］.Physical Review B，2011，84（22）：220101.

［14］ An Q，Goddard III W A，Zybin S V，et al.Highly shocked polymer bonded explosives at a nonplanar interface：hot-spot formation leading to detonation［J］.Journal of Physical Chemistry C，2013，117（50）：26551-26561.

［15］ Ken-Ichi N，Kalia R K，Aiichiro N，et al.Dynamic transition in the structure of an energetic crystal during chemical reactions at shock front prior to detonation［J］.Physical Review Letters，2007，99（14）：12243-12254.

［16］ Ge N N，Wei Y K，Ji G F，et al.Initial decomposition of the condensed-phase β-HMX under shock waves：molecular dynamics simulations［J］.Journal of Physical Chemistry B，2012，116（46）：13696-13704.

［17］ Landenberger K B.Cocrystallization of energetic materials［D］. Michigan：The University of Michigan.2013.

［18］ Bolton O，Matzger J A.Improved stability and smart-material functionality realized in an energetic cocrystal［J］.Angewandte Chemie International Edition，2011，50（38）：8960-8963.

［19］ Yang Z W，Li H Z，Huang H，et al.Preparation and performance of a HNIW/TNT cocrystal explosive［J］.Propellants ExplosivesPyrotechnics，2013，38（4）：495-510.

［20］ Bolton O，Simke L R，Pagoria P F，et al.High power explosive with good sensitivity：a 2：1 cocrystal of CL-20：HMX［J］.Crystal Growth&Design，2012，12（9）：4311-4314.

［21］ Yang Z W，Li H Z，Zhou X Q，et al.Characterization and properties of a novel energetic-energetic cocrystal explosive composed of HNIW and BTF crystal［J］.Crystal Growth&Design，2012，12（11）：5155-5158.

［22］ Wei C X，Huang H，Duan X H，et al.Structures and properties prediction of HMX/TATB co-crystal［J］.Propellants Explosives Pyrotechnics，2011，36（5）：416-423.

［23］ Landenberger K B，Matzger A J.Cocrystal engineering of a prototype energetic material：supramolecular chemistry of 2，4，6-trinitrotoluene［J］.Crystal Growth&Design，2010，10（12）：5341-5347.

［24］ 劉海，李啟楷，何遠(yuǎn)航.CL20-TNT共晶高溫?zé)峤獾腞eaxFF/lg反應(yīng)力場分子動(dòng)力學(xué)模擬［J］.物理學(xué)報(bào)，2013，62（20）：208202. LIU Hai，LI Qi-kai，HE Yuan-hang.Pyrolysis of CL20-TNT cocrystal from ReaxFF/lg reactive molecular dynamics simulations［J］.Acta Physica Sinica，2013，62（20）：208202.（in Chinese）

［25］ Guo D Z，An Q，Goddard W A，et al.Compressive shear reactive molecular dynamics studies indicating that cocrystals of TNT/ CL-20 decrease sensitivity［J］.Journal of Physical Chemistry C，2014，118（51）：30202-30208

［26］ Sun T，Xiao J J，Liu Q，et al.Comparative study on structure，energetic and mechanical properties of aε-CL-20/HMX cocrystal and its composite with molecular dynamics simulation［J］.Journal of Materials Chemistry A，2014，34（2）：13898-13904.

［27］ Groom C.The Cambridge Crystallographic data centre（CCDC）［DB/OL］.［2015-8-22］.http：∥www.ccdc.cam.ac.uk.

［28］ Plimpton S.Fast parallel algorithms for short-range molecular dynamics［J］.Journal of Computational Physics，1995，117（1）：1-19.

［29］ Plimpton S.LAMMPS Molecular dynamics simulator［CP/OL］. ［2015-7-26］.http：∥lammps.sandia.gov.

［30］ Liu L C，Liu Y，Zybin S V，et al.ReaxFF-lg：correction of the ReaxFF reactive force field for London dispersion，with applications to the equations of state for energetic materials［J］.Journal of Physical Chemistry A，2011，115（40）：11016-11022.

［31］ 宋華杰，周婷婷，黃風(fēng)雷，等.低壓長脈沖載荷下β-HMX單晶滑移系的微觀物理化學(xué)響應(yīng)［J］.物理化學(xué)學(xué)報(bào)，2014，30（11）：2024-2034. SONG Hua-jie，ZHOU Ting-ting，HUANG Feng-lei，et al.Microscopic physical and chemical responses of slip systems in the β-HMX single crystal under low pressure and long pulse loading ［J］.Acta Physico-Chimica Sinica，2014，30（11）：2024-2034.（in Chinese）

［32］ Zhou T T，Song H J，Liu Y，et al.Shock initiated thermal and chemical responses of HMX crystal from ReaxFF molecular dynamics simulation［J］.Physical Chemistry Chemical Physics，2014，16（27）：13914-13931.

［33］ 劉海，李啟楷，何遠(yuǎn)航.高速?zèng)_擊壓縮梯恩梯的分子動(dòng)力學(xué)模擬［J］.力學(xué)學(xué)報(bào)，2015，47（1）：174-179. LIU Hai，LI Qi-kai，HE Yuan-hang.Molecular dynamics simulations of high velocity shock compressed TNT［J］.Chinese Journal of Theoretical and Applied Mechani，2015，47（1）：174-179.（in Chinese）

［34］ Furman D，Kosloff R，Dubnikova F，et al.Decomposition of condensed phase energetic materials：interplay between uni-and bimolecular mechanisms［J］.Journal of the American Chemical Society，2014，136（11）：4192-4200.

［35］ Robertson A J B.The decomposition，boiling and explosion of trinitrotoluene at high temperatures［J］.Transactions of the Faraday Society，1948，44：977-983.

［36］ Zinn J，Rogers R N.Thermal initiation of explosivers［J］.Journal of Physical Chemistry 1962，66（12）：2646-2653.

［37］ Beckmann J W，Wilkes J S，McGuire R R.2，4，6-Trinitrotoluene thermal decomposition：kinetic parameters determined by the isothermal differential scanning calorimetry technique［J］.Thermochimica Acta，1977，19（1）：111-118.

［38］ Cohen R，Zeiri Y，Wurzberg E，et al.Mechanism of thermal unimolecular decomposition of TNT（2，4，6-trinitrotoluene）：a DFT study［J］.Journal of Physical Chemistry A，2007，111（43）：11074-11083.

［39］ Rom N，Hirshberg B，Zeiri Y，et al.First-principles-based reaction kinetics for decomposition of hot，dense liquid TNT from ReaxFF multiscale reactive dynamics simulations［J］.Journal of Physical Chemistry C，2013，117（41）：21043-21054.

［40］ Ornellas D L.Calorimetric determination of the heat and products of detonation for explosives，AD-A409329［R］.Livermore，CA：Lawrence Livermore National Laboratory，1982.

［41］ Guo D Z，An Q，Zybin S V，et al.The co-crystal of TNT/CL-20 leads to decreased sensitivity toward thermal decomposition from first principles based reactive molecular dynamics［J］.Journal of Materials Chemistry A，2015，10（3）：5409-5419.

Initial Thermal Decomposition of TNT/C10H8Cocrystal：Reactive Molecular Dynamics Simulations

LIU Hai，HE Yuan-hang
（State Key Laboratory of Explosion Science and Technology，Beijing Institute of Technology，Beijing 100081，China）

The discussion of shock initiation and detonation with ReaxFF reactive force field has centered nearly around conventional energetic materials.The comprehensive molecular dynamics of thermal decomposition of condensed phase TNT/C10H8cocrystal with the ReaxFF/lg potential function are studied，which adds a long range correction term to the total energy expression of original ReaxFF.The thermal decomposition of TNT single crystal is performed similarly to determine the main effect of C10H8on the overall thermal characteristics of cocrystal.Furthermore，the activation energies associated with endothermic and exothermic stages of cocrystal thermal decomposition are found to be 35.7 kcal/mol and 56.1 kcal/mol，respectively.Its activation energy at primary endothermic stage is compatible with that of TNT single crystal，but the activation energy at secondary exothermic stage is higher than the value calculated for TNT.In addition，the heat output calculated during the exothermic step of TNT/C10H8cocrystal is less than that of TNT at the same temperature.The decay rate of reactants by curve fitting shows thatC10H8inhibits the decomposition of TNT in cocrystal.NO2，NO and HONO are major primary decomposition products by the way of identification analysis.The distribution of NO2/TNT—NO2，NO/TNT—NO，HONO/TNT—HONOshows that the decomposition proceeds via bimolecular reaction.Main final products of cocrystal are N2，H2O，CO2，and CO，and the initial decomposition routes of C10H8is C—H bond dissociation as a result of C—C bond dissociation energy of C10H8in cocrystal is higher than C—H bond，and this process may promote the formation of H2O.

ordnance science and technology；cocrystal；pyrolysis；reaction dynamics；ReaxFF reactive force field

O64

A

1000-1093（2016）03-0414-10

10.3969/j.issn.1000-1093.2016.03.005

2015-04-03

劉海（1985—），男，博士研究生。E-mail：tristan_l@bit.edu.cn；何遠(yuǎn)航（1964—），男，教授，博士生導(dǎo)師。E-mail：heyuanhang@bit.edu.cn