梯恩梯/萘共晶初始熱分解的反應分子動力學模擬

劉海，何遠航

（北京理工大學爆炸科學與技術國家重點實驗室，北京100081）

梯恩梯/萘共晶初始熱分解的反應分子動力學模擬

劉海，何遠航

（北京理工大學爆炸科學與技術國家重點實驗室，北京100081）

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

兵器科學與技術；共晶；高溫熱分解；反應動力學；ReaxFF反應力場

0 引言

極端條件下固相含能材料的物理化學響應問題是爆轟學的核心內容之一，并且，近幾十年的實驗測試和理論分析多圍繞常規含能材料沖擊起爆和爆轟波的宏觀問題，如爆壓、爆速、爆轟波結構以及穩定性等［1］。極端溫度/壓力條件下固相含能材料的初始反應路徑是研究沖擊起爆問題的基礎，同樣也與含能材料爆轟波能量輸出結構直接關聯。由于固相含能材料的沖擊起爆及爆轟問題涉及極短的時間和極小的空間尺度，這給現有的實驗技術帶來了很大的挑戰。近期，各國學者利用美國加州理工學院Goddard小組提出的ReaxFF反應力場分子動力學方法［2］已經就常規含能材料高溫熱分解的初始反應路徑［3-5］以及沖擊起爆機制［6-16］進行了廣泛而深入地研究。

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

1 模型及計算細節

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

ReaxFF/lg是在ReaxFF勢函數的基礎上增加了長程修正項Elg.ReaxFF/lg勢函數在含能材料晶體結構以及狀態方程的預測相比于ReaxFF勢函數更加接近實驗值［30］。當前，ReaxFF/lg勢函數已經在含能材料感度分析［31］、熱分解［24］、熱沖擊［32］、沖擊壓縮［33］、沖擊起爆［6-7］等問題得到應用。

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

2 結果及分析

2.1勢能演化及曲線擬合

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

含能材料高溫熱分解過程中，勢能的變化和系統內化學反應進程緊密關聯。各溫度條件下，兩種含能材料高溫熱分解過程中系統勢能隨時間的演化分布表明，系統首先吸收熱量，勢能快速達到峰值。此階段對應分子由平動模態轉換為振動模態，同時伴隨著初始化學反應［34］，此階段可視為“平衡-誘導期（tE-I）”［3］，且為吸熱反應階段。隨后，伴隨著次級放熱化學反應的進行，勢能開始衰減并逐漸達到平衡，并且系統勢能衰減曲線可用下述指數函數［3-5］進行擬合：

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

式中：U（t）為t時刻系統勢能（kcal/mol）；U∞為勢能漸進平衡值（kcal/mol）；ΔUexo為次級反應放熱量（kcal/mol）；tE-I為平衡-誘導期（ps）；τ為特征反應時間（ps，反應速率的倒數）。部分擬合曲線如圖3黑線所示，擬合得到的各動力學以及熱力學參數如表1.

擬合結果顯示，隨著溫度的升高，含能材料高溫熱分解的“平衡-誘導期”和特征反應時間越短，也就是初級吸熱反應，勢能增加至峰值的時間以及次級放熱反應達到平衡所需的時間越少。并且，同溫度條件下，TNT高溫熱解的放熱量均明顯高于TNT/ C10H8共晶。

系統勢能峰值為高溫熱分解過程中含能材料吸熱和放熱反應的分界點。為了比較TNT/C10H8共晶和TNT高溫熱分解吸熱和放熱反應的活化能，這里采用經典 Arrhenius方程，同時利用表1中 TNT/ C10H8共晶和TNT高溫熱解勢能衰減曲線（TNT/C10H8共晶：3000 K、3 500 K、4 000 K、4 500 K；TNT：3 000 K、3 500 K、4 000 K）擬合得到“平衡-誘導期”、特征反應時間進行擬合：

式中：Ea為活化能；A為指前因子；R為理想氣體常數。變換（2）式得到：lnτ=（Ea/RT）-lnA，可見lnτ 與1/T呈線性關系。擬合結果為TNT/C10H8共晶和TNT高溫熱解的吸熱和放熱反應階段的活化能分別為

圖3　兩種含能材料高溫熱分解過程吸熱和放熱反應階段的活化能Fig.3 Activation energy of endothermic and exothermic stages of thermal decomposition for cocrystal and TNT single crystal

表1　利用指數函數擬合得到的各動力學以及熱力學參數Tab.1 Kinetic and thermodynamic parameters by exponential function fitting

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

2.2反應物衰減速率比較

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

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

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

表2　反應物衰減速率比較Tab.2 Attenuation rates of reactants

圖4　反應物分子數量衰減曲線及擬合Fig.4 Decay curves of reactant molecules and curve fitting

2.3初始反應路徑

