彭亞晶 孫爽 劉偉娜 劉宇輝

(渤海大學物理科學與技術學院, 錦州 121000)

1 引 言

環(huán)三亞甲基三硝胺(RDX)是一種高性能的含能化合物, 廣泛應用于多種聚合物鍵合炸藥(PBXs)和推進劑配方中[1].由于晶體生長過程中的各種不利條件, 通常會在晶體內(nèi)部形成一些缺陷, 包括空位、雜質(zhì)和位錯等, 這對其初始反應性能有一定的影響[2-4].當含能材料中含有大量缺陷時, 對外界熱或沖擊的感度將顯著下降[5].Wang等[6]和Kuklja等[7]運用第一性原理的方法對含空位和位錯缺陷的RDX晶體的電子結構性質(zhì)進行了研究, 結果表明空位缺陷的存在會導致晶體光學帶隙變窄, 并降低了其絕緣體相態(tài)向金屬相態(tài)轉變的臨界壓力.我們先前的工作對含有分子空位的RDX晶體的幾何結構、電子結構及振動特性進行了詳細地分析,表明了分子空位缺陷的存在使體系能帶隙變小, 并使缺陷附近的分子結構松弛, 電荷分布增多, 反應活性增強[8].Chakraborty等[9]利用密度泛函理論計算了HMX和RDX含能材料的分解途徑, 表明了N—NO2消去和氫轉移形成HONO的能壘較低(約39 kcal/mol); Ge等[10,11]利用多尺度沖擊技術(MSST)對凝聚相β-HMX 和RDX進行分子動力學模擬, 給出N—NO2鍵的解離和C—N鍵斷裂是RDX的主要分解途徑, 并指出單一RDX分子在11 km/s速度沖擊下, 主要產(chǎn)物包括NO2, NO,N2O, CO和N2.鄭朝陽和趙紀軍[12]利用第一性原理分子動力學方法研究了RDX的熱解初始反應機制, 表明了分解方式主要包括質(zhì)子轉移、C—N鍵斷裂和N—NO2鍵斷裂等方式, 并提出沖擊加載可能帶來C—H鍵斷裂的反應機制.然而, 在沖擊加載下, 含空位的RDX晶體內(nèi)部微觀結構如何變化? 反應機理怎樣? 還沒有被澄清.由于晶體內(nèi)部微觀現(xiàn)象較難通過實驗分辨, 因此, 第一性原理和分子動力學方法被廣泛地應用來研究含能材料的微觀動態(tài)過程[13-16].

本文主要利用基于ReaxFF (reactive force field)的分子動力學方法結合多尺度沖擊技術(MSST)研究沖擊加載下RDX的完美晶體和分子空位晶體的初始動態(tài)響應特性, 進而預測兩種RDX晶體可能的初始反應機制, 為提高含能材料的安全性能提供必要的微觀信息.

2 模型構建與計算方法

ReaxFF是一種以第一性原理為基礎的反應力場[17-22], 它能夠較好地模擬含能材料在熱或沖擊作用下的化學反應過程.ReaxFF的參數(shù)是用量子力學計算或?qū)嶒灚@得的, 具有量子力學的準確性[19,20].因此, 對于較小的模擬尺度, ReaxFF分子動力學仍然能夠給出較好的模擬結果[22].使用Material Studio建立完美RDX超胞(2 × 1 × 1)和含有1個分子空位的超胞, 如圖1(a)和圖1(b)所示.A, B和C代表3個晶格方向.采用了密度泛函理論中的廣義梯度近似方法及Perdew-Burke-Ernzerhof泛函[23], 對RDX完美超晶胞以及分子空位超晶胞進行幾何結構優(yōu)化.由于RDX超晶胞的原子數(shù)較多, 計算量較大, 所以選取了Gammapoint, 即k = 1 × 1 × 1.優(yōu)化完美RDX超胞的晶格參數(shù)為a = 26.731756 ?, b = 11.500793 ?,c = 10.786047 ?, α = β = γ = 90°, 與文獻[24]報道的數(shù)據(jù)基本相符.在未加載沖擊波情況下, 采用NVE系蹤并結合Berendsen 控溫方式模擬RDX平衡過程1 ps得到穩(wěn)定的RDX晶體結構, 并以此作為沖擊加載的初始結構, 初始溫度為300 K, 初始壓力設為1個大氣壓, 沿晶體A方向進行多尺度沖擊壓縮模擬.時間步長為0.01 fs, 原子截斷半徑取1.5 ?.分析了沖擊速度為10, 11和13 km/s的沖擊波作用下完美RDX晶體和空位晶體的微觀動態(tài)過程, 計算可能參與反應的原子間的徑向分布函數(shù), 揭示兩種晶體結構的可能的初始反應機制.

