葉寶云，范家珂，安崇偉,3，王晶禹，3，趙鳳起，秦 釗

(1.西安近代化學研究所，陜西 西安 710065；2.中北大學 環(huán)境與安全工程學院，山西 太原 030051；3.山西省超細粉體工程技術(shù)研究中心，山西 太原 030051)

引 言

3,3′-二氨基-4，4′-氧化偶氮呋咱(3,3′-diamino-4,4′-azofurazan, DAAF)是一種呋咱類高氮雜環(huán)化合物，在常溫下的晶體密度為1.745g/cm3，生成焓為444kJ/mol，氧平衡為-52.8%，爆速為8300m/s[1]，且具有較低的機械感度和良好的熱穩(wěn)定性，其綜合性能與TATB相似，但合成過程卻比TATB經(jīng)濟環(huán)保，是一種有望替代TATB的高能鈍感含能化合物。

研究表明，炸藥的晶體形貌對其感度和能量性能影響很大。由于合成過程中結(jié)晶條件的復雜性，導致直接合成得到的DAAF晶體缺陷較多，嚴重影響了DAAF的使用性能[2-3]。因此，控制DAAF的結(jié)晶工藝、提高其晶體品質(zhì)是DAAF走向工程應用的前提。晶體的形貌不僅與晶體的內(nèi)部化學結(jié)構(gòu)有關(guān)，還取決于晶體生長的外部條件，如溶劑、過飽和度和雜質(zhì)。

目前，DAAF晶體形貌調(diào)控的方法有溶劑-非溶劑法、蒸發(fā)法。李旭陽[4]采用二甲基亞砜(20℃)/水(0℃)體系重結(jié)晶獲得平均粒徑500nm的類球形DAAF。姚笑璐[5]使用蒸發(fā)法制備DAAF聚晶,研究了溶劑極性對DAAF形貌的影響。結(jié)果表明，在弱極性溶劑(乙酸乙酯、丙酮)中DAAF晶形呈長片狀。解佳妮[6-7]以DMSO為溶劑，水為非溶劑，在溶劑、非溶劑比為1∶40、攪拌速度為150r/min條件下，得到球形DAAF。當前研究主要集中在含能材料領(lǐng)域中常見溶劑對DAAF形貌的調(diào)控，而探究新溶劑對DAAF形貌的影響則相對較少。

隨著計算機科學的發(fā)展，量子化學、分子動力學及蒙特卡羅等方法被廣泛應用于計算材料的相關(guān)性質(zhì)。其中，修正的附著能模型被廣泛應用于含能材料等材料的晶體形貌模擬[8-9]。目前，對DAAF在不同溶劑中晶體形貌的模擬還未被系統(tǒng)研究。

為探究DAAF在不同溶劑中的晶體形貌，同時也為擴大溶劑選擇范圍，基于COSMO-RS理論使用COSMOquick軟件對DAAF的溶劑進行篩選。而后基于分子動力學模擬，使用修正附著能模型，計算了DAAF在二甲基亞砜、1, 2-二氨基乙烷、六甲基磷酰胺、N, N-二甲基乙酰胺、N, N-二甲基甲酰胺、N-甲基吡咯烷酮、哌啶7種單溶劑(由COSMOquick篩選得到)中的形貌。通過理論計算的方式探究溶劑對DAAF晶習的影響，以期對DAAF結(jié)晶的實驗研究提供理論指導。

1 理論與計算方法

1.1 修正COSMO-RS理論

COSMO-RS(Conductor like Screening Model for Real Solvents)理論于1993年由Klamt等[10]提出并發(fā)表。該理論在COSMOthermX代碼中實現(xiàn)。在化學、制藥等領(lǐng)域的應用中，COSMOthermX被證明在為給定的溶質(zhì)篩選溶劑時，具有實際應用價值?；诒砻骐姾擅芏确植?σ-profile)文件或其他包含COSMO表面電荷信息的文件，COSMO-RS可以輕松計算出溶質(zhì)X在純?nèi)軇┗蚧旌先軇﹊中的化學勢。基于計算得到的化學勢，摩爾分數(shù)溶解度xi可以根據(jù)式(1)計算得到：

