含能材料晶形預(yù)測(cè)方法：附著能模型及其發(fā)展

劉英哲

(西安近代化學(xué)研究所，陜西 西安 710065)

引 言

晶形控制是目前含能材料改性研究的重要途徑之一[1]。對(duì)于相同的含能晶體，晶形能夠在很大程度上影響感度、能量、流散性等性質(zhì)[2-3]。例如，在相同尺寸下，球形晶體通常比其他晶形具有更高的裝填密度和更低的機(jī)械感度[4-5]。然而，含能材料的晶形控制仍然較為困難，因?yàn)榫紊L(zhǎng)不僅與晶體結(jié)構(gòu)有關(guān)，還受到溶劑、溫度、添加劑、結(jié)晶工藝等多種條件的影響。溶液結(jié)晶是含能材料最常用的晶形生長(zhǎng)方式，但采用不同的溶劑時(shí)，得到的晶形差異很大[6]。

與基于“試錯(cuò)法”的結(jié)晶實(shí)驗(yàn)相比，通過(guò)理論方法預(yù)測(cè)晶形，有助于從微觀上深入理解晶體生長(zhǎng)機(jī)制，能夠快速篩選溶劑、溫度等結(jié)晶條件，從而為含能材料的晶形控制提供理論指導(dǎo)。當(dāng)前，含能材料的晶形預(yù)測(cè)方法主要包括：Bravais-Friedel-Donnay-Harker (BFDH) 準(zhǔn)則[7]、附著能模型[8-9]、占有率模型[10-11]、蒙特卡洛模擬方法[12-14]、界面結(jié)構(gòu)分析模型[15]、螺旋生長(zhǎng)模型和二維成核模型[16-19]等，這些方法的原理逐漸從結(jié)構(gòu)因素、能量因素向生長(zhǎng)機(jī)制演化，預(yù)測(cè)精度也逐漸增加。其中，附著能模型在考慮晶體結(jié)構(gòu)因素的基礎(chǔ)上，又引入了能量因素，預(yù)測(cè)結(jié)果相對(duì)準(zhǔn)確，計(jì)算簡(jiǎn)單，又有很好的軟件支持，成為含能材料領(lǐng)域中應(yīng)用最廣泛的晶形預(yù)測(cè)方法。

為此，本文綜述了國(guó)內(nèi)外近十年來(lái)附著能模型在含能材料中的應(yīng)用情況，詳細(xì)闡述了附著能的理論基礎(chǔ)與計(jì)算方法，討論了影響晶形預(yù)測(cè)結(jié)果的重要因素，并介紹了附著能模型的最新研究進(jìn)展。

1 附著能模型理論與方法

1.1 基礎(chǔ)理論

附著能模型基于恒定生長(zhǎng)速率的簡(jiǎn)單假定[8-9]，即晶體中不同晶面的相對(duì)生長(zhǎng)速率(Rhkl)與其附著能(Eatt)的絕對(duì)值成正比關(guān)系：

Rhkl∝|Eatt|

(1)

晶面的附著能絕對(duì)值越大，其相對(duì)生長(zhǎng)速率越快，在生長(zhǎng)過(guò)程中越容易消失。附著能定義為一層厚度為dhkl的晶片附著在晶面(h k l)上所釋放的能量，可由式(2)進(jìn)行計(jì)算：

Eatt=Elatt-Eslice

(2)

式中：Elatt為晶格能量；Eslice為晶片能量。對(duì)于穩(wěn)定生長(zhǎng)的晶面，附著能為負(fù)值。需要注意的是，不同晶面的Eatt需要在相同的化學(xué)計(jì)量下(如每分子、每晶胞)比較才有意義，目前在絕大部分文獻(xiàn)中，Eatt的單位簡(jiǎn)寫(xiě)為kcal/mol，實(shí)際上完整的單位應(yīng)為kcal/(mol·molecule)或kcal/(mol·unitcell)。

基于給定的晶體結(jié)構(gòu)，通過(guò)計(jì)算不同生長(zhǎng)晶面的附著能，可預(yù)測(cè)真空條件下的生長(zhǎng)形貌。在溶劑條件中，溶劑將附著在生長(zhǎng)晶面上，溶質(zhì)在晶面上生長(zhǎng)時(shí)要克服溶劑層的阻力，導(dǎo)致晶面生長(zhǎng)速率受到抑制。為此，需要對(duì)附著能進(jìn)行修正：

