稀有氣體二聚體基準集NgD×15構建及其對DFT方法的評價

王裕平, 魏孝珍，王一波

(貴州大學化學系 貴州省高性能計算化學重點實驗室, 貴陽 550025)

1 引 言

目前在分子間相互作用領域，用于訓練和評價各種近似理論方法的基準集主要有Hobza研究組發展的A24[1]，S22[2]，S22×5[3]，S66，S66×8[4]以及X40，X40×10[5]等，這些基準集未包含稀有氣體二聚體；Grimme研究組發展了GMTKN24[6]，GMTKN30[7]，GMTKN55[8]等基準集數據庫，其中GMTKN30僅給出了Ne2至 Rn2平衡位置數據，未給出勢能曲線；Head Gordon在MGCD84數據庫的RG10基準集中仍沿用了Tang-Toennies經驗勢[9].

1998年Partridge等用CCSD(T)方法計算了He2解離能，同時考察鍵函數的作用[10].其后Hellmann，Sheng等計算并擬合了Ne2[11]，Ar2[12]，Kr2[13]，Xe2[14]和Rn2[15]的勢能曲線.但是由于這些研究結果所采用的方法、基函數不一致，難以構成一套完整、系統和精確的基準集.以至于2017年Kovacs等評測13種DFT-D3方法時，是以自己的CCSD(T) /aug-cc-pVTZ計算結果為評測標準[16]，Austin 等在發展APF-PFD方法時也是采用類似的做法[17].以CCSD(T) /aug-cc-pVTZ結果作為標準，基函數還遠遠沒有收斂，對DFT方法的評價不嚴謹.作者在發展B972-PFD[18]DFT方法時也是缺乏可靠的Ng2參考基準.因此建立完整、精確的Ng2基準集非常必要.

根據前人的大量計算實踐，在基函數極限下用CCSD(T)計算結果構建Ng2體系相互作用勢是合適的，但該體系作用勢很小，對基函數非常依賴，且收斂很慢，關鍵在于使用收斂的基函數.中點鍵函數能夠有效地加速基函數收斂，本文用aug-cc-pV5Z核中心基，加上陶福明等所推薦的中點鍵函數[19, 20]，仔細考察了鍵函數的影響，精確地給出Ng2體系的勢能曲線，建立了一個較為準確、完整的Ng2基準集，并重新測試、評價了目前常用的DFT色散及非局域校正方法.

2 NgD×15基準集的構建

本文所建立的NgD×15基準集，對He2，Ne2，Ar2，Kr2，Xe2和Rn2每個體系各取15個測試點，共計90個點.我們確定核間距時，以平衡核間距Re×0.9為起始點，以0.2 ?或0.3 ?及以上的間距作為步長取點，核間距具體數值見表2.

3 研究方法

首先在CCSD(T)[21]計算水平下，以He2為例分別使用aug-cc-pVXZ(X=T, Q, 5, 6)以及d-aug-cc-pVXZ(X=7, 8)基組，并在aug-cc-pV5Z基組下分別加上{6s6p6d3f2g1h}，{6s6p6d3f2g2h}，{6s6p6d6f3g2h1i}以及{6s6p6d6f3g3h2i1k}等高角量子數鍵函數集(以下分別標記為BF1-BF4)，研究相互作用勢對基函數的收斂性.所有計算均使用Boys和Bernard[22]所提出的均衡校正法(Counterpoise Procedure, CP) 消除BSSE.

以CCSD(T)/ aug-cc-pV5Z-{6s6p6d3f2g1h}方法建立基準集后；對B97M-V[23]，ωB97M-V[24]，ωB97X-V[25]，B97M-rV，ωB97X-rV，ωB97M-rV[26]，ωB97M(2)[27]等7種DFT-NL方法以及Kovacs等所測試的DFT-D3方法[16]，在def2- QZVPPD基組，(250, 974)積分網格下進行了測試，并統計分析了各種方法的誤差.

全部研究工作采用Gaussian16 RevA.03[28]、Q-Chem5.1[29]程序，在貴州大學云計算平臺上完成計算.

4 結果與討論

4.1 NgD×15基準集

圖1是以He2體系為例的計算結果.由圖1可知：CCSD(T)/d-aug-cc-pV7Z計算水平下，He2的勢能曲線計算沒有得到收斂；直到d-aug-cc-pV8Z，He2勢能曲線才接近收斂.比較CCSD(T)/aug-cc-pV5Z- BF1的計算結果發現，此時的計算精度已經達到CCSD(T)/d-aug-cc-pV8Z水平；當鍵函數再進一步增大到{6s6p6d6f3g3h2i1k}時，其勢能面的最低點并無明顯的變化，最大誤差只有0.015cm-1.因此，CCSD(T)/aug-cc-pV5Z- BF1已趨于收斂.