1）NO2的出現預示著硝基型含能材料高溫熱分解過程的能量勢壘被克服［3］。氣相硝基型含能材料熱分解的活化能主要由NO2的脫落決定［34］。而固相硝基型含能材料在脫去NO2的同時還將發生一系列雙分子鏈式反應。Cohen［38］利用可測性設計（DFT）方法分析了氣相TNT高溫熱分解過程，并得出可能的3種初始反應路徑：①C—NO2鍵斷裂形成NO2（g）以及C7H5O4N2；② 甲基上C—H鍵斷裂和鄰位NO2上的O結合形成H2O和C7H3O5N3（DNAn）；③TNT分子上—NO2異構形成O—NO，最后形成NO（g）.并且當溫度高于1 700 K時，氣相TNT分解的主要初始反應路徑為—NO2鍵斷裂，其他的路徑可以忽略。

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

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

2.4最終氣相產物分析

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

比較至模擬結束時刻，兩種含能材料高溫熱分解形成的主要最終產物的數量，TNT單晶內，平均每個TNT分子產生的CO、CO2和N2的數量高于共晶內每個TNT分子產生的數量，而H2O的產量則相反。這主要是由于共晶內C10H8分子上C—H鍵斷裂所需的能量為 C—H（約100 kcal/mol），小于C—C鍵斷裂需要的能量（約100～150 kcal/mol）［38］，因此，共晶高溫熱分解過程中，C10H8分子開環斷裂前，C—H鍵斷裂，脫落形成的H原子，將促進共晶內H2O的形成。另外，兩種含能材料高溫熱分解最終產物（N2、H2O、CO和CO2）的演化分布同樣滿足指數函數特征［4］：

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

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

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

擬合結果顯示，TNT/C10H8共晶及TNT單晶高溫熱分解最終產物的產生速率均隨著溫度的升高而加快。當溫度為2000 K、2500 K時，至模擬結束時，TNT/C10H8共晶高溫熱解產物未出現N2和CO2.當溫度為3 000 K～4 500 K，H2O分子的產生速率＞N2＞CO2.對于TNT單晶，H2O分子的產生速率＞N2≈CO2.并且，各溫度條件下共晶內H2O的產生速率均高于TNT單晶。此現象可用上述共晶內C10H8分子上C—H鍵斷裂形成的H原子促進了H2O的產生速率加以解釋。

3 結論

ReaxFF分子動力學方法已經在含能材料熱分解、沖擊起爆、爆轟穩定性等問題得到廣泛應用，且較實驗技術得到了更多的微觀細節，但當前的研究對象多為常規含能材料。本文利用增加了長程修正項Elg的反應力場分析了TNT/C10H8共晶初始高溫熱分解情況，并通過TNT單晶高溫熱分解進行了對比分析。研究結果表明，同溫度條件下，TNT單晶放熱量高于TNT/C10H8共晶，并且利用經典反應速率方程擬合得到共晶高溫熱分解初級吸熱和次級放熱反應的活化能分別為35.7 kcal/mol和56.1 kcal/mol.共晶熱分解初級吸熱階段的活化能和TNT單晶基本相同，這主要是由于共晶內TNT分子上—NO2鍵斷裂，共晶的熱分解能量勢壘已經越過。次級放熱反應的活化能明顯高于TNT單晶，這主要是由于固相含能材料高溫熱分解的初始反應路徑為雙分子反應機制，C10H8的加入抑制了共晶熱分解的放熱反應。此結果可以與擬合得到的共晶和TNT單晶的反應放熱以及共晶內TNT和TNT單晶內分子數量的衰減速率相互驗證。共晶含能材料的熱分解特征主要受晶體內材料自身物性以及二者相互作用機制的影響，并且具有類三明治結構的共晶與純物理混合結構的含能材料熱分解特征具有較明顯的差別［41］。

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

當前，利用ReaxFF分子動力學通過施加Chapman-Jouget點溫度模擬含能材料爆轟 Chapman-Jouget點后方的初始反應路徑和最終產物，但相對爆轟來說，此過程并沒有涉及壓力條件，另外，還需要含能材料爆轟Chapman-Jouget點密度也是一個重要因素。因此，高溫熱分解模擬僅能對含能材料爆轟波內的反應情況提供參考。除此之外，利用ReaxFF分子動力學模擬缺乏自由邊界條件，并且對系統體積加以固定，這將與實驗結果造成差異，因此，未來需要建立真實的初始和邊界條件，同時不斷完善勢函數以模擬含能材料爆轟波結構內的反應情況。

（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］ 劉海，李啟楷，何遠航.CL20-TNT共晶高溫熱解的ReaxFF/lg反應力場分子動力學模擬［J］.物理學報，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］ 宋華杰，周婷婷，黃風雷，等.低壓長脈沖載荷下β-HMX單晶滑移系的微觀物理化學響應［J］.物理化學學報，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］ 劉海，李啟楷，何遠航.高速沖擊壓縮梯恩梯的分子動力學模擬［J］.力學學報，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；何遠航（1964—），男，教授，博士生導師。E-mail：heyuanhang@bit.edu.cn