(3)

式中：Es為溶劑吸附所引入的能量修正項(xiàng)，表示溶劑對(duì)晶面生長(zhǎng)的抑制作用；S為修正因子，體現(xiàn)了生長(zhǎng)晶面的結(jié)構(gòu)特征，可理解為不同晶面結(jié)構(gòu)對(duì)Es的權(quán)重，其定義為:

(4)

式中：Aacc為晶面單元(h k l)的溶劑可及表面積；Ahkl為晶面單元的橫截面積。顯然，S值越大，晶面越粗糙；S值越小，晶面越平坦。Es可通過(guò)溶劑與晶面的相互作用能(Eint)進(jìn)行計(jì)算：

(5)

式中：Abox為模擬盒子中生長(zhǎng)晶面(h k l)的橫截面積；Eint為晶面-溶劑體系的總能量(Etot)與晶面能量(Ecry)和溶劑能量(Esol)之差，即：

Eint=Etot-(Ecry+Esol)

(6)

在溶劑條件下，假定晶面的相對(duì)生長(zhǎng)速率與修正后的附著能仍然成正比關(guān)系，即：

(7)

因此，通過(guò)計(jì)算溶劑中的修正附著能，即可預(yù)測(cè)晶體在溶劑中的生長(zhǎng)形貌。為了便于理解溶劑條件下附著能公式中的相關(guān)計(jì)算參數(shù)，將其物理含義形象地表示出來(lái)，如圖1所示。另外，基于類似的思想，也可將附著能模型進(jìn)一步拓展至添加劑條件下[20-22]，此處不再贅述。

圖1 溶劑中附著能計(jì)算公式示意圖Fig.1 Schematic diagram of attachment energy calculation in solvent

1.2 計(jì)算模型

為了預(yù)測(cè)溶劑中的生長(zhǎng)形貌，需要構(gòu)建合理的晶面-溶劑界面模型，從而準(zhǔn)確地計(jì)算附著能。圖2為晶面-溶劑界面模型的構(gòu)建示意圖。首先，基于給定的晶胞結(jié)構(gòu)，計(jì)算低米勒指數(shù)晶面的附著能，得到真空中的預(yù)測(cè)晶形和穩(wěn)定生長(zhǎng)晶面。然后，沿著每個(gè)生長(zhǎng)晶面的(hkl)指數(shù)，對(duì)晶胞進(jìn)行切割獲得生長(zhǎng)晶面單元，并構(gòu)建超晶胞，從而獲得三維的晶面模型。對(duì)于溶劑模型，可根據(jù)給定密度將溶劑分子隨機(jī)地填充到三維盒子中。最后，沿著生長(zhǎng)晶面的法向量方向，將溶劑模型“拼接”到晶面模型上方，形成晶面-溶劑界面的初始結(jié)構(gòu)。在添加劑條件下，也可通過(guò)類似方法構(gòu)建計(jì)算模型，但一般用單個(gè)添加劑分子替代溶劑層結(jié)構(gòu)[23]。

圖2 晶面-溶劑界面模型搭建示意圖 Fig.2 Scheme in construction of crystal-solvent interfacial model

隨后，在NVT系綜中對(duì)上述界面結(jié)構(gòu)進(jìn)行分子動(dòng)力學(xué)模擬。由于施加了三維周期性邊界條件，為了盡可能地避免溶劑與晶面的鏡像發(fā)生相互作用，還需要在溶劑層上方設(shè)置一定厚度的真空層。真空層的厚度要大于非鍵作用的截?cái)喟霃剑ǔ閹资?。另外，附著能模型假定生長(zhǎng)晶面上不發(fā)生結(jié)構(gòu)重排[8]，一般可在分子動(dòng)力學(xué)模擬中固定晶面的位置。待界面體系到達(dá)平衡后，對(duì)模擬軌跡進(jìn)行平均求得晶面-溶劑相互作用能，從而進(jìn)一步計(jì)算溶劑中的附著能。

2 附著能模型在含能材料中的應(yīng)用

2.1 在晶形預(yù)測(cè)中的應(yīng)用