圖1 構建的完美RDX超胞(2 × 1 × 1)(a)和含1個分子空位的RDX超胞(b), 紅、藍、灰、白分別表示氧、氮、碳和氫原子Fig.1.Constructed perfect RDX supercell (2 × 1 × 1) (a)and RDX supercell with a molecular vacancy (b).Red,blue, gray and white represent oxygen, nitrogen, carbon,and hydrogen atom, respectively.

3 結果與討論

3.1 晶體內(nèi)部結構的動態(tài)變化

圖2(a)和圖2(b)分別為完美RDX晶胞和空位晶胞在沖擊速度為11 km/s的沖擊波作用下的動態(tài)變化過程.可見, 隨著沖擊波作用時間的增加,晶胞在A方向上逐漸被壓縮, 尺寸減小.晶體中的分子逐漸從有序排列到無序狀態(tài).從圖2(b)可見,在沖擊作用下, 空位附近的分子發(fā)生較大移動, 并逐漸填充到空位部分.圖3(a)和圖3(b)分別為10 km/s和11 km/s沖擊速度加載下RDX完美晶胞和空位晶胞的內(nèi)能、體積比、壓力和溫度隨時間的變化情況.顯然, 隨著時間的增加, 兩種情況下體系的熱力學參數(shù)逐漸趨于穩(wěn)定.反應溫度均穩(wěn)定在5000 K左右, 壓強大約維持在70 GPa左右, 沖擊速度越大, 體積比減小的越快.這些與文獻[25,26]報道基本符合.此外, 圖4給出了沖擊波加載條件下的沖擊速度和粒子速度關系與先前實驗和理論數(shù)據(jù)的比較[27,28].可見, 完美晶體的沖擊速度-粒子速度關系與目前實驗和理論數(shù)據(jù)符合較好.空位晶體的粒子速度比完美晶體的稍高, 但相差不多, 表明了計算結果的有效性.

圖2 沖擊速度為11 km/s加載下, 完美RDX超胞(a)和空位RDX超胞(b)在不同時刻的體系結構, 圖中紅色代表氧原子、黑色代表碳原子、藍色代表氮原子、綠色代表氫原子Fig.2.System structure of perfect RDX supercell (a) and RDX supercell with a molecular vacancy (b) at different time under shock velocity of 11 km/s.The red for oxygen,the black for carbon, the blue for nitrogen and the green for hydrogen.

圖3 沖擊加載下, RDX完美超胞(a)和含1個分子空位超胞(b)的總能量、體積比、壓強和溫度隨時間的變化Fig.3.Changes of total energy, volume ratio, pressure and temperature vs.time for RDX perfect supercell (a) and RDX supercell with a molecular vacancy (b) under shock loading.

圖4 RDX沖擊速度與粒子速度對應關系, 其中, 黑色方形數(shù)據(jù)為起爆前(沖擊速度低于8.6 km/s)實驗數(shù)據(jù)和起爆后的理論計算數(shù)據(jù)[27,28], 紅色圓形和藍色三角形數(shù)據(jù)為本計算的完美晶體和空位晶體的相應數(shù)據(jù)Fig.4.Shock velocity vs.particle velocity for RDX.Here,the black square data are the experimental data before detonation (below 8.6 km/s in shock velocity) and the theoretical calculation data after detonation[27,28], and the red circle and blue triangle data are respectively that of the perfect crystal and the vacancy crystal in this calculation.

