優化長程校正泛函w參數的簡便工具optDFTw
優化長程校正泛函w參數的簡便工具optDFTw
文/Sobereva @北京科音
First release: 2016-Sep-20 Last update: 2019-Jul-12
1 簡介
這幾年優化長程校正泛函的w參數(其實是ω,為打字方便就寫成w)的做法很火,文章接連不斷,有興趣者可以看看鐘成等人的綜述《最優化“調控”區間分離密度泛函理論的研究進展》(DOI: 10.3866/PKU.WHXB201605301)。對w進行調節的一種較好方法是使當前w下計算的HOMO軌道能量盡可能接近電離能。這么做的物理思想是對于精確的交換相關泛函,HOMO能量精確等于電離能,即Koopmans定理完全滿足。這么調節w之后,長程校正泛函可以很好地計算激發能、(超)極化率、fundamental gap、單-三重態激發態能量差等問題(但也并非萬能,比如有大小一致性問題、JPCB,119,1202發現有的體系的超極化率即便w經過優化還是算不好)。w對體系依賴性大,針對一個體系優化的w,對于另一個體系就往往很不適合,所以對每個被計算的體系都需要優化w,導致比普通泛函計算要多花不少代價。
這里提供筆者寫的極其簡單易用的優化w參數的工具optDFTw,以及附帶的對w做掃描的工具scanDFTw。
用于Gaussian09版的下載地址:http://www.shanxitv.org/soft/optDFTw_1.0_g09.zip
用于Gaussian16版的下載地址:http://www.shanxitv.org/soft/optDFTw_1.0_g16.zip
其中.exe的是Windows版可執行文件,沒后綴的是Linux版可執行文件。
若文章使用這兩個程序,請引用:Tian Lu, optDFTw program v1.0, webpage: http://www.shanxitv.org/346。
這兩個程序目前只支持Gaussian程序。只支持中性體系的計算。Windows下運行之前需要在系統中添加GAUSS_EXEDIR環境變量使之指向g09.exe或g16.exe所在目錄,并且在PATH環境變量里也添加這個目錄,使得能通過命令行方式順利調用Gaussian。
2 optDFTw程序
此程序涉及兩個量,J和J^2,都是衡量N電子態以及N+1電子態時HOMO能量與電離能的偏差之和。注意J^2不是J直接取平方

令J或J^2函數最小化就找到了最優的w。由于J^2對w更敏感,所以optDFTw優化的是J^2。這個程序是基于Brent算法來最小化J^2的。優化過程是迭代過程,令w收斂到0.0001就已經足夠精確了,這一般只需要十來圈,如果收斂到0.001一般也就<=十圈。
使用optDFTw程序前首先要編輯一個Gaussian的長程校正泛函的單點任務文件作為模板,存到當前目錄下template.gjf中,比如
%mem=60GB
%nproc=16
# LC-wPBE/6-311+G**
test
0 1
C 0.00000000 0.00000000 -0.52710800
H 0.00000000 0.93885600 -1.11413900
H 0.00000000 -0.93885600 -1.11413900
O 0.00000000 0.00000000 0.67386600
基組可以是混合基組,照常寫即可。泛函可以是比如LC-wPBE、LC-BLYP等標準長程校正泛函,wB97、wB97XD等近程HF成分不為0的泛函可能也能用但結果未必可靠(數據自行負責)。之后在optDFTw運行過程中,就會基于這個文件產生對應不同電子數的N.gjf、N-1.gjf、N+1.gjf,并調用Gaussian進行運算,然后從輸出文件中讀取計算J^2所需的HOMO軌道能量和體系總能量。在Gaussian中用長程校正泛函計算時,將IOp 3/107和3/108都設為MMMMM00000就相當于用了w參數為MMMMM/10000的范圍分離泛函,所以每一步optDFTw都是靠這倆IOp來實現不同w下計算的。
在template.gjf準備好后直接啟動optDFTw就可以進行對w的優化。Brent優化算法需要給定初始的w范圍以及初猜,給得越合適越可能用較少步數收斂,范圍一定要能夠把實際的w值囊括在此范圍中。默認的w下限是0.05(不能寫0,否則Gaussian沒法運行),默認上限是0.6(一般足夠大了),默認的w收斂限是0.0001,默認的初猜值是上下限的中間值。如果要自己設這些參數,需要以命令行方式運行,即:
下面是實際運行的輸出例子(隨便選的分子,和上面的示例輸入文件不對應),可見每一輪對N、N+1、N-1體系分別算一次單點。經過14輪,最終優化的w是0.373547 Bohr^-1,之后在Gaussian中使用此范圍分離泛函時就應當寫IOp(3/107=0373500000,3/108=0373500000)了。
Lower limit: 0.050 Upper limit: 0.600 Init w: 0.325 Tol: 0.00010
The initial point:
Running: g09 N.gjf N.out
Running: g09 N-1.gjf N-1.out
Running: g09 N+1.gjf N+1.out
w: 0.325000 J: 0.0158995406 J^2: 0.0001837031
Iteration: 1
Running: g09 N.gjf N.out
Running: g09 N-1.gjf N-1.out
Running: g09 N+1.gjf N+1.out
w: 0.430041 J: 0.0179941437 J^2: 0.0001839521
Iteration: 2
Running: g09 N.gjf N.out
Running: g09 N-1.gjf N-1.out
Running: g09 N+1.gjf N+1.out
w: 0.219959 J: 0.0661930709 J^2: 0.0026101120
...略
Iteration: 13
Running: g09 N.gjf N.out
Running: g09 N-1.gjf N-1.out
Running: g09 N+1.gjf N+1.out
w: 0.373584 J: 0.0024404394 J^2: 0.0000039539
Iteration: 14
Converged!
The final w: 0.373547 Bohr^-1 J^2: 0.0000039466
3 scanDFTw程序
scanDFTw程序是在指定范圍內按照指定步長對w進行掃描,對每個w會輸出J和J^2值。運行前需要以和optDFTw同樣的方式編寫template.gjf放到當前目錄下。默認從0.05掃到1.0,步長是0.05。如果自行調節設定,用命令行方式運行:
iverb默認為0,如果想同時輸出每個w值的HOMO軌道能量和體系總能量則設為1。比如scanDFTw 0.3 0.5 0.02 1。下面是輸出例子,可見最優的w在w=0.35附近。
w: 0.050000 J: 0.20647538 J^2: 0.02246148
w: 0.100000 J: 0.15601700 J^2: 0.01319407
w: 0.150000 J: 0.11338736 J^2: 0.00721878
w: 0.200000 J: 0.07836289 J^2: 0.00359102
w: 0.250000 J: 0.04967293 J^2: 0.00151705
w: 0.300000 J: 0.02609101 J^2: 0.00045375
w: 0.350000 J: 0.00661196 J^2: 0.00004245
w: 0.400000 J: 0.00958336 J^2: 0.00004770
w: 0.450000 J: 0.02312111 J^2: 0.00031436
w: 0.500000 J: 0.03450162 J^2: 0.00073904
w: 0.550000 J: 0.04413207 J^2: 0.00125304
w: 0.600000 J: 0.05233264 J^2: 0.00181038
w: 0.650000 J: 0.05936051 J^2: 0.00238061
w: 0.700000 J: 0.06541519 J^2: 0.00294496
w: 0.750000 J: 0.07066110 J^2: 0.00349063
w: 0.800000 J: 0.07522949 J^2: 0.00400930
w: 0.850000 J: 0.07921135 J^2: 0.00449682
w: 0.900000 J: 0.08271471 J^2: 0.00495013
w: 0.950000 J: 0.08577732 J^2: 0.00536856