(1)

因為準確預測熔融吉布斯自由能需要對晶體的結(jié)構(gòu)進行預測，所以熔融吉布斯自由能通常無法直接計算得到。基于篩選溶劑的目的，Klamt[11]提出一種基于多種參考溶劑來預測熔融吉布斯自由能的方法：

(2)

對于每一種新溶劑j，其校正項可以通過權(quán)重因子cij對溶劑j與每一個參考溶劑i的相似性進行加權(quán)，進而計算出溶劑j的化學勢校正項。

(3)

(4)

(5)

(6)

最終，可以計算得到溶質(zhì)X在溶劑j中的溶解度：

(7)

1.2 MAE模型

附著能模型(Attachment Energy Model, AE模型)中，附著能(Eatt)的定義為厚度為dhkl的晶片附著到(hkl)面上釋放的能量，而晶體生長出一層厚度為dhkl的晶片所釋放的能量為晶片能。

Eatt=Elatt-Eslice

(8)

式中：Eatt為附著能；Elatt為晶格能；Eslice為晶片能。

當考慮溶劑與晶面層之間的相互作用時，可得到修正的附著能計算公式:

(9)

式中：Aacc為晶面單元(hkl)上溶劑可及面積；Ahkl為晶面單元(hkl)的橫截面積。Es代表溶劑對晶面的影響能力，其定義為：

(10)

式中：Abox為模擬盒子的橫截面積；Eint為晶面層與溶劑之間的相互作用能。

Eint計算公式見式(11):

Eint=Etot-(Ecry+Esol)

(11)

式中：Eint為晶面層與溶劑之間的相互作用能；Etot為晶面與溶劑混合體系的總能量；Ecry為計算模型中晶面的能量；Esol為計算模型中溶劑的能量。同時，溶劑條件下與真空下一致，晶體的相對生長速率與(修正后)附著能成正比關(guān)系，即：

(12)

1.3 計算方法

以Qiu[12]文章中給出的298.15K下DAAF在水、丙酮、二氯甲烷、乙酸乙酯、甲醇、乙腈6種溶劑中溶解度作為參考，使用COSMOquick2020溶解度預測功能，計算了DAAF在軟件內(nèi)114種常用溶劑中的溶解度，根據(jù)溶解度大小，篩選出了7種溶解度較大的溶劑，分別為二甲基亞砜、1, 2-二氨基乙烷、六甲基磷酰胺、N, N-二甲基乙酰胺、N, N-二甲基甲酰胺、N-甲基吡咯烷酮和哌啶。上述7種物質(zhì)的分子式、簡化分子線性輸入規(guī)范(Simplified molecular input line entry specification，SMILES)、相對分子質(zhì)量及DAAF在其中的溶解度見表1。

表1 7種溶劑的分子式、SMILES、相對分子質(zhì)量和溶解度

計算使用的DAAF晶體的單胞結(jié)構(gòu)來自劍橋晶體數(shù)據(jù)庫。DAAF的初始晶胞參數(shù)a=9.3212?,b=9.6326?,c=8.9004?,α=γ=90°,β=91.3434°,空間群為P21/c，屬于單斜晶系。DAAF的分子結(jié)構(gòu)及單胞結(jié)構(gòu)如圖1所示。

圖1 DAAF分子結(jié)構(gòu)與單胞結(jié)構(gòu)

COMPASSII[13]力場作為一種高階版本的COMPASS力場，不僅增加了更多的力場類型，還可以更好地模擬凝聚態(tài)體系，且COMPASSII力場對于DAAF具有較好的適用性[14-15]。因此后續(xù)的經(jīng)典分子動力學計算過程中，計算使用的力場選擇COMPASSII力場。后續(xù)的計算均在Materials Studio 2020中完成。