為了明確沖擊過程中RDX分子結構的具體變化情況, 分別給出了沖擊速度為11 km/s作用下完美晶胞和分子空位晶胞在不同時刻的分子結構的情況, 如圖5和圖6所示.從圖5(a)可見,在1.015 ps時刻, 完美RDX晶體中的N—NO2鍵發(fā)生了斷裂, 而在1.125 ps時刻, C—N鍵發(fā)生了斷裂(如圖5(b)所示).這表明, 在沖擊加載下, 完美RDX晶體中N—NO2鍵首先發(fā)生斷裂, 隨后C—N鍵發(fā)生斷裂.從圖6(a)可以看出, 含分子空位缺陷的RDX晶體中N—NO2鍵斷裂時刻為1.010 ps,而C—N鍵斷裂時刻為1.043 ps (如圖6(b)).顯然,空位晶體中斷鍵時間要比完美晶體中的斷鍵稍早.這表明了分子空位缺陷的存在使含能材料對沖擊加載更敏感, 更容易發(fā)生化學反應.這與目前文獻[4, 5, 8]報道的結果是一致的.

圖5 在沖擊速度為11 km/s作用下, RDX完美晶胞中N—NO2鍵隨時間演化情況(a)和C—N鍵隨時間演化情況(b), 其中藍綠色代表氮、紫色代表氧、綠色代表氫、橘紅色代表碳Fig.5.Evolution of N—NO2 bonds with time (a) and that of C—N bonds with time (b) in RDX perfect cell under shock velocity of 11 km/s.The blue-green represents nitrogen atom, the purple represents oxygen atom, the green represents hydrogen atom and the tangerine represents carbon atom.

圖6 在沖擊速度為11 km/s作用下, 含分子空位的RDX晶胞中N—NO2鍵隨時間的演化情況(a)和C—N鍵隨時間的演化情況(b), 其中藍綠色代表氮、紫色代表氧、綠色代表氫、橘紅色代表碳Fig.6.Evolution of N—NO2 bonds with time (a) and that of C—N bonds with time (b) in RDX cell with a molecular vacancy under shock velocity of 11 km/s.The blue-green represents nitrogen atom, the purple represents oxygen atom, the green represents hydrogen atom and the tangerine represents carbon atom.

3.2 沖擊加載下完美RDX晶體的初始反應分析

從RDX分子本身結構和可能產(chǎn)生的產(chǎn)物碎片結構入手[11,15,16,29,30], 計算可能發(fā)生反應的原子間的徑向分布函數(shù), 如N—N, C—N, C—H, H—O等.圖7為計算獲得的完美RDX晶胞在10, 11和13 km/s沖擊速度作用下, 主要原子間徑向分布函數(shù)情況.為了分析原子間可能的成鍵概率情況, 考慮了徑向分布函數(shù)相關峰的展寬效應, 即對幾個主峰所包圍的面積進行了計算, 結果見表1.可見, 在不同沖擊速度下, 由于徑向分布函數(shù)的峰寬變化較小, 其所包圍的面積基本由峰值決定, 峰值高的所包圍的面積較大, 代表了該處原子成鍵數(shù)量較多.從圖7(a)可以看出, 對于氮氮原子, 當沖擊速度為10 km/s時, 在徑向半徑為1.55 ?附近, 其徑向分布函數(shù)gN—N(r)出現(xiàn)最大峰值.該徑向長度符合RDX晶體中的N—NO2鍵的成鍵范圍, 表明該狀態(tài)下, 氮-氮原子主要以N—NO2形式存在.隨著沖擊速度的增加, 1.55 ?附近徑向分布函數(shù)的峰值逐漸降低, 而表征氮-氮三鍵的1.18 ?附近徑向分布函數(shù)的峰值逐漸升高.這表明在沖擊加載下,RDX中的N—NO2鍵將發(fā)生部分斷裂, 形成產(chǎn)物氮氣, 并隨著沖擊速度的增加而斷裂程度增大, 最終形成氮氣的量也增多; 從圖7(b)可見, 在 1.45 ?附近, 碳-氮原子間的徑向分布函數(shù)gC—N(r)出現(xiàn)峰值, 該距離屬于RDX 環(huán)鏈中C—N 鍵的成鍵范圍.隨著沖擊速度的增加, gC—N(r)峰值下降, 表明了C—N鍵將發(fā)生部分分解.這些與對晶體中微觀結構的觀察結果是一致的.從圖7(c)可見, 在1.04 ?附近碳氫原子間的徑向分布函數(shù)gC—H(r)具有最高峰, 表明了C和H原子多數(shù)以環(huán)鏈上的C—H鍵存在.隨著沖擊速度的增加, 該位置徑向分布函數(shù)峰值降低, 表明了隨著沖擊速度的增大, C—H鍵可能會發(fā)生形變或斷裂致使部分碳氫原子間距離改變.從圖7(d)可以看出, 當沖擊速度為10 km/s時, 氫-氧原子的徑向分布函數(shù)gH—O(r)的峰值在2.45 ?附近, 該位置為RDX晶體中分子間的H和O原子距離范圍.而HONO和H2O中H—O鍵的合理范圍均為1.00 ?左右.這表明該狀態(tài)下RDX晶體中基本未發(fā)生C—H鍵斷裂, 即未出現(xiàn)H原子轉移到硝基中的氧原子上形成HONO的情況.當沖擊速度由10 km/s增大到13 km/s時, 2.45 ?附近的gH—O(r)的峰值降低, 而1.52 ?附近的徑向分布函數(shù)峰值明顯增大, 這表明了隨著沖擊速度的增大, 分子間H和O形成了分子間氫鍵, 并可能會有較少部分C-H鍵斷裂, 使H原子轉移到硝基的O原子上形成HONO, 進而使H—ONO鍵的數(shù)量增多[22,25,31,32].因此, 對于完美RDX晶體, 在沖擊加載下, 隨著沖擊速度的增加, 表征N—NO2鍵、C—N鍵和C—H鍵的原子間徑向分布函數(shù)峰值均降低, 而分子間的H和O的距離先減小形成氫鍵,之后在氫鍵作用下會有部分C—H鍵斷裂形成H—ONO.表明了RDX可能的初始分解反應為N—NO2斷裂或C—N鍵斷裂.當沖擊速度增大到一定值時(如13 km/s), 會在H和O之間形成分子間氫鍵, 并有部分C—H鍵斷裂, 進而H轉移到硝基上形成HONO.