圖1 CCSD(T)方法在aug-cc-pVXZ(X=5, 6)，d-aug-cc-pV XZ(X=7, 8)以及 aug-cc-pV5Z-BF(1-4)基組下對He2勢能面的計算結果.Fig.1 The results of potential energy surface calculated for the He2 at the CCSD(T)/ aug-cc-pVXZ(X=5, 6)，CCSD(T)/ d-aug-cc-pVXZ(X=7, 8) and CCSD(T)/ aug-cc-pV5Z-BF(1-4) levels.

我們繼續在CCSD(T)/ aug-cc-pV5Z- BF1水平計算了Ne2至Rn2體系，并且與aug-cc-pVXZ(X=T, Q, 5)基組結果進行了對比.

表2為CCSD(T)/ aug-cc-pV5Z-BF1計算水平下，所有同核Ng2體系的相互作用勢在不同核間距下的具體數值.我們發現Head-Gordon等所使用的Tang-Toennies經驗勢結果與NgD×15計算結果相一致.對于He2體系誤差只有0.1 cm-1，絕對誤差最大的Kr2，也只有5.5 cm-1，說明稀有氣體的Tang-Toennies經驗勢模型是可信的.

表1 NgD×15的平衡核間距(?)與解離能(cm-1)文獻參考值與本次研究相關的計算數值

Table 1 Literature values and calculated values of this study for the equilibrium internuclear distancesR(?) and dissociation energiesΔE(cm-1) of NgD×15

DimerRefThis workRΔERΔEHe22.980[10]7.4502.9767.430Ne23.10[11]28.753.10628.74Ar23.80[12]97.273.79097.31Kr24.0[13]134.44.067134.8Xe24.4[14]192.24.392194.9Rn24.483245.4

Ref[10-14]為分別在CCSD(T)/aug-cc-pVTZ-{6s6p6d3f3g3h}，CCSD(T)/t-aug-cc-pV6Z-{4s4p3d3f2g}，CCSD(T)/d-aug-cc-pV6dZ-{4s4p3d3f2g}，CCSD(T)/t-aug-cc-pV6Z-{4s4p3d3f2g}，CCSD(T)/d-aug-cc-pV6Z-{4s4p3d3f2g}計算水平下的結果；this work為CCSD(T)/ aug-cc-pV5Z- BF1的計算結果.

表2 NgD×15核間距(?)與對應的相互作用勢(cm-1)

4.2 DFT-D3和DFT-NL方法的評價

以下表3是以表2數據為基準集，對DFT-NL以及DFT-D3等方法計算結果的絕對均方根誤差RMSD以及相對均方根誤差RMSRD的具體數值.

由表3可知，在Kovacs等人的測試結果中表現較好的B3LYP-D3(BJ)，CAM-B3LYP-D3，B3PW91-D3，PBE0-D3以及B97-D3等方法在此測試中與標準參考值的總體誤差相對還是比較大的；尤其是對于He2體系的計算，絕對均方根誤差RMSD數值在4-22 cm-1之間；PBE-D3，TPSS-D3泛函，其絕對均方根誤差RMSD已經超出15 cm-1.DFT-NL方法中表現較好的ωB97X-V以及ωB97X-rV方法對He2體系計算的絕對誤差雖然只有0.9 cm-1左右，但是相對誤差RMSRD依然較大.所以對于He2這種極弱的相互作用體系，DFT-D3以及DFT-NL均不能提供較好的描述.從He2至Rn2絕對誤差增大，而Rn2的絕對誤差大都在20-60 cm-1以上；在這些DFT-NL方法中，除了B97M-V以及ωB97M-V，ωB97M(2)方法的誤差略大，其余的均在25cm-1之間，DFT-NL方法相對于DFT-D3優勢明顯，特別是ωB97X-V，ωB97X-rV兩種DFT-NL方法.

表3 DFT-NL以及DFT-D3(BJ)等方法的測試誤差分析結果(cm-1)

5 結 論

本文在CCSD(T)/aug-cc-pV5Z-{6s6p6d3f2g1h}水平下，經BSSE均衡校正后，獲得了He2至Rn2全部同核稀有氣體二聚體的15點相互作用勢，構建了NgD×15基準集.并以此為標準重新評價了DFT-D3和DFT-NL方法對色散作用的描述，發現DFT-NL總體上要強于DFT-D3方法.Kovacs等的測試結果中表現最好的B3LYP-D3， B3PW91-D3，PBE0-D3以及B97-D3等方法絕對均方根誤差RMSD相對于DFT-NL方法較大.在所有的DFT-NL方法中，ωB97X-V及ωB97X-rV兩種方法是該體系最優性價比的方法.