當(dāng)前，附著能模型已經(jīng)廣泛應(yīng)用于含能材料的晶形預(yù)測(cè)，研究體系包括常見(jiàn)的分子晶體、離子晶體和含能共晶，分別見(jiàn)表1、表2和表3。這些研究預(yù)測(cè)了含能材料在真空、(混合)溶劑、添加劑等不同生長(zhǎng)條件的晶形，并討論了溫度效應(yīng)對(duì)晶形生長(zhǎng)的影響。此外，將附著能結(jié)果進(jìn)一步結(jié)合徑向分布函數(shù)、擴(kuò)散系數(shù)等計(jì)算分析方法，從微觀上揭示了溶劑與生長(zhǎng)晶面的相互作用，為含能材料的晶形控制，甚至是共晶形成提供了理論支撐[24-26]。

表1 含能分子晶體的晶形預(yù)測(cè)文獻(xiàn)匯總Table 1 Summary on the morphology predictions of energetic molecular crystals

表2 含能離子晶體的晶形預(yù)測(cè)文獻(xiàn)匯總Table 2 Summary on the morphology predictions of energetic ionic crystals

表3 含能共晶的晶形預(yù)測(cè)文獻(xiàn)匯總Table 3 Summary on the morphology predictions of energetic cocrystals

2.2 不同計(jì)算結(jié)果間的差異

對(duì)于大部分含能晶體，附著能模型能夠從定性上給出與實(shí)驗(yàn)數(shù)據(jù)相吻合的預(yù)測(cè)結(jié)果。然而，從相同體系的預(yù)測(cè)晶形預(yù)測(cè)對(duì)比圖(見(jiàn)圖3)可以看出，有時(shí)附著模型預(yù)測(cè)的晶形也存在一定差異。盡管附著能模型的理論假設(shè)非常明晰，但附著能的計(jì)算結(jié)果受很多因素影響，如力場(chǎng)參數(shù)，模型尺寸，晶面-溶劑相互作用能等，這都將直接決定晶形的預(yù)測(cè)結(jié)果。因此，如何準(zhǔn)確地計(jì)算附著能是一個(gè)關(guān)鍵問(wèn)題。

圖3 相同體系下的晶形預(yù)測(cè)結(jié)果對(duì)比Fig.3 Comparison on the morphology predictions of same systems

3 附著能的準(zhǔn)確計(jì)算

3.1 Eatt的計(jì)算

從附著能的定義上看，附著能的本質(zhì)是相互作用能。因此，在實(shí)際計(jì)算過(guò)程中，真空中的附著能可通過(guò)式(2)求得：

(8)

為了準(zhǔn)確計(jì)算附著能，需要考慮其隨著晶面厚度的收斂情況。如圖4所示，可計(jì)算一系列不同厚度的晶面能量，當(dāng)其能量差趨于穩(wěn)定時(shí)，則附著能達(dá)到收斂。此處以β-HMX晶體為例，晶胞的Z值為2，5個(gè)生長(zhǎng)晶面(0 1 1)、(1 1 0)、(0 2 0)、(1 0 1)和(1 0 -1)的Z值分別為2、2、1、2和2。

圖4 真空中附著能計(jì)算示意圖 Fig.4 Schematic diagram of attachment energy calculation in vacuum

表4給出了真空中HMX生長(zhǎng)晶面的附著能隨晶面厚度的變化情況[31]。不同晶面的附著能對(duì)晶面厚度的收斂情況是不同的。顯然，當(dāng)晶層厚度不足時(shí)，附著能計(jì)算偏差較大。對(duì)于HMX晶體，當(dāng)n≥6時(shí)，所有晶面的附著能都已收斂，與Materials Studio軟件[68]給出的結(jié)果基本一致。另外，考察附著能隨晶面厚度的收斂情況，也可為界面模型中晶面層厚度的構(gòu)建提供參考。

表4 真空中HMX生長(zhǎng)晶面的附著能隨晶面厚度的變化情況Table 4 Attachment energies of HMX growth faces in a vacuum as a function of crystal layer thickness (unit: kcal/mol/unitcell)

3.2 力場(chǎng)參數(shù)

附著能的計(jì)算是基于分子力場(chǎng)完成的。因此，在開(kāi)始計(jì)算附著能之前，要對(duì)力場(chǎng)參數(shù)的可靠性進(jìn)行評(píng)價(jià)。選擇合適的分子力場(chǎng)是準(zhǔn)確計(jì)算附著能的基礎(chǔ)。通常，可將擬選用的分子力場(chǎng)用于晶體結(jié)構(gòu)優(yōu)化，根據(jù)優(yōu)化后的晶格參數(shù)與實(shí)驗(yàn)數(shù)據(jù)的偏差，對(duì)力場(chǎng)的準(zhǔn)確性進(jìn)行驗(yàn)證。此外，也可基于分子力場(chǎng)計(jì)算晶格能，與其實(shí)驗(yàn)值進(jìn)行對(duì)比。晶格能可通過(guò)實(shí)驗(yàn)測(cè)定的升華焓(ΔHsub)數(shù)據(jù)進(jìn)行求解[8,29]：