圖7 完美RDX晶體在沖擊速度為10, 11和13 km/s作用下N—N (a), C—N (b), C—H (c)和H—O (d)原子間的徑向分布函數(shù)Fig.7.Radial distribution function between atoms N—N (a), C—N bond (b), C—H bond (c), H—O bond (d) in the prefect RDX crystal at shock velocity of 10, 11 and 13 km/s.

表1 完美晶體中原子間的徑向分布函數(shù)相關峰所包圍的面積Table 1.Area enclosed by a certain peak of inter-atoms radial distribution function in the perfect crystal.

3.3 沖擊加載下含1個分子空位的RDX晶體的初始反應分析

圖8(a)—(d)為含1個分子空位的RDX晶胞在沖擊速度10, 11和13 km/s作用下的N—N, C—N,C—H和H—O原子間的徑向分布函數(shù).表2為計算獲得的主要徑向分布函數(shù)峰所包圍的面積.顯然, 對于主要位置的N—N, C—N, C—H和H—O原子之間成鍵數(shù)量基本還是受峰值影響較大, 但對13 km/s沖擊速度情況, 在2.45 ?附近的H—O原子的徑向分布函數(shù)面積比11 km/s的略大一些,這可能是統(tǒng)計誤差所致.從圖8(a)可見, N—N原子間的徑向分布函數(shù)gN—N(r)的最高峰仍然出現(xiàn)在1.55 ?附近, 隨著沖擊速度的增加, 該處峰值逐漸降低, 而1.55 ?附近峰值逐漸升高.表明了沖擊加載下, 空位RDX晶體中也將發(fā)生N—NO2鍵斷裂, 并且隨著沖擊速度的增加, N—NO2斷裂數(shù)量增加, 而最終形成的氮氣數(shù)量也增多.從圖8(b)和圖8(c)可看出, RDX中的C—N和C—H原子間的徑向分布函數(shù)峰值也隨著沖擊速度的增加而逐漸降低, 表明了沖擊加載下, 空位RDX也可能發(fā)生C—N鍵和C—H鍵斷裂.圖8(d)表明, 隨著沖擊速度的增加, 2.45 ?附近的氫-氧原子間徑向分布函數(shù)gH—O(r)峰值同樣也有所減小, 而在1.52 ?附近的峰值明顯升高.這說明隨著沖擊速度的增加, RDX空位晶體受壓縮程度逐漸增大, 使分子間H和O更多地以分子間氫鍵形式存在, 并可能引發(fā)C—H鍵斷裂, 使氫原子轉移到NO2的氧原子上形成HONO.因此, 對于含1個分子空位的RDX晶體, 其初始分解方式基本與完美RDX分解方式相同, 即可能發(fā)生N—NO2鍵斷裂、C—N鍵斷裂和C—H鍵斷裂, 并且H會轉移到硝基中的O原子上.