首先，對DAAF單胞結(jié)構(gòu)進行優(yōu)化，得到DAAF單胞結(jié)構(gòu)的穩(wěn)定狀態(tài)，優(yōu)化算法選擇軟件內(nèi)置的Smart算法。范德華力求和方法選擇Atom based，截斷半徑設為12.5?；靜電力求和方法選擇Ewald求和方法，Ewald精度設為0.042kJ/mol；電荷分配方式選擇Forcefield assigned。使用Morphology模塊對DAAF晶體在真空條件下的形貌進行了計算，獲得DAAF在真空條件下的重要生長面。然后，對優(yōu)化的DAAF單胞結(jié)構(gòu)進行切面處理，將相應的重要生長面分別擴展為5×2×4、3×3×4和5×2×4的超晶胞結(jié)構(gòu)。使用Amorphous Cell模塊構(gòu)建包含300個溶劑分子的溶劑層。使用Build Layers功能，構(gòu)建“炸藥/溶劑”雙層模型，溶劑層上方設置50?的真空層。以DMSO溶劑為例，DAAF/DMSO雙層模型如圖2所示。對雙層模型進行結(jié)構(gòu)優(yōu)化，然后進行分子動力學計算。進行分子動力學計算時，選擇NVT系綜，模擬時間設定為300ps，步長為1fs，共執(zhí)行300000步計算，每1000步輸出一次結(jié)構(gòu)信息，控溫方式采用Andersen控溫方法，粒子的初始速度按照高斯分布進行隨機分配。最后，對計算得到的軌跡文件的穩(wěn)定部分，使用Perl腳本對相應的炸藥與溶劑之間的相互作用能進行統(tǒng)計平均。

圖2 DAAF/DMSO雙層模型

為詳細分析DAAF與溶劑分子之間的相互作用，基于密度泛函理論，使用B3LYP-D3方法，def-TZVP基組對DAAF分子以及7種溶劑分子進行結(jié)構(gòu)優(yōu)化及振動分析。而后，對優(yōu)化得到的結(jié)構(gòu)進行分子表面靜電勢計算，同時，參照盧天[16]對分子極性指數(shù)(molecular polarity index, MPI)的定義，使用Multiwfn[17]計算了7種溶劑與DAAF的MPI指數(shù)。根據(jù)分子表面靜電勢計算得到的結(jié)果，基于分子動力學計算得到的分子運動軌跡，計算了DAAF氨基中的氮原子與溶劑分子給電子基團中的主要原子間的徑向分布函數(shù)。

2 結(jié)果與討論

2.1 真空條件下的形貌

基于附著能模型(AE模型)計算得到DAAF在真空中的形貌，見圖3。真空條件下，DAAF晶體的重要生長面及面積占比見表2。根據(jù)圖3可知，真空條件下，DAAF的形貌呈立方體狀，晶體長徑比為2.094。DAAF的重要生長面為(0 1 1)、(1 0 0)和(1 1 0)。由表2可知，DAAF在真空條件下面積最大的為(0 1 1)晶面，其次是(1 0 0)晶面，最小的是(1 1 0)晶面。由于AE模型中附著能絕對值的大小與對應晶面的相對生長速率成正比，DAAF各個晶面的相對生長速率由大到小的順序為：(1 1 0)>(0 1 1)>(1 0 0)。

圖3 DAAF真空下的形貌(AE模型)

表2 DAAF真空下晶習參數(shù)(AE模型)

在附著能模型中，晶面單元(hkl)的溶劑可及面積與晶面單元橫截面積之比表示晶面(hkl)的表面粗糙度[18]，即：

(13)

通常，S數(shù)值越大，對應的(hkl)晶面表面越粗糙，越有利于溶劑的吸附。通過計算可知，(0 1 1)、(1 0 0)和(1 1 0)晶面對應的S數(shù)值分別為2.130、1.338和1.750。(0 1 1)晶面粗糙度最大，說明(0 1 1)晶面最有利于溶劑吸附，這將抑制(0 1 1)晶面上溶質(zhì)的吸附，不利于該晶面的生長。圖4中給出了DAAF的3個重要生長面(0 1 1)、(1 0 0)及(1 1 0)晶面的堆積結(jié)構(gòu)以及相應的Connolly表面。