Elatt=-ΔHsub-2RT

(9)

式中：R和T分別為氣體常數(shù)和溫度；2RT為校正因子。晶格參數(shù)與晶格能考察的主要是分子間作用參數(shù)，同時(shí)也要注意晶體優(yōu)化后分子結(jié)構(gòu)是否發(fā)生較大變化。另外，可通過(guò)溶劑密度計(jì)算來(lái)驗(yàn)證分子力場(chǎng)對(duì)溶劑的兼容性。

COMPASS力場(chǎng)[69]是目前含能材料領(lǐng)域中應(yīng)用最廣泛的力場(chǎng)，可適用于絕大部分CHON含能分子。然而，隨著近年來(lái)新型含能化合物的陸續(xù)合成，有一些特殊的成鍵類型，難以用COMPASS力場(chǎng)描述，比如TKX-50中的N→O鍵。表5對(duì)比了COMPASS、PCFF[70-72]、CVFF[73]3種力場(chǎng)對(duì)TKX-50晶格參數(shù)的優(yōu)化結(jié)果[23]。從表5可知，COMPASS力場(chǎng)優(yōu)化的晶格參數(shù)與實(shí)驗(yàn)數(shù)據(jù)相差較大。CVFF力場(chǎng)優(yōu)化的b軸長(zhǎng)度與實(shí)驗(yàn)值也存在較大偏差。相比之下，只有PCFF力場(chǎng)的優(yōu)化結(jié)果與實(shí)驗(yàn)結(jié)果較為接近，因此被用于描述TKX-50晶體[21,59,60]。

表5 不同力場(chǎng)優(yōu)化下的TKX-50晶格參數(shù)Table 5 Lattice parameters of TKX-50 crystal optimized by different force fields

需要注意的是，TKX-50是一種含能離子鹽，正負(fù)離子之間存在很強(qiáng)的靜電相互作用，對(duì)生長(zhǎng)晶面的附著能影響很大。在PCFF力場(chǎng)中，TKX-50的正、負(fù)離子分別帶有+1e和-2e的電量，與量化計(jì)算結(jié)果+0.82e和-1.64e不相符[74]。因此，可通過(guò)密度泛函理論GGA-PBE/DNP對(duì)TKX-50晶胞進(jìn)行電荷布居分析，將Mulliken電荷作為原子電荷，進(jìn)一步修正TKX-50晶體內(nèi)的靜電作用。經(jīng)該參數(shù)優(yōu)化后，TKX-50的晶格參數(shù)變?yōu)閍=5.45?,b=10.89?,c=6.23?,α=90.00°,β= 8.44°,γ=90.00°，與實(shí)驗(yàn)數(shù)據(jù)也較為相近。

表6匯總了真空下TKX-50的晶形預(yù)測(cè)結(jié)果。同樣采用PCFF力場(chǎng)描述TKX-50時(shí)，晶形預(yù)測(cè)結(jié)果也略有不同，但主要生長(zhǎng)晶面和長(zhǎng)徑比基本相同。對(duì)PCFF力場(chǎng)的電荷參數(shù)修正后，晶形與生長(zhǎng)晶面發(fā)生了較大變化。另外，COMPASS力場(chǎng)計(jì)算的晶形類似于紡錘狀。由此可見(jiàn)，晶形預(yù)測(cè)結(jié)果對(duì)力場(chǎng)參數(shù)較為敏感。

除了TKX-50，劉寧等[44]對(duì)FOX-7的Dreiding力場(chǎng)參數(shù)[75]進(jìn)行了部分修正，使晶格參數(shù)的相對(duì)計(jì)算誤差在5%之內(nèi)；Song等[40]結(jié)合COMPASS力場(chǎng)和Gasteiger電荷，進(jìn)一步將晶格參數(shù)的相對(duì)計(jì)算誤差控制在3%以內(nèi)，均比COMPASS力場(chǎng)(7%)有了一定提升[41]。另外，Song等[76]結(jié)合COMPASS力場(chǎng)參數(shù)和RESP電荷描述DNTF晶體，晶胞參數(shù)的最大相對(duì)誤差僅為-2.86%，相比之下，采用COMPASS力場(chǎng)中的默認(rèn)電荷時(shí)，晶胞參數(shù)的最大相對(duì)誤差達(dá)-6.00%[77]。

