• 談談該從Gaussian輸出文件中的什么地方讀電子能量

    談談該從Gaussian輸出文件中的什么地方讀電子能量

    文/Sobereva@北京科音   2019-Jun-1


    1 前言

    在網上答疑時,經常看到有Gaussian初學者問諸如“我要算xxx,是不是用HF能量”、“我應該用SCF Done能量對比么?”。這種問法被內行人看了很煩,完全是根本沒入門、沒常識的問法,回答他本來想問的問題前還得先給對方科普一下從輸出文件里該讀什么數據。鑒于每次都科普一遍太麻煩、甚至都懶得回復,索性這里專門寫個小文章說一下對于常見類型的計算,該在什么地方讀電子能量。但在此之前,筆者覺得有必要先給初學者簡單科普一下什么叫電子能量。


    2 什么叫電子能量、怎么得到

    首先要知道量子化學程序能算出的描述體系的能量的量包括電子能量(electronic energy)以及一些熱力學數據,即內能、焓、自由能。

    一般量子化學書籍、文獻、程序里所謂的“電子能量”包括四項:(1)電子的動能(2)電子與電子間的庫侖互斥能(3)核與核之間的庫侖互斥能(4)電子與核之間的庫侖吸引能。所以一般所說的“電子能量”其實并不完全是它的字面意思,因為也把核之間的互斥能包含了進去。

    電子能量的能量零點是假設核與電子都沒有動能,體系中所有電子和原子核都被分離到無窮遠的情況的能量。這個能量也可以被視為是體系的絕對能量。其數值本身并沒有什么化學意義,只有通過求差來得到與物理/化學上研究的問題對應的量的時候才能體現出意義,如反應能、電離能。除了對微型體系使用極高精度的方法,否則電子能量的計算是不可能達到定量準確的,但由于求差的時候大部分誤差都可以被抵消,因此目前常用的方法給出的有化學意義的數據的精度還是不錯的。

    任何量子化學計算任務都可以給出電子能量。量子化學任務中最簡單的就是單點(Single point)計算任務,指的就是對于一個給定的幾何結構去計算電子能量,這也叫單點能。所以“電子能量”和“單點能”這兩個詞有時候不加區分地使用。如果沒有額外任務設定,量子化學程序默認就是算單點任務。

    單點計算以外的任務也可以給出電子能量,只不過它們的目的不是專門算電子能量,電子能量只是作為一個副產物順帶產生的。比如做幾何結構優化,最后一次輸出的電子能量對應的就是優化后的幾何結構下的電子能量。如果你用優化后的結構手動再做一次單點能計算任務,若此時用的計算級別以及各種影響能量的設定都與做幾何優化時相同,那么此時得到的電子能量也和之前做優化時最后一次輸出的相同(但也有不同的可能,因為Gaussian的幾何優化每一步用的波函數初猜都是上一步收斂的波函數。直接算單點任務時,自動產生的初猜不一定總能收斂到與優化任務最后一步相同的波函數)。

    做振動分析、NMR計算、極化率計算等等考察基態特征的任務,順帶給出的電子能量和計算單點任務給出的是完全相同的(前提是計算條件完全相同)。

    如果想計算內能、焓、自由能這些熱力學上的能量,必須要做振動分析才能得到(或者做熱力學組合方法計算,如CBS-QB3)。它們與電子能量的一個關鍵差異就在于這些熱力學量里也體現了核運動的貢獻,而我們算電子能量的時候是完全忽略掉核的運動的。一般說的零點能(ZPE)對應的是0K下體系的內能/焓/自由能(此時這三者數值相同)與電子能量的差值,來自于原子核的振動運動。顯然ZPE也必須通過做振動分析才能得到,電子能量里是不包含ZPE的。

    關于上述熱力學量的具體計算方法,看《使用Shermo結合量子化學程序方便地計算分子的各種熱力學數據》(http://www.shanxitv.org/552)里面的例子。在文中介紹的Shermo程序手冊的附錄部分有詳細的熱力學數據計算的基礎知識介紹,仔細看了就能徹底理清關系了。


    3 單點任務該在哪里讀電子能量?

    老有初學者搞不清楚單點任務的輸出文件中哪項才是自己需要的當前級別下的電子能量,這里專門說一下,這和你用的理論方法有直接關系。Gaussian輸出文件里給出電子能量的時候單位都是Hartree,即原子單位(a.u.)制下的能量單位。不同版本的輸出信息可能有異,下面的說法至少對于Gaussian 16 A.03是適用的。

    3.1 Hartree-Fock

    Hartree-Fock(HF)方法需要做SCF(自洽場)迭代,迭代過程中電子能量會一點一點逐漸逼近HF方法的最終結果,迭代收斂后會顯示SCF Done,這后面的值顯然就是這個結構下HF級別下給出的電子能量了。

    另外,在輸出文件末尾會有一大堆密密麻麻的輸出,這叫archive段落,是當前計算過程中重要信息的匯總。其中HF=后面的數值和之前輸出的SCF Done后面顯示的完全相同,因此你直接從輸出文件末尾來讀電子能量也可以。

    PS:HF已經徹底過時了,這年頭還用HF絕對文章發不出去,絕對沒有任何理由使用HF方法。這點我在《簡談量子化學計算中DFT泛函的選擇》(http://www.shanxitv.org/272)里專門強調了。


    3.2 普通泛函的DFT計算

    對于雙雜化泛函以外的DFT泛函的計算,如B3LYP、M06-2X、PBE0等,讀電子能量也是讀SCF Done。

    注意,在末尾的archive段落中輸出電子能量時此時還是用HF=這個標簽來輸出,即此時HF=后面的值其實是DFT的結果,和SCF Done后面的值是一樣的,而非指的是Hartree-Fock的結果。這是為什么我特別反感有人用諸如“是不是用HF能量”這種問法,誰知道他是真的做了巨垃圾的Hartree-Fock計算,還是指他做的是DFT計算但是指的實際上是HF=后面的那個DFT結果?


    3.3 MP2

    MP2是后HF方法,計算分為兩個過程,先做常規的HF計算,然后基于HF波函數計算MP2能量。因此,SCF Done,以及archive段落的HF=后面的值,指的都是HF級別的電子能量。MP2級別的電子能量是輸出文件里EUMP2 =后面的值,這個值和archive段落的MP2=后面的值是相同的。


    3.4 雙雜化泛函

    雙雜化泛函里面包含了基于DFT軌道按照MP2公式計算的一部分能量,因此雙雜化泛函也分為兩步,第一步和普通泛函DFT計算一樣要進行SCF迭代,然后輸出SCF Done,這個值和archive段落的HF=后面的值是相同的。然后程序會基于收斂的DFT軌道,再做類似MP2的那一部分。比如用的是B2PLYP雙雜化泛函,則最終雙雜化泛函的能量就是輸出文件中E(B2PLYP) =后面的值,這個值在archive段落里也有,但是標簽用的是MP2=。

    因此,對于雙雜化泛函計算,最省事的讀取電子能量的做法就是直接從archive段落里讀MP2=后頭的值。有的初學者誤將SCF Done后面的值當做雙雜化泛函能量,此時相當于花了雙雜化泛函的耗時,然而得到的卻是一個垃圾數據。

    另外,老有初學者非要用GaussView的Results - Summary界面來讀取結果,我強烈不建議大家這樣做,因為GaussView從Gaussian輸出文件中解析出來的數據可能是有極度誤導性的!!!比如用B2PLYP泛函計算,Summary界面里就給出一個能量數據E(RB2PLYP),這個數據是前面提到的SCF Done后面的數據,根本不是B2PLYP級別的真正的電子能量!這個問題起碼對于GaussView 5.0.9和6.0.16都有。顯然,初學者盲目相信Summary這個很小兒科界面里顯示的數據很容易被坑得頭破血流,而只有從輸出文件里自己讀恰當的數據,才能萬無一失。


    3.5 CCSD(T)

    CCSD(T)是如今最常用的高精度后HF方法。和MP2一樣也是先做HF的SCF迭代,然后基于HF軌道再計算CCSD(T)的電子相關部分,最后給出CCSD(T)電子能量。CCSD(T)級別的電子能量是輸出文件中CCSD(T)=后面的值,直接從archive段落里讀取CCSD(T)=后面的值也是一樣的。

    CCSD(T)是高級別后HF方法,計算過程中也會利用中間數據順帶著把一些更低級別的后HF能量給算出來。因此在archive段落中你也會看到MP2能量、不同形式的MP4能量、CCSD能量。


    3.6 CASSCF方法

    CASSCF是個SCF迭代過程,用Gaussian做CASSCF計算在迭代收斂后,會輸出各個被計算的電子態的能量和組態系數。可以直接搜EIGENVALUES AND  EIGENVECTORS OF CI MATRIX,下面的EIGENVALUE后面的值就是相應的電子態的電子能量了。archive段落的HF=后面的值也是CASSCF的電子能量,如果同時算了多個態的話,這個值是序號最高的那個態的值。


    3.7 TDDFT計算

    TDDFT的計算在《Gaussian中用TDDFT計算激發態和吸收、熒光、磷光光譜的方法》(http://www.shanxitv.org/314)有非常詳細的介紹,其中也專門說了root的設置,它默認為1。輸出文件中,root設成了第幾激發態,在輸出激發態信息的段落中,第幾激發態下面就會有一行Total Energy, E(TD-HF/TD-DFT) =,后面的值就是這個激發態的電子能量,相當于基態電子能量加上這個態的激發能。在archive段落里HF=后面的值是基態的電子能量,當前激發態的電子能量在archive段落里并沒有輸出。


    3.8 半經驗方法

    半經驗方法如AM1、PM7等,計算過程和HF一樣要做SCF迭代,之后輸出的SCF Done的值就是半經驗方法的能量了,在archive段落里依然是以HF=作為標簽來輸出。

    特別要注意的是,半經驗方法在原作者擬合參數的時候大多都是朝著生成焓來擬合的,因此Gaussian做半經驗方法計算給出的能量并不是一般意義的電子能量,而是當前體系在標況下的生成焓。由于生成焓是個相對能量,這也是為什么半經驗方法計算得到的能量比起用前面那些方法的數量級都要小得多得多。


    3.9 分子力學計算

    Gaussian雖然是量子化學程序,但也支持基于AMBER、UFF等分子力場做分子力學計算。從輸出文件中搜Energy=,后面的值就是分子力學能量,在archive段落以HF=標簽來輸出。

    注意分子力學計算給出的能量也不是一般意義的電子能量,分子力學計算時根本都沒有把原子核與電子單獨進行描述,而是每個原子被整個當做一個質點考慮。分子力學的能量零點是每個幾何變量都恰處于分子力場定義的能量最低位置。

    由于半經驗方法和分子力學方法給出的能量都不是體系的絕對能量,雖然你說成“電子能量”別人也知道你指的是什么,但是說成“單點能”或簡單說成“能量”從語義上更合適。


    4 總結

    希望看過本文后,讀者能清楚知道當前用的理論方對應的電子能量應該從哪里讀。也絕對不要再在提問時不給出任何前提就直接說什么“HF能量”、“SCF Done值”這種令內行人看著頭疼、糾結的含糊不清的表達。正確、沒有絲毫歧義的表達應當是比如“雙雜化泛函xxx下算的電子能量”、“PM6方法下算的單點能”。另外,沒有前提時,直接說“能量”的時候,內行人都一律會默認當做你說的是電子能量/單點能,所以如果你要問熱力學量相關問題的時候必須明確說清楚是哪種熱力學量。

    如果你做的是幾何優化任務,輸出文件末尾archive段落輸出的能量數據都是最后一步結構的。

    久久精品国产99久久香蕉