在Multiwfn中分析比CCSD更高級別波函數的方法
在Multiwfn中分析比CCSD更高級別波函數的方法
文/Sobereva @北京科音
First release: 2017-Nov-28 Last update: 2021-May-19
1 前言
Gaussian程序能產生的最高級別的波函數是CCSD,雖然這已經非常精確了,比如在Science, 355, 49 (2017)中作者以CCSD密度為金標準,考察了不同泛函對電子密度的重現性,但是有時候由于特殊原因,就是希望在更高級別下做波函數分析。例如在J. Chem. Theory Comput., 13, 4753 (2017)中,作者對CCSD密度還不知足,于是以CCSDTQ密度為參考考察了不同理論方法對電子密度及衍生量的計算精度;再比如,在J. Chem. Theory Comput., 13, 2705 (2017)中,作者提出了基于自然軌道占據數將動態和靜態相關圖形化展現的方法(此方法在Multiwfn中支持,看3.4.1及之后版本手冊4.A.6節的例子),文中用到了FCI級別自然軌道,以使得對文中所考慮的電子結構復雜的情況,靜態和動態相關都能完美展現。能夠產生比CCSD更高級別波函數的程序較少。PSI4可以產生最高到CCSD(T)級別波函數,MRCC可以產生無窮高階耦合簇和無窮高階CI的波函數(包括FCI。但是產生不了微擾近似版本的耦合簇的波函數)。Multiwfn程序能夠結合這兩個程序在這些比CCSD更高的理論級別做波函數分析,下面就說一下具體方法。本文使用的PSI4為1.2.1版,MRCC為Sept 25, 2017版。Multiwfn為官網上的最新版本,絕對不要用更老版本!
對Multiwfn不了解者請參考《Multiwfn入門tips》(http://www.shanxitv.org/167)、《Multiwfn波函數分析程序的意義、功能與用途》(http://www.shanxitv.org/184)。在閱讀本文之前應當先閱讀此文《詳談Multiwfn支持的輸入文件類型、產生方法以及相互轉換》(http://www.shanxitv.org/379)了解一些基本信息。
2 結合PSI4在CCSD(T)級別下做波函數分析
PSI4網址為http://www.psicode.org。Linux下安裝方法:去官網下載比如psi4conda-1.2.1-py36-Linux-x86_6.sh,然后運行之,按照提示操作即可,會自動下載所有要用的小程序和庫。最后按照提示在用戶目錄下的.bashrc里添加環境變量,如
export PATH=/sob/psi4/bin:$PATH
export PSI_SCRATCH=/sob/psi4scr
下面示例輸入文件用于產生氟化氫的CCSD(T)級別的密度矩陣,并在當前目錄下輸出HF_CCSDpT.fch文件。PSI4默認是不凍核的,如果要求凍核的話是沒法產生CCSD(T)級別的密度矩陣的。
molecule HF {
H 0.0 0.0 -0.831975
F 0.0 0.0 0.092442
}
set basis cc-pVTZ
grad, wfn = gradient('CCSD(T)', return_wfn=True)
fchk_writer = psi4.FCHKWriter(wfn)
fchk_writer.write('HF_CCSDpT.fchk')
注:如果你用的PSI4版本>=1.4,上例最后兩行應當被替換為fchk(wfn,'HF_CCSDpT.fchk')。
將上面輸入文件的內容存到比如test.inp里,然后運行psi4 test.inp test.out,即開始計算并產生輸出文件test.out。
得到的HF_CCSDpT.fchk里面記錄的軌道是參考態HF級別的,而記錄的密度矩陣則是CCSD(T)的。Multiwfn是基于軌道做波函數分析的,從輸入文件里載入的也是軌道而非密度矩陣。因此,若以正常方式令Multiwfn載入這個.fchk,程序只會載入HF軌道,因此之后分析的也都是HF級別的情況。為了能讓Multiwfn分析CCSD(T)波函數,需要做以下步驟:
(1)將HF_CCSDpT.fchk載入Multiwfn
(2)進入主功能200的子功能16
(3)輸入CCSD。此時Multiwfn就會載入此fchk文件里的Total CCSD Density字段的矩陣,由于當前PSI4做的是CCSD(T)計算,所以這個字段對應的是CCSD(T)密度。之后程序立刻輸出了將CCSD(T)密度矩陣對角化得到的自然軌道占據數
(4)輸入y,這將在當前目錄下導出new.mwfn,其中包含了CCSD(T)級別的自然軌道。然后Multiwfn會自動再載入之,此時內存里的密度矩陣、基函數的系數矩陣、GTF信息都是CCSD(T)自然軌道的了。之后任何分析也都是CCSD(T)級別的了
如果.fchk文件是對應開殼層體系,程序還會問你產生哪種類型的自然軌道,可以產生總密度對應的自然軌道、自旋密度對應的自然軌道、alpha/beta各自的自然軌道。更多相關信息可參考《在Multiwfn中基于fch產生自然軌道的方法與激發態波函數、自旋自然軌道分析實例》(http://www.shanxitv.org/403)。
3 結合MRCC在高階耦合簇/CI級別下做波函數分析
MRCC網址為https://www.mrcc.hu,下載預編譯的包,解壓到本機,將目錄添加到PATH環境變量后即可直接使用。輸入文件名必須叫MINP。運行諸如dmrcc |tee out.txt,就會啟動MRCC,讀取當前目錄下的MINP,并將信息同時輸出到屏幕和out.txt文件中。下面的輸入文件內容用于在CCSDT/cc-pVTZ下計算氟化氫單點,并且產生相應級別的密度矩陣(若不寫dens=1則只計算能量)。
basis=cc-pvtz
calc=CCSDT
mem=2500MB
dens=1
geom=xyz
2
H 0.0 0.0 -0.831975
F 0.0 0.0 0.092442
計算完畢后,會在當前目錄下產生MOLDEN文件,這是Molden輸入文件,里面記錄的是參考態HF軌道。還輸出了CCDENSITIES文件,可以用文本編輯器打開,第一部分記錄了二階約化密度矩陣(2RDM),有四個標號;之后記錄的是一階約化密度矩陣(1RDM),有兩個標號(另外兩個標號都是0)。這些密度矩陣以MO為基,標號是從第一個考慮電子相關的MO開始計的。我們也得像上一節那樣,把1RDM轉化成自然軌道存到.molden文件里,才能被Multiwfn利用并在當前級別下做波函數分析。
啟動Multiwfn,輸入MRCC產生的MOLDEN文件路徑,然后進入主功能1000的子功能97,輸入CCDENSITIES文件的路徑,然后輸入有多少個軌道被凍結。MRCC默認是凍核的,當前體系有兩個內核電子(輸出文件中可看到Number of core electrons: 2),當前體系又是閉殼層的,每個占據軌道有倆電子,所以凍結的是1個軌道,因此輸入1。之后Multiwfn會計算自然軌道,然后輸出自然軌道占據數,會看到那個凍核的軌道占據數是精確的2.0,其它都是非整數。如屏幕上的提示所示,Multiwfn已將自然軌道輸出到了一個.mwfn文件里。之后如果輸入y,程序會直接載入這個.mwfn文件,之后做的各種波函數分析就都是CCSDT級別的了。
MRCC做各種CI,以及FCI,也都是把密度矩陣輸出到CCDENSITIES里。基于CI密度進行分析和上述完全相同。MRCC里做FCI很簡單,把calc=后面寫FCI即可。下面是在FCI/aug-cc-pVDZ下計算拉長的LiH的輸入文件,不做凍核近似
basis=aug-cc-pvdz
calc=fci
mem=2500MB
dens=1
core=0
geom=xyz
2
H 0.0 0.0 0.0
Li 0.0 0.0 3.0