3.3 模型尺寸

在分子動(dòng)力學(xué)模擬中，模型尺寸也是一個(gè)影響附著能計(jì)算結(jié)果的重要因素。顯然，模型越大，附著能計(jì)算越精確，耗時(shí)越長(zhǎng)。相反，模型尺寸也不能過(guò)小，至少要大于非鍵作用截?cái)喟霃降膬杀叮员苊庵芷谛赃吔鐥l件下分子與其鏡像發(fā)生相互作用[78]。因此，如何搭建合適的模型體系，在準(zhǔn)確性與耗時(shí)性上取得平衡，是一個(gè)需要考慮的問(wèn)題。在構(gòu)建晶面-溶劑界面模型時(shí)，研究人員通常采用3×3×3超晶胞的晶面模型和幾百個(gè)分子的溶劑模型，并沒(méi)有一個(gè)通用的構(gòu)建規(guī)則。

為此，Li lijie等[30,36,37,41,46]系統(tǒng)研究了模型尺寸對(duì)附著能計(jì)算結(jié)果的影響，包括晶面層、溶劑層、真空層的厚度以及模擬盒子的長(zhǎng)度和寬度。首先，考察了非鍵作用能隨截?cái)喟霃降氖諗壳闆r，認(rèn)為非鍵作用的截?cái)喟霃?dc)選取15.5?較為合適。然后，分別計(jì)算了溶液中附著能隨晶面層、溶劑層、真空層厚度和模擬盒子長(zhǎng)度、寬度的變化。最后，通過(guò)對(duì)β-HMX、FOX-7、ε-CL-20、NTO多種含能材料的測(cè)試，給出了一個(gè)合適尺寸的模型構(gòu)建規(guī)則：晶面層的長(zhǎng)度、寬度和厚度應(yīng)大于2dc、2dc和dc，溶劑層和真空層的厚度應(yīng)大于dc，即模型尺寸應(yīng)在圖5填色區(qū)域之外?；谠撃Ｐ统叽?，溶液中生長(zhǎng)形貌的預(yù)測(cè)結(jié)果與電鏡數(shù)據(jù)極為吻合。

圖5 模擬盒子的建議尺寸 Fig.5 Suggested size of simulation box

3.4 Eint的計(jì)算

探索溶劑在晶面上的微觀吸附行為，有助于理解溶劑效應(yīng)對(duì)晶面生長(zhǎng)的影響，從而更加準(zhǔn)確地計(jì)算晶面-溶劑相互作用。研究表明[24,31,63,67,79-82]，在晶面-溶劑作用的影響下，溶劑在晶面上呈現(xiàn)出了典型的結(jié)構(gòu)特征——層狀分布。此處，以HMX和丙酮溶劑為例，進(jìn)一步闡述溶劑的吸附特點(diǎn)。

圖6展示了丙酮溶劑(紅色)和HMX(黑色)沿生長(zhǎng)晶面(0 0 1)法向的質(zhì)量密度分布情況[31]。顯然，溶劑分子在HMX晶面法向上并不是均勻分布，而是呈層狀分布。隨著溶劑分子與HMX晶面法向距離的增加，溶劑峰的波動(dòng)逐漸減小，直至局域密度趨于體密度，說(shuō)明此處溶劑分子已經(jīng)不受晶面作用的影響，與體相中的分布一致。因此，若要完整體現(xiàn)出溶劑在生長(zhǎng)晶面上的結(jié)構(gòu)排布，溶劑層需要具有一定的厚度。盡管不同晶面結(jié)構(gòu)會(huì)導(dǎo)致局域密度有所不同，但均呈現(xiàn)出層狀分布的特征。

圖6 沿晶面法向的質(zhì)量密度分布Fig.6 Mass density profiles along the normal to the crystal face