圖4 DAAF各重要生長面的分子堆積結(jié)構(gòu)(左)及對應Connolly表面(右)

根據(jù)圖4可知，DAAF在3個生長面上均有4個DAAF分子裸露。3個晶面中DAAF裸露部分由多到少順序為：(0 1 1)>(1 1 0)>(1 0 0)。這與計算得到的S的大小順序基本一致。

2.2 溶劑中的形貌

2.2.1 分子表面靜電勢

分子表面靜電勢的電子密度等值面圖可以直觀顯示分子表面的電荷分布，通過其表面的顏色可以判斷分子中指定區(qū)域的靜電勢?；贐3LYP-D3(BJ)/def-TZVP計算得到了DAAF分子及7種溶劑分子的表面靜電勢圖，如圖5所示。

圖5 DAAF與各溶劑分子的分子表面靜電勢：(a)3,3′-二氨基-4.4′-氧化偶氮呋咱;(b)二甲基亞砜;(c)1,2-二氨基乙烷;(d)N-甲基吡咯烷酮;(e)哌啶;(f)N,N-二甲基乙酰胺;(g)N,N-二甲基甲酰胺;(h)六甲基磷酰胺

紅色區(qū)域表示靜電勢為負，說明該區(qū)域具有給電子特性；藍色區(qū)域靜電勢為正，說明該區(qū)域具有吸電子特性。根據(jù)圖5可知，DAAF的氨基附近的靜電勢為正，具有吸電子特性；二甲基亞砜中氧原子附近的靜電勢為負，具有給電子特性；1,2-二氨基乙烷中氮原子附近的靜電勢數(shù)值為負，具有給電子特性，而與氮原子相連的氫原子附近的靜電勢為正，說明該區(qū)域具有較弱的吸電子特性；N-甲基吡咯烷酮、N,N-二甲基乙酰胺、N,N-二甲基甲酰胺、六甲基磷酰胺中具有給電子特性的原子均為氧原子，其附近的靜電勢均為負；哌啶中氮原子附近的靜電勢為負值，具有一定的吸電子特性，與氮原子相連的氫原子的給電子特性相對較弱。因此，DAAF中氨基的氫原子與溶劑中具有給電子特性的氧原子或氮原子之間具有較強的吸引力，具有形成分子間氫鍵的潛力。

2.2.2 相互作用能

溶劑與溶質(zhì)之間的相互作用越強，對相應晶面生長的抑制作用越強，導致該晶面面積占比越小。由于不同晶面的擴展倍數(shù)不同，DAAF超晶胞內(nèi)DAAF分子的數(shù)目不一致。在此對計算得到的相互作用能進行歸一化處理。計算公式見式(14)：

(14)

式中：E為歸一化之后計算得到的能量；Eint為計算得到的初始相互作用能；xa為超晶胞在a方向的擴展倍數(shù)；xb為超晶胞在b方向的擴展倍數(shù)。

DAAF的3個晶面與不同溶劑間的相互作用如表3所示。

表3 不同溶劑中DAAF與溶劑的相互作用能

根據(jù)表3可知，歸一化處理之后，7種溶劑與DAAF之間的相互作用能均為負值，則對應的結(jié)合能為正值，說明兩者之間的相互作用以吸引作用為主。DAAF的3個晶面與溶劑間的結(jié)合能由強到弱均為(1 1 0)>(0 1 1)>(1 0 0)。7種溶劑中，1,2-二氨基乙烷與DAAF之間的結(jié)合作用最強，二甲基亞砜與DAAF之間的結(jié)合作用最弱。在(0 1 1)面上，7種溶劑與DAAF之間的結(jié)合作用由強到弱順序為：1,2-二氨基乙烷>N,N-二甲基甲酰胺>N-甲基吡咯烷酮>N,N-二甲基乙酰胺>哌啶>六甲基磷酰胺>二甲基亞砜；在(1 0 0)面上，7種溶劑與DAAF之間的結(jié)合作用由強到弱的順序為：1,2-二氨基乙烷>N-甲基吡咯烷酮>哌啶>N,N-二甲基甲酰胺>N,N-二甲基乙酰胺>六甲基磷酰胺>二甲基亞砜；在(1 1 0)面上，DAAF與7種溶劑分子之間的結(jié)合作用由強到弱的順序為：1,2-二氨基乙烷>N,N-二甲基甲酰胺>N,N-二甲基乙酰胺>N-甲基吡咯烷酮>六甲基磷酰胺>哌啶>二甲基亞砜。7種溶劑在DAAF不同晶面上相互作用能的差異會影響DAAF在不同溶劑中的形貌。

