計算化學鍵鍵能時考慮BSSE不僅是多余的甚至是有害的
1 前言
BSSE在很早以前寫的《談談BSSE校正與Gaussian對它的處理》(http://www.shanxitv.org/46)中已經介紹過了。也正如文中提到的以及在無數文獻中大家看到的,考慮BSSE影響的場合幾乎都是在計算弱相互作用的時候。雖然BSSE可以說無處不在,只要使用原子中心基組時就必定有BSSE效應,但實際上,通過整體減片段能量來得到化學鍵鍵能的時候,只要基組用得基本靠譜,那么再用常規的counterpoise方式考慮BSSE問題就純屬多余,白浪費時間。很多初學者對BSSE的本質一知半解,容易胡思亂想,對BSSE可謂敏感得至極,稍一牽扯整體減片段的計算,不管是什么具體情況,腦子中都立馬瞬間閃現“BSSE”四個字母(大量初學者對于自旋污染也是這個毛病,凡是開殼層,言必曰自旋污染)。另外,很多審稿人也是一知半解就瞎出餿主意,明明文章算的是化學鍵鍵能,基組也合理,某些低水平審稿人偏偏非要讓考察BSSE影響(筆者發的一篇JCC就碰到過,被筆者義正言辭懟回去了,之后對方也沒再說什么),搞得一些本來持觀念正確的作者,也開始懷疑自己的立場。一些IF不低的期刊文章中算鍵能用counterpoise的例子也不少,導致一些初學者錯誤地盲目模仿。
恰逢昨日公社群里又有人談起自己文章中算化學鍵鍵能時,審稿人讓自己考慮BSSE。為了以正視聽,這里通過兩個極為簡單的計算化學鍵鍵能的例子,通過考察不同基組時的結果,證明算化學鍵鍵能時只要基組合理,考慮BSSE問題不僅純屬多余,甚至還可能令精度不增反降。如果不知道基組怎么樣選才算合理,看此文《談談量子化學中基組的選擇》(http://www.shanxitv.org/336)。第一個體系是乙烷的C-C鍵均裂成兩個甲基,第二個體系是NaCl異裂成Na+和Cl-,這就把最典型共價鍵和離子鍵的情況都考慮了,至于極性共價鍵會在文末再做簡單提及。本文理論方法使用HF。參考結果是在jul-cc-pV5Z下計算得到的,此時基組完備性極高,結果已達到HF極限,BSSE的影響可以認為已經為0。計算通過G09 D.01完成。
校正BSSE問題的方式不止Boys提出的counterpoise一種,Grimme提出的gCP也挺知名(這里介紹了:http://www.shanxitv.org/214),但為了敘述起來簡單,下面凡是說BSSE的時候都是指使用counterpoise方式考慮的。要注意正確理解本文傳達的意思:在基組選用合適的情況下,計算鍵能時不適合用counterpoise方式試圖刻意修正此時本來就沒必要考慮的BSSE問題,刻意這么做甚至還有害處。而帶著gCP的話,對這種情況既沒額外的好處也沒壞處,因為對于短距離的作用gCP對結果不起明顯影響,而對于中、長距離的真正可能需要解決BSSE對弱相互作用不良影響的時候gCP才會起到作用。
2 乙烷均裂
乙烷均裂成兩個甲基自由基過程的能量變化如下,可以視為是此體系C-C鍵鍵能。計算時沒有考慮斷鍵后甲基自由基的結構弛豫。raw一列對應沒考慮BSSE的結果,corrected是考慮BSSE后的結果,E(BSSE)就是BSSE校正能(對整體能量的修正量)。
我們先只看圖中虛線上方的數據,可以看出以下幾點:
(1)平時我們試圖計算定量準確數據的時候,基組都得3-zeta或以上。從數據看到在cc-pVTZ這個檔的時候,E(BSSE)已經非常小了,考不考慮BSSE無所謂,本身理論方法的誤差比這大得多。到了cc-pVQZ時,BSSE則幾乎不產生任何影響。
(2)在基組基本合理的情況下,具體來說從cc-pVDZ開始,BSSE校正后結果并非變得更好!與參考值相對比,可見不考慮BSSE時高估了鍵能,而考慮BSSE后又低估了鍵能。仔細看的話,從cc-pVTZ開始,或者從cc-pVDZ加彌散開始,考慮BSSE后結果的誤差反倒更大!這反映出counterpoise方法本身的問題,往往矯枉過正,為什么會這樣在《談談BSSE校正與Gaussian對它的處理》中提到過。
(3)STO-3G這種破基組誤差極大,算出來的鍵能誤差高達25kcal/mol。考慮BSSE后結果誤差依然很大,但改進不小,達到14kcal/mol。改進程度這么大是容易理解的,畢竟極小基完備性非常低,考慮BSSE校正勢必會帶來顯著補償效果。然而,誰在實際中也不會用極小基,因此完全不能用來說明考慮BSSE對于計算鍵能是有必要的。
(4)什么時候該加彌散在《談談彌散函數和“月份”基組》(http://www.shanxitv.org/119)已經充分談過了,當前的計算顯然不屬于此文提及的情況。因此正如預期的,給cc-pVDZ加彌散前后結果的差異甚微,才不到0.2kcal/mol,這還是對于小基組的情況。而比較表中cc-pV5Z和參考值jul-cc-pV5Z的結果,則看到加彌散對結果沒有絲毫影響,因為cc-pV5Z對此體系已經充分達到HF極限了。
(5)從數據可見,HF計算(也包括DFT計算),cc-pVTZ已經足夠大了,提升到cc-pVQZ根本沒什么意義,對結果影響的程度已經遠低于方法本身的誤差。只有高級別電子相關方法才有必要上到cc-pVQZ。
筆者還順手測了幾個其它基組,即上圖中虛線下面的數據。從數據中可見,給6-311G**增加更多極化函數,即增大到6-311G(2df,2p)并不會給精度帶來什么改進。比TZVP(具體指def-TZVP)極化函數更多的def2-TZVP也并沒顯示出精度比TZVP有什么優勢。這也正一定程度反映出《談談量子化學中基組的選擇》里提到的對于HF/DFT加很多極化函數是白浪費時間(當然,極化函數起到的改進和體系、研究的問題有關,不能以偏概全)。從數據中還看出,TZVP的結果就已經十分接近參考值了,比大小相仿佛的6-311G**和比它更大的6-311(2df,2p)結果還好,這也是我為什么鼓勵使用def/def2基組代替常用的Pople系列基組的原因。def2-TZVP雖然對此體系的raw結果和TZVP相同,但并非比TZVP沒有本質上改進,從它更小的BSSE校正能上可以間接反映出def2-TZVP的完備性比TZVP更高。注意TZVP、def2-TZVP在BSSE校正后誤差比原先更大,再次反映出只要基組質量不錯,去額外花時間做counterpoise來試圖校正BSSE只會自取其辱。
3 NaCl異裂
NaCl解離成Na+和Cl-過程的能量變化,也即NaCl的異裂鍵能如下。
由于無論是在NaCl中,還是Cl-狀態,Cl都帶明顯負電荷,因此測試中對每個基組都考慮了帶和不帶彌散函數時的結果,因為眾所周知,彌散函數對于計算陰離子體系或含有顯著帶負電荷的原子的體系的能量必須的,前面列的博文里都強調過。
先看不考慮BSSE時的數值,即raw那一列。從數據可見,彌散函數對此體系是非常重要的。雖然隨著cc-pVnZ的序數n增大,結果也能逐漸向參考值收斂,但是有彌散函數時誤差小得多。比如aug-cc-pVTZ的誤差都已經明顯小于cc-pV5Z。哪怕只加最低限度的彌散函數,也能大幅改進結果。比如對cc-pVDZ只加上一層sp彌散函數而忽略d彌散函數成為maug-cc-pVDZ,誤差從8.8kcal/mol降低到3.0kcal/mol。如果把d彌散函數也考慮則還會有進一步改進,誤差就只有1.6kcal/mol了。
表格中(aug)-cc-pVTZ指的是只對帶顯著負電荷的Cl用aug-cc-pVTZ,而Na還用cc-pVTZ。可見結果和aug-cc-pVTZ幾乎沒任何差別。這充分表明,算離子鍵鍵能時,只對帶明顯負電的原子用彌散函數就足夠了。
然后我們再看BSSE如何影響結果。從數據中可見只要帶了彌散函數,BSSE的影響就微乎其微了,完全不需要考慮之。而不加彌散函數的時候BSSE校正能挺大,且考慮BSSE時結果會明顯與參考值更相符。另外還看出,隨著cc-pVnZ的序數增加,由于基組逐漸變得完備,BSSE校正的影響逐漸降低。
帶彌散函數(只給帶顯著負電的原子帶就夠了)和通過counterpoise考慮BSSE,都能使cc-pVnZ計算離子鍵鍵能的結果有顯著改進,而又帶彌散函數(哪怕只是最低程度的單層sp彌散)又考慮BSSE則毫無必要。而帶著彌散函數時,不僅鍵能的描述得到改進,當前體系其它各方面性質的描述都會得到很大改進(如靜電勢、極化率),而用counterpoise考慮BSSE的話改進的只有鍵能,而且在優化、振動分析的時候則幾乎沒法用,因為此時沒有解析導數(而且往往還沒法利用對稱性加速計算)。
4 總結與討論
雖然本文只通過最簡單的非極性共價鍵和離子鍵體系的數據進行討論,但可以以小見大。最主要結論是,只要基組選用合理,那么計算化學鍵鍵能根本沒必要考慮BSSE!用counterpoise來考慮的話,甚至結果精度反倒可能更差。何謂基組選擇合理?只要基于《談談量子化學中基組的選擇》文中的討論就行了。對于算鍵能來說,基組這樣設是合理且劃算的:算共價鍵體系鍵能時用3-zeta級別基組(如果是后HF方法,建議4-zeta級別),彌散函數不需要加;算離子鍵鍵能時,彌散函數一定要有,但只給帶顯著負電的原子加上就行了。對于算鍵能,考慮BSSE并非完全無用,但僅對于基組明顯太小,或者該加彌散時候沒加的時候才對結果有顯著改進。相信除了菜鳥以外,誰也不會把基組選得那么爛。
當你的文中計算了鍵能,而且基組使用合理的時候,卻遇到審稿人亂提要求讓你考慮BSSE的時候,大可直接引用本文的數據給懟回去。
作為本文例子的乙烷均裂和NaCl異裂屬于很理想的情況,而對于實際經常碰到的極性鍵,counterpoise方法甚至從原理上就根本沒法用。counterpoise是通過下面的方式計算對AB整體能量的修正量的:
CH3當陽離子、F當陰離子,E(BSSE)=8.30kcal/mol
CH3和F都當中性自由基,E(BSSE)=0.52kcal/mol
可見,用counterpoise方式考慮BSSE時,計算的鍵能會直接依賴于CH3和F狀態的選取,結果差異高達7.8kcal/mol。然而無論選取哪種狀態都不合理,因為實際中CH3F分子中甲基只是轉移了零點幾個電子給F。故強行使用counterpoise必會導致最終計算的鍵能誤差增加,完全是費力不討好!
綜上所述,計算鍵能時基組選好了就完了,不要管BSSE。希望此文能夠以正視聽。