為了分析層狀分布對(duì)晶面-溶劑相互作用能的影響，根據(jù)溶劑峰的個(gè)數(shù)將整個(gè)溶劑層進(jìn)行分層，如圖6中為L(zhǎng)1、L2和L3的色塊，然后分別計(jì)算區(qū)域L1、L1+L2、L1+L2+L3和整個(gè)溶劑層與生長(zhǎng)晶面的相互作用能。對(duì)于HMX和丙酮體系，HMX晶面與溶劑層L1的作用能約占總相互作用能的90%以上?；诓煌軇拥腅int，最終預(yù)測(cè)的晶形結(jié)果較為相似，只是長(zhǎng)徑比略有不同[31]。

需要注意的是，在動(dòng)力學(xué)模擬過(guò)程中，個(gè)別溶劑分子由于自由擴(kuò)散作用可能會(huì)吸附在生長(zhǎng)晶面的鏡像上，為了確保計(jì)算精度，在計(jì)算Eint時(shí)要將這些溶劑分子剔除。如果采用了“三明治”計(jì)算模型，即生長(zhǎng)晶面上下均有溶劑層，則不存在這個(gè)問(wèn)題[76]。同時(shí)，還應(yīng)注意附著能的化學(xué)計(jì)量。在真空中計(jì)算附著能時(shí)，若以晶胞為計(jì)量單位，需要將能量校正項(xiàng)Es轉(zhuǎn)換為與Eatt相同的單位。

4 附著能模型的發(fā)展

4.1 溶劑效應(yīng)的校正

如何更加合理準(zhǔn)確地描述溶劑效應(yīng)對(duì)生長(zhǎng)形貌的影響是附著能模型發(fā)展的一個(gè)方向。從式(2)可知，能量校正項(xiàng)Es體現(xiàn)了晶面與溶劑的相互作用，校正因子S相當(dāng)于生長(zhǎng)晶面結(jié)構(gòu)帶來(lái)的權(quán)重。值得注意的是，Eatt和Es的本質(zhì)都是相互作用能，區(qū)別在于：前者是晶面與生長(zhǎng)晶片的相互作用，而后者反映了晶面與溶劑的相互作用。然而，在真空中計(jì)算Eatt時(shí)，并沒(méi)有用S進(jìn)行校正。因此，在溶液中計(jì)算附著能時(shí)，也不應(yīng)采用S對(duì)Es進(jìn)行校正。這是因?yàn)椋涸谟?jì)算相互作用能時(shí)，已經(jīng)隱含了晶面結(jié)構(gòu)的影響，即晶面越粗糙，相互作用的接觸面積也就越大，同樣條件下的相互作用能也就越強(qiáng)，無(wú)需進(jìn)一步修正。

為了更好地解釋說(shuō)明，定義了接觸面積：

(10)

其中，Atot、A1和A2分別為體系和所包含組分的溶劑可及表面積。以HMX為例，計(jì)算表明，在真空中和溶劑中，晶面分別與生長(zhǎng)晶片和溶劑的接觸面積幾乎是一樣的[31]。因此，在同樣的接觸面積下，生長(zhǎng)晶面在真空中的相互作用能Eatt和在溶劑中的相互作用能Es的計(jì)算標(biāo)準(zhǔn)應(yīng)是一致的。實(shí)際上，校正因子S對(duì)晶面-溶劑相互作用額外加了權(quán)重，將對(duì)溶劑中的附著能帶來(lái)一定偏差，特別是結(jié)構(gòu)相對(duì)粗糙的生長(zhǎng)晶面，有時(shí)可能使得附著能為正值，致使晶面不穩(wěn)定。最近的研究基于這一新的校正公式，對(duì)溶劑中HMX的晶形進(jìn)行了預(yù)測(cè)，計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合較好[30-31]。

4.2 過(guò)飽和度的探索

附著能模型假設(shè)晶面的生長(zhǎng)速率與附著能成正比，即每個(gè)晶面的生長(zhǎng)速率是恒定的。因此，從原理上講，附著能模型難以處理生長(zhǎng)速率變化顯著的生長(zhǎng)條件，如過(guò)飽和度等[7]。

最近，Li lijie等[30,37]在構(gòu)建溶劑層時(shí)，將一定數(shù)量的溶質(zhì)分子也添加在內(nèi)，研究了過(guò)飽和度對(duì)β-HMX和ε-CL-20在溶劑中生長(zhǎng)形貌的影響。過(guò)飽和度(δ)可通過(guò)式(11)在模型構(gòu)建時(shí)進(jìn)行計(jì)算：

(11)

