revTPSS泛函結合pcSseg-1基組是計算NMR很好的選擇
revTPSS泛函結合pcSseg-1基組是計算NMR很好的選擇
文/Sobereva@北京科音 2021-Oct-24
1 前言
計算NMR化學位移是量子化學常見而且重要的應用之一，之前筆者也寫過一篇文章《談談如何又好又快地計算NMR化學位移》（http://www.shanxitv.org/354）做了相關討論，不了解這方面計算者強烈建議先看看。《使用Multiwfn繪制NMR譜》（http://www.shanxitv.org/565）里也有許多相關知識講解。
最近在JCTC期刊上出現了一篇文章https://doi.org/10.1021/acs.jctc.1c00604，里面對一批泛函計算氫的化學位移的精度進行了橫測，使用的是一批有機分子在氯仿中的氫的實驗化學位移作為參照。他們是用的最普通的NMR計算方式，即優化極小點，算各項同性磁屏蔽值，將TMS與當前體系的值相減折算成化學位移然后和實驗對照。此文聲稱revTPSS泛函（meta-GGA檔次）算氫的化學位移很好，不輸于依賴擬合參數的標度法，并發現比普通泛函昂貴得多的雙雜化泛函、MP2算氫的NMR的精度平平，比B3LYP等普通泛函并無優勢，更是顯著輸于revTPSS。這些測試結論和之前很多橫測不同方法算NMR精度的文章的結論頗有差異，例如在JCTC, 10, 572 (2014)的測試中MP2算氫的NMR精度強于所有所測試的普通泛函，在JCTC, 14, 4756 (2018)的測試中發現雙雜化泛函DSD-PBEP86算各種元素化學位移的精度顯著好于普通泛函，也略好于MP2。結論有這么明顯差異的原因有多方面，一方面是測試集明顯不同（后兩篇里一般有機體系很少），還有一大原因是后兩篇文章都是用的高精度的CCSD(T)在平衡結構的數據作為參照值，但平衡結構下的CCSD(T)的數據和實驗值相比也是有明顯差異的，除了方法本身誤差外還來自于振動校正、對溶劑效應的考慮、用的幾何結構等方面的誤差。對于搞應用性計算的研究者，真正感興趣的無疑是在某個便于使用的級別下，且在不考慮昂貴、復雜的振動校正和顯式溶劑等情況下，計算的結果和實驗相差多大；他們需要的是一個整體誤差小、計算省事而且又普適的級別，哪怕表現好的實質上是由于和其它因素誤差抵消得較充分也無妨。所以上面提到的第一篇JCTC文中直接以實驗數據作為參照衡量誤差的測試結論明顯更具有實際意義。
筆者下面將對三個有代表性的有機小分子計算NMR，有以下幾個目的
(1)實際檢驗一下revTPSS算氫的NMR是否如前面JCTC文章里說得表現很好，和文章里該測卻沒測的B97-2、KT1、KT2泛函相比如何、和標度法相比如何
(2)前面JCTC文章里沒有測試revTPSS算碳譜的NMR的精度，但這也是很重要的問題，這里做一下粗略測試
(3)看看revTPSS結合哪種基組表現最好。JCTC文中的圖中體現此泛函結合pcSseg-1精度最好而且又便宜，但文中推薦的卻是revTPSS/cc-pVTZ，非常蹊蹺
(4)看看幾何優化級別對revTPSS算NMR結果的影響，以確定一個比較劃算的protocol
(5)看看雙雜化泛函算NMR和revTPSS相比如何，值不值得花費比普通泛函更多的時間
本文的根本目的是確定算普通有機體系NMR理想且劃算的級別。
雖說做計算級別的精度測試必須用很大的樣本才有說服力，但筆者沒有那么多精力，拿三個分子也能多多少少獲得一定有用的結論，對于驗證上述JCTC文章里的測試結論是否靠譜更是足夠了。
2 計算結果和分析
本文測試的三個分子如下所示，其中惡丁環和異惡唑在《談談如何又好又快地計算NMR化學位移》（http://www.shanxitv.org/354）中也被當過測試例子。實驗的氫1和碳13的化學位移數據是氯仿下測的。
實驗值和不同方式計算的以上三個分子的化學位移（ppm），以及相對于實驗值的平均絕對誤差(MAE)結果如下。除標度法外，用的幾何結構都是B3LYP/def2-TZVPP在真空下優化得到的。其中revTPSS、B97-2、B3LYP、MP2計算是Gaussian 16在IEFPCM溶劑模型表現的氯仿環境下算的。DSD-PBEP86雙雜化泛函是用ORCA 5.0.1算的，分別在真空和CPCM表現的氯仿環境下計算。“標度法”是《談談如何又好又快地計算NMR化學位移》里的“標度法1”，即B3LYP/6-31G*真空下優化然后在B3LYP/6-31G*結合SMD溶劑模型表現的氯仿環境利用前人擬合的標度法參數得到的（為了和前人擬合標度法參數時用的Gaussian版本一致，刻意用的Gaussian 09做計算）。雖然也有很多其它計算級別的標度法，但沒有比當前用的這個明顯更好的。
通過以上數據可以看出（注意必須三個體系的數據同時綜合考慮，不要光盯著個別數據說事，后同），revTPSS結合恰當的基組算氫的化學位移相當出色。pcSseg系列是專門給計算NMR搞的基組，在計算NMR方面性價比普通基組高得多。其中pcSseg-1大約是6-31G**的大小，revTPSS/pcSseg-1既便宜結果又好，故在此強烈推薦拿它算氫的NMR。這個級別和標度法的耗時基本一樣，精度互有勝負，反正都很好。若你不想用標度法（嫌經驗性強等原因），那么revTPSS/pcSseg-1是首選，而用標度法的時候也可以同時用revTPSS/pcSseg-1算一算，令結果相互檢驗，更放心。如果你用其中一個方法算出來的結果和實驗值相比存在較大誤差，此時可以試另一個。
revTPSS/pcSseg-1算碳的化學位移的精度明顯不如算氫的時候那么準，但至少比其它普通泛函結合較貴的def2-TZVPP時還整體更好，不過整體略遜于標度法。所以，算碳的NMR我優先建議用標度法，而必須用標準方式計算的時候則推薦使用revTPSS/pcSseg-1。
上面表格里也簡單測試了revTPSS計算結果受基組的影響。用更貴的pcSseg-2反倒結果更差，這體現出revTPSS/pcSseg-1結果好有一部分因素是泛函和基組誤差抵消得較好，而且抵消得好還不是個例，而是比較普遍的，這在前述的JCTC文章的圖3的對比測試中也有所體現。revTPSS結合普通的2-Zeta檔次基組def2-SVP遠不如結合與之大小差不多的pcSseg-1，這是因為def2-SVP實在太小，又不是專門給NMR計算設計的屬性基組。revTPSS結合常用且較好的def2-TZVPP基組的時候遠強于def2-SVP，但比pcSseg-1沒優勢，再加上又貴得多，因此不值得考慮。由于revTPSS在結合較好的通用基組時表現也很不錯，說明revTPSS本身算NMR就是不錯，而非必須明顯依賴于與基組誤差的抵消。
和氯仿下的實驗NMR對照時，用隱式溶劑模型表現氯仿環境不一定會使得結果更好，一方面是氯仿的極性很弱，溶劑效應對結果產生的影響有限，另一方面還有誤差抵消的因素，所以考慮溶劑模型是否會改進精度有一定不確定性。而有意思的是，雙雜化泛函DSD-PBEP86在CPCM表現的氯仿環境下算的NMR明顯不如在真空下的，因此強烈不建議此時帶溶劑模型。我有所懷疑是ORCA 5.0.1在雙雜化泛函算NMR的過程中對溶劑模型產生的影響考慮得不充分。在真空下，DSD-PBEP86結合較好的def2-TZVPP基組時比B97-2、B3LYP、MP2整體都更好，但是比revTPSS/pcSseg-1并沒啥優勢。所以，我認為碰到某些情況，如果發現標度法、revTPSS/pcSseg-1和實驗對比時表現得都不佳的話，才值得再去考慮用明顯貴得多的DSD-PBEP86/def2-TZVPP，而且應當在真空下算（假定實驗是氯仿環境）。
值得一提的是，Gaussian在MP2下做NMR是非常昂貴的，主流36核服務器在MP2/def2-TZVPP下算吡啶的NMR都得花六分鐘。ORCA里利用RI近似做MP2或雙雜化泛函的NMR計算則相對來說快得多得多的，在i7-10870H這種主流Intel 8核CPU上用DSD-PBEP86/def2-TZVPP算吡啶的NMR才花了一分半鐘就算完了。
在上面的測試做完后，筆者還有些懸念，遂又做了幾個額外的測試，見下表。KT系列泛函（包括KT1、KT2、KT3）是極少數專門考慮算磁屬性的泛函，在此泛函原文里KT1算NMR是KT系列里最好的，但很多NMR測試文章里只考慮了KT2，并發現KT2完勝其它測試的泛函。因此這里用ORCA 5.0.1對KT1、KT2也做了測試，并且對前者在真空（下表里標注了vac）和CPCM表現的氯仿下（下表里標注了solv）都做了計算。另外也考察了在revTPSS/pcSseg-1在真空中的表現。最后還考察了DSD-PBEP86結合很大基組pcSseg-3的情況，此基組對于計算NMR來說可以算CBS了。
由上面數據以及和之前表格里的數據對比可見，在真空下，KT2其實比KT1好一些，KT2結合def2-TZVPP算氫的精度和revTPSS/pcSseg-1基本可以打成平手，但鑒于KT泛函被量子化學程序支持得較少（Gaussian不支持，ORCA得通過libxc界面才能用），所以有了revTPSS就沒必要用它了。考慮溶劑模型并沒給KT1計算結果帶來什么改進，對于算氫還有所劣化。KT1和KT2結合def2-TZVPP算碳的化學位移比revTPSS/pcSseg-1差得很遠。
revTPSS/pcSseg-1在真空下計算結果和在IEFPCM表現的氯仿環境下互有勝負，刻意用真空并不會令整體變得更好。用revTPSS/pcSseg-1的時候還是建議考慮氯仿溶劑環境，畢竟原理上更嚴格，耗時也多不了多少。
上表中的DSD-PBEP86/pcSseg-3在計算的時候用的是ORCA 5.0發布會的ppt上列的關鍵詞（其實RI-以及D3BJ都是多余的)
! RI-DSD-PBEP86/2013 D3BJ pcSseg-3 cc-pwCVQZ/C tightSCF NMR defgrid3
%basis AuxJ "def2/JK" end
從上表的結果來看，雖然算氫的精度相比于用def2-TZVPP時有一絲改進，但算碳的精度則下降很多。這也體現出DSD-PBEP86算NMR方面比revTPSS沒優勢是泛函自身層面的事，用很大的基組也并不改變結論。
最后，筆者測試了一下基于B3LYP結合不同基組（6-31G*、def2-SVP、def2-TZVPP）在真空中優化的結構，在上面推薦的revTPSS/pcSseg-1結合IEFPCM表現的氯仿情況下算的NMR的精度，如下所示。B3LYP優化普通有機分子的精度眾所周知足夠好，見《談談量子化學研究中什么時候用B3LYP泛函優化幾何結構是適當的》（http://www.shanxitv.org/557），因此不考慮用別的泛函優化。
可見結合不同基組優化結構，對最終的NMR計算結果有些許影響，但影響不大，也沒有什么明顯規律性。一個很好的發現是優化只需要用很便宜的6-31G*就完全夠了，用貴得多的def2-TZVPP不會帶來什么改進。也可以說，如果是為了算NMR，只需要用B3LYP/6-31G*在真空下優化即可；而如果是為了其它某些目的在其它恰當級別下優化了結構，也可以直接用這結構來算NMR，不用為了算NMR而再額外特意用某級別來優化。
3 結論
本文利用了三個典型的有機小分子，通過測試證明了一開始提到的JCTC上測試文章的結論合理性，即revTPSS特別適合算氯仿環境中的有機體系的氫的化學位移，而且本文進一步明確指出revTPSS/pcSseg-1（IEFPCM表現的氯仿環境下）和標度法（基于B3LYP/6-31G*的）算氫的化學位移在耗時和精度上都相仿佛，前者的好處是經驗性低，不依賴于針對元素擬合的參數，因此是更普適的方法。在算氯仿中有機體系的碳的化學位移方面，revTPSS/pcSseg-1也不明顯輸于標度法，比其它普通泛函以及MP2結合較大基組時明顯更好，昂貴得多的雙雜化泛函相對于它也沒優勢。在revTPSS/pcSseg-1算NMR之前的幾何優化方面，用B3LYP/6-31G*在真空下優化即可。
相對于revTPSS/pcSseg-1和標度法，至少對于算氯仿中有機體系的氫、碳的NMR方面雙雜化泛函沒什么優點，不值得為此花高得多的計算代價，但不排除它在計算其它類型體系、其它元素的NMR的時候會有相對于revTPSS泛函在內的普通泛函更好的表現。
由于revTPSS/pcSseg-1和標度法表現出色，ORCA又支持雙雜化泛函算NMR，導致MP2在NMR計算上沒有什么實用性了。KT2算氫的NMR很不錯，但算碳的NMR不推薦，有了revTPSS/pcSseg-1就沒必要考慮它了。
需要強調的是，對于精確計算能夠與實驗相對照的NMR化學位移來說，振動校正不可忽略，但本文始終沒有考慮這一點，因為這考慮起來對于中小體系耗時都超極高，而且Gaussian還不直接支持（但ORCA 5.0開始支持VPT2非諧振模型下考慮這種校正）。所以記住，本文關于精度的結論實際上是與忽略振動校正等因素帶來的誤差相抵消之后的，并不是直接反映計算級別內在的計算磁屏蔽張量的理論精度。由于平時我們做NMR計算一般都不考慮振動校正，所以本文的結論對于實際應用性研究是有直接指導意義的。
鑒于用B3LYP/6-31G*在真空下優化，然后用revTPSS/pcSseg-1在IEFPCM描述的氯仿算氫和碳的NMR非常有價值，筆者已經將Gaussian 16下這個情況的參考物質TMS的碳和氫的磁屏蔽值納入到《使用Multiwfn繪制NMR譜》（http://www.shanxitv.org/565）介紹的Multiwfn的功能里了（從2021-Oct-23更新的版本開始），在繪制NMR的界面里的選項7 Set how to determine chemical shifts里選1 Set reference shielding value to determine chemical shift，再根據屏幕提示選擇相應選項就可以把讀入的此級別的磁屏蔽值轉化為化學位移，然后可以直接作圖，省得自行計算TMS的值并且每次都手動輸入了。
4 附：計算用的關鍵詞和使用細節
這里對上文涉及的計算方法的關鍵詞使用做一些說明。
revTPSS在Gaussian中要寫成revTPSSrevTPSS，代表交換泛函和相關泛函都是revTPSS文中定義的。
pcSseg基組在Gaussian 16里（至少截止到目前最新的C.01）還沒有內置，可以從BSE基組數據庫拷定義來使用，參考《詳解Gaussian中混合基組、自定義基組和贗勢基組的輸入》（http://www.shanxitv.org/60）、《在線基組和贗勢數據庫一覽》（http://www.shanxitv.org/309）。pcSseg在ORCA里是內置的，可以直接寫比如pcSseg-1關鍵詞使用pcSseg-1基組。
B97-2在Gaussian里寫作B972。
ORCA里，在CPCM表現的氯仿環境下用DSD-PBEP86/def2-TZVPP做NMR計算的關鍵詞可以這么寫
! DSD-PBEP86 def2-TZVPP def2-TZVPP/C tightSCF NMR miniprint nopop defgrid3 CPCM(chloroform)
雙雜化算NMR時RI近似是強制用的，因此必須用庫侖擬合輔助基組，即這里的def2-TZVPP/C。defgrid3要求使用比默認的defgrid2更好的DFT積分格點和COSX格點，手冊里建議對于雙雜化算NMR用這個。如果要用pcSseg系列基組做雙雜化泛函的NMR計算，沒有標配的庫侖擬合輔助基組可用，但可以用cc-pwCVnZ/C系列，用的pcSseg基組越大則輔助基組應當檔次越高，比如如果是高檔次的pcSseg-3應當用cc-pwCVQZ/C。之所以pcSseg可以用這種輔助基組，在于pcSseg有比較緊的基函數，cc-pwCVQZ/C本來搭配的cc-pwCVQZ為了描述核-價相關也有比較緊的基函數，因此可以搭配pcSseg。若使用autoaux關鍵詞讓ORCA自動構建適合相應pcSseg的輔助基組也可以，但在質量相同的情況下比cc-pwCVQZ/C大不少。對比測試見JCTC, 14, 4756 (2018)的圖1。
ORCA目前自身的代碼不支持KT系列泛函，但可以通過libxc界面使用。例如用KT2算NMR可以這么寫關鍵詞
! def2-TZVPP tightSCF NMR miniprint nopop defgrid3
%method
method dft
functional gga_xc_kt2
end
這里gga_xc_kt2是libxc泛函庫的泛函列表里的名字。orca -libxcfunctionals命令可以顯示支持的所有泛函的寫法。
對于ORCA自身支持的泛函，是否利用libxc界面來實現，實測對NMR計算結果的影響可忽略不計，但我發現用libxc界面計算耗時會高一些。