2.2.3 徑向分布函數(shù)

徑向分布函數(shù)表示指定原子A的距離r處B原子出現(xiàn)的概率，可以用來對原子對之間的非鍵相互作用進行分析。根據(jù)DAAF及7種溶劑的分子表面靜電勢分析，對DAAF氨基中氮原子與溶劑中給電子原子之間的徑向分布函數(shù)進行了計算，結(jié)果如圖6所示。一般說來，非鍵相互作用包括氫鍵(小于3.1?)，范德華力(3.1～5.0?)以及靜電力(大于5.0?)[19]。

圖6 DAAF氨基的氮原子與溶劑中給電子原子之間的徑向分布函數(shù)

以二甲基亞砜為例，紅色曲線對應的是DAAF(0 1 1)面上二甲基亞砜中氧原子與DAAF氨基中氮原子之間的徑向分布函數(shù)，其在r為3.1?處出現(xiàn)第一個峰值，在r=4.21?處出現(xiàn)第二個峰值，在r=3.1～5.0?(不包含3.1)范圍內(nèi)出現(xiàn)兩個峰值；類似的，(1 0 0)和(1 1 0)晶面在3.1～5.0?范圍內(nèi)均有峰值出現(xiàn)。這表明，在DAAF的3個重要生長面上，DAAF與二甲基亞砜之間存在范德華相互作用，沒有氫鍵相互作用。根據(jù)徑向分布函數(shù)峰值位置及數(shù)量可知，在DAAF的(1 1 0)面上，DAAF與二甲基亞砜之間的相互作用最強，其次是(0 1 1)，最后是(1 0 0)面。徑向分布函數(shù)的分析結(jié)果與相互作用能的分析結(jié)果一致。

根據(jù)圖6可知，DAAF與1，2-二氨基乙烷在3個重要晶面上徑向分布函數(shù)第一個峰值出現(xiàn)的位置均在r=3.21?之后，故DAAF與1，2-二氨基乙烷之間均只有范德華相互作用，無氫鍵相互作用；類似的，在3個重要晶面上，DAAF與六甲基磷酰胺、N,N-二甲基乙酰胺、N,N-二甲基甲酰胺、N-甲基吡咯烷酮、哌啶之間的徑向分布函數(shù)的第一個峰值出現(xiàn)的位置均在r=3.1?之前，故DAAF與這5種溶劑分子之間氫鍵相互作用與范德華相互作用同時存在。

2.2.4 溶劑中DAAF的晶體形貌與溶劑極性

根據(jù)修正附著能模型計算得到的DAAF不同生長晶面與不同溶劑之間的修正附著能及縱橫比見表4。計算得到的相應溶劑條件下的DAAF晶體形貌以及298K下DAAF在二甲基亞砜中的重結(jié)晶形貌見圖7。從圖7中可以看出，在二甲基亞砜中實驗重結(jié)晶得到的DAAF形貌與計算得到的DAAF在DMSO中的形貌基本一致，均近似棱柱狀。

表4 DAAF在溶劑作用下的修正附著能及相關(guān)參數(shù)

圖7 不同溶劑條件下DAAF的晶體形貌及298K下DAAF在二甲基亞砜中的形貌