式中：C為計(jì)算模型中溶質(zhì)的實(shí)際濃度，由溶劑層中溶質(zhì)和溶劑的數(shù)量決定；C0是實(shí)驗(yàn)測(cè)定的飽和濃度。晶形預(yù)測(cè)結(jié)果表明，不同的過(guò)飽和度對(duì)形貌的長(zhǎng)徑比影響較大。由于結(jié)晶實(shí)驗(yàn)中的過(guò)飽和度很難保持不變，因此晶形預(yù)測(cè)結(jié)果尚不能與實(shí)驗(yàn)結(jié)果進(jìn)行對(duì)比。

4.3 晶形預(yù)測(cè)新策略

圖6表明溶劑沿生長(zhǎng)晶面法向呈層狀分布，當(dāng)溶劑分子遠(yuǎn)離晶面一定距離后，其與晶面的作用可忽略不計(jì)，此時(shí)溶劑可看作自由擴(kuò)散。當(dāng)溶劑距離晶面較近時(shí)，溶劑吸附在晶面上，對(duì)溶質(zhì)在晶面上的生長(zhǎng)起到一定阻礙作用。

最近的分子動(dòng)力學(xué)模擬表明，溶劑分子能夠長(zhǎng)時(shí)間地停留在晶面上的某些特定位置(吸附位點(diǎn))[82]，具有最強(qiáng)的結(jié)合能力，對(duì)溶質(zhì)生長(zhǎng)的阻礙作用最大，這說(shuō)明溶劑與晶面的相互作用也是不均勻的。然而，在計(jì)算溶液中的附著能時(shí)，溶劑效應(yīng)是通過(guò)“平均化”的晶面-溶劑相互作用能來(lái)描述的(見(jiàn)式(6))。這種溶劑層與晶面層的相互作用并不能真正體現(xiàn)出對(duì)溶質(zhì)生長(zhǎng)最為重要的抑制效應(yīng)。

實(shí)際上，吸附位點(diǎn)處的結(jié)合能可以看作是“木桶效應(yīng)”中的最短板，是決定晶體生長(zhǎng)速率的關(guān)鍵因素。圖7給出了晶體表面生長(zhǎng)過(guò)程的示意圖。由于溶劑在吸附位點(diǎn)處的結(jié)合能力最強(qiáng)，溶質(zhì)將吸附位點(diǎn)處的溶劑替換則可以看作是晶面生長(zhǎng)的“限速”步驟。為此，基于附著能模型的理論框架，提出了晶形預(yù)測(cè)新策略，即采用吸附位點(diǎn)處溶劑與溶質(zhì)的結(jié)合能代替晶面-溶劑相互作用能，采用吸附位點(diǎn)處溶質(zhì)與溶質(zhì)的結(jié)合能代替真空中的附著能，具體計(jì)算細(xì)節(jié)可參考文獻(xiàn)[82]。顯然，相比于溶劑中原有附著能的計(jì)算，從晶體生長(zhǎng)過(guò)程上看，這種策略更加合理。基于該策略，預(yù)測(cè)了HMX在丙酮溶劑中的晶形，計(jì)算結(jié)果與實(shí)驗(yàn)結(jié)果吻合良好[82]。

圖7 晶體表面生長(zhǎng)過(guò)程示意圖 Fig.7 Scheme in the crystal growth process in solvent

5 結(jié)束語(yǔ)

附著能模型是目前含能材料領(lǐng)域中應(yīng)用最為廣泛的晶形預(yù)測(cè)方法，能夠定性地給出溶劑、溫度等不同條件下的生長(zhǎng)形貌，可為含能材料的晶形控制提供理論參考。在使用附著能模型時(shí)，需要充分考慮力場(chǎng)選擇、模型構(gòu)建、溶劑效應(yīng)等因素，從而準(zhǔn)確地計(jì)算附著能。否則，相同體系的預(yù)測(cè)結(jié)果也可能會(huì)存在差異。值得注意的是，附著能模型僅僅基于恒定生長(zhǎng)速率的簡(jiǎn)單假設(shè)，尚不能考慮結(jié)晶實(shí)驗(yàn)中的諸多因素，預(yù)測(cè)晶形有時(shí)與實(shí)驗(yàn)結(jié)果偏差較大也是合理的。為此，仍需進(jìn)一步發(fā)展基于生長(zhǎng)機(jī)制的晶形預(yù)測(cè)方法。