表2 空位晶體中原子間的徑向分布函數(shù)峰所包圍的面積Table 2.Area enclosed by a certain peak of inter-atoms radial distribution function in the vacancy crystal.

圖8 含1個分子空位RDX超胞在沖擊速度為10, 11和13 km/s作用下原子間(N—N (a), C—N (b), C—H (c)和H—O (d))的徑向分布函數(shù)Fig.8.Radial distribution between atoms (N—N (a), C—N (b), C—H (c), H—O (d)) of the vacancy RDX supercell at 10, 11 and 13 km/s.

3.4 沖擊加載下分子空位缺陷對RDX初始反應的影響

為了明確分子空位缺陷對RDX晶體分解的影響, 圖9(a)—(d)給出了完美RDX晶體和空位RDX晶體在11 km/s沖擊速度下主要原子間的徑向分布函數(shù).表3為計算的主峰所包圍的面積.可見, 此關系仍為峰值越高面積越大.在相同沖擊速度作用下, 空位晶體中表征N—NO2鍵的氮-氮原子徑向分布函數(shù)(1.55 ?)峰值減小, 而表征氮-氮三鍵的1.18 ?位置處的徑向分布函數(shù)的峰值增大,表明了N—NO2鍵斷裂程度增大, 而形成產(chǎn)物氮氣的數(shù)量增多.可見, 分子空位的存在增加了N—NO2鍵的活性, 促進了其斷鍵反應.C—N和C—H原子間的徑向分布函數(shù)的峰值均與完美晶體的峰值相差不多.對于氫氧原子, 空位晶體在2.15 ?和1.52 ?附近的徑向分布函數(shù)峰值要比相應完美晶體的徑向分布函數(shù)峰值低, 而在1.00 ?附近, 空位晶體比完美晶體的峰值稍高, 但不明顯(面積差不到0.001).這表明, 空位缺陷對HONO的形成影響較小.因此, 在相同沖擊速度加載下, 分子空位缺陷的存在, 主要增加了N—NO2基團的反應活性, 使空位晶體更容易發(fā)生N—NO2鍵斷裂的化學反應,這與目前文獻[5]報道的缺陷的存在增加了含能化合物的反應活性的結果是一致的.

圖9 完美RDX晶體和空位RDX晶體在沖擊速度11 km/s作用下原子間(N—N (a), C—N (b), C—H (c)和H—O (d))的徑向分布函數(shù)Fig.9.Radial distribution function between atoms (N—N (a), C—N (b), C—H (c) and H—O (d)) of the perfect RDX crystal and molecular vacancy RDX crystal at shock velocity of 11 km/s.

表3 11 km/s沖擊速度加載下完美晶體和空位晶體中原子間的徑向分布函數(shù)峰所包圍的面積Table 3.Area enclosed by a certain peak of inter-atoms radial distribution function in the perfect and vacancy crystals at the shock velocity of 11 km/s.

4 結 論

本文運用ReaxFF分子動力學方法結合MSST研究了沖擊加載完美RDX晶體和含有1個分子空位的RDX晶體的動態(tài)響應過程.通過對兩種微觀結構的動力學分析, 表明了在沖擊加載下, 兩種RDX晶體的初始分解機理均為首先發(fā)生N—NO2斷裂, 然后是C—N斷裂.含有空位的晶體的初始反應時刻要早于完美晶體, 表明了空位晶體對沖擊更敏感, 更容易發(fā)生分解反應.此外, 本文計算了主要的可能參與反應的原子間的徑向分布函數(shù), 分析了不同沖擊速度及分子空位缺陷對分解反應的影響.結果表明, 隨著沖擊速度的增加, 兩種RDX晶體中N—NO2鍵及C—N鍵斷裂數(shù)目增多, 還可能出現(xiàn)C—H鍵斷裂, 并有氫轉移到硝基中的氧原子上形成HONO.在相同沖擊速度作用下, 分子空位晶體中N—NO2鍵的徑向分布函數(shù)峰值比相應的完美晶體的徑向分布函數(shù)峰值明顯降低, 表明了空位缺陷的存在增加了N—NO2的反應活性, 使其更易發(fā)生斷裂, 形成產(chǎn)物氮氣.