根據(jù)表4計算結(jié)果可知，溶劑為二甲基亞砜時，DAAF的3個重要晶面顯露面積最大的為(1 1 0)面，其次是(0 1 1)和(1 0 0)。溶劑為1,2-二氨基乙烷、六甲基磷酰胺、N,N-二甲基乙酰胺、N,N-二甲基甲酰胺、N-甲基吡咯烷酮時，DAAF的3個重要晶面面積占比由大到小的順序均為：(1 0 0)>(1 1 0)>(0 1 1)，形貌呈厚度不一的六棱片狀；溶劑為哌啶時，DAAF的3個重要生長面中，(1 0 0)面消失，(1 1 0)面的面積最大達到99.9%，(0 1 1)面相對面積占比僅占0.1%，其長徑比達到696.366，形貌呈針狀。DAAF在上述7種溶劑中的晶體形貌與真空條件下計算得到的DAAF晶體形貌差距較大，表明溶劑對晶體形貌影響較大。

計算得到的分子極性指數(shù)結(jié)果為：DAAF的分子極性指數(shù)為66.19kJ/mol；二甲基亞砜的分子極性指數(shù)最大，為78.13kJ/mol；哌啶的分子極性指數(shù)最小，為29.06kJ/mol；1，2-二氨基乙烷、六甲基磷酰胺、N，N-二甲基乙酰胺、N，N-二甲基甲酰胺、N-甲基吡咯烷酮的分子極性指數(shù)分別為：56.35、41.83、58.53、67.11、59.24kJ/mol。7種溶劑中，僅二甲基亞砜及N, 5N-二甲基甲酰胺的極性比DAAF的極性強。7種溶劑的極性由大到小順序為：二甲基亞砜>N,N-二甲基甲酰胺>N-甲基吡咯烷酮>N,N-二甲基乙酰胺>1,2-二氨基乙烷>六甲基磷酰胺>哌啶。

DAAF作為一種極性較強的分子，通過結(jié)晶的方法使用單溶劑對DAAF晶體形貌進行改性時，使用的溶劑的極性對DAAF的形貌有較大的影響。當溶劑極性小于DAAF且兩者的分子極性指數(shù)差值大于33.49kJ/mol時，得到的DAAF的晶習呈針狀；當使用的一元溶劑的極性小于DAAF的極性時，DAAF的晶習呈六棱片狀；當使用的一元溶劑的極性大于DAAF且兩者分子極性指數(shù)差值接近12.56kJ/mol時，DAAF的晶習接近棱柱狀。通過計算DAAF在一元溶劑中的晶體形貌可知，通過結(jié)晶方法改善DAAF的晶體形貌時，不能僅考慮DAAF在溶劑中的溶解度還需要考慮溶劑分子與DAAF之間的分子極性指數(shù)的差值。

3 結(jié) 論

(1)DAAF的晶體形貌在真空中呈塊狀，3個主要顯露面分別為：(0 1 1)、(1 0 0)以及(1 1 0)面。其中，(0 1 1)面的面積占比最大，為62.276%；(1 1 0)面的面積占比最小，為0.041%。

(2)7種溶劑與DAAF之間的結(jié)合能均為正值，7種溶劑中1,2-二氨基乙烷與DAAF之間的結(jié)合作用最強，二甲基亞砜與DAAF之間的結(jié)合作用最弱。以二甲基亞砜為例，DAAF與其之間的相互作用以范德華力占主導。

(3)分子表面靜電勢與徑向分布函數(shù)分析表明，DAAF中氨基中的氮原子與二甲基亞砜和1，2-二氨基乙烷中的給電子原子之間僅存在范德華相互作用；DAAF氨基中的氮原子與其余5種溶劑中的給電子原子之間存在氫鍵及范德華相互作用。DAAF(1 0 0)晶面表面具有較多的氨基，有利于與溶劑分子內(nèi)的給電子原子之間形成氫鍵及范德華相互作用。

(4)不同極性的溶劑對DAAF的晶體形貌具有較大的影響。溶劑的極性小于DAAF的極性時，DAAF在其中的形貌主要呈現(xiàn)六棱片狀；當溶劑的極性遠小于DAAF的極性時，DAAF的形貌呈針狀；當溶劑的極性大于DAAF且兩者分子極性指數(shù)的差值接近12.56kJ/mol時，DAAF的晶習接近棱柱狀。