计算化学公社

标题: 优化长程校正泛函w参数的简便工具optDFTw [打印本页]

作者
Author:
sobereva    时间: 2016-9-21 00:02
标题: 优化长程校正泛函w参数的简便工具optDFTw
:在北京科音高级量子化学培训班(http://www.keinsci.com/workshop/KAQC_content.html)中笔者用50页幻灯片对交换-相关泛函的调控做十分全面和深入的讲解,信息量是本文的5倍以上,非常推荐对这部分内容感兴趣的读者参加!

优化长程校正泛函w参数的简便工具optDFTw
OptDFTw: A utility for optimizing the w parameters of long-range correction functionals

文/Sobereva @北京科音
First release: 2016-Sep-20  Last update: 2024-Nov-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://sobereva.com/soft/optDFTw_1.0_g09.zip
用于Gaussian16版的下载地址:http://sobereva.com/soft/optDFTw_1.0_g16.zip
其中.exe的是Windows版可执行文件,没后缀的是Linux版可执行文件。

若文章使用这两个程序,请引用:Tian Lu, optDFTw program v1.0, webpage: http://sobereva.com/346

这两个程序目前只支持Gaussian程序。只支持中性体系的计算。Windows下运行之前需要在系统中添加GAUSS_EXEDIR环境变量使之指向g09.exe或g16.exe所在目录,并且在PATH环境变量里也添加这个目录,使得能通过命令行方式顺利调用Gaussian。


2 optDFTw程序

此程序涉及两个量,J和J^2,都是衡量N电子态以及N+1电子态时HOMO能量与电离能的偏差之和。注意J^2不是J直接取平方
(, 下载次数 Times of downloads: 219)

(, 下载次数 Times of downloads: 220)

令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,默认的初猜值是上下限的中间值。如果要自己设这些参数,需要以命令行方式运行,即:
optDFTw [下限] [上限] [初猜] [收敛限]

没设的参数会自动用默认值。迭代次数上限是100,如果要改的话需要改源代码里的maxit参数。

下面是实际运行的输出例子(随便选的分子,和上面的示例输入文件不对应),可见每一轮对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。如果自行调节设定,用命令行方式运行:
scanDFTw [下限] [上限] [步长] [iverb]
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


若optDFTw或scanDFTw调用Gaussian计算时意外中断,显然是Gaussian计算出错,自行打开屏幕上提示的Gaussian的out文件,根据末尾的报错信息结合Gaussian使用常识自行判断怎么解决,不要在群里或论坛里问我。>90%的概率是SCF不收敛,参考《解决SCF不收敛问题的方法》(http://sobereva.com/61)试图在模板文件里加帮助收敛的关键词,或者尝试其它基组、几何结构等方式解决。也可能是净电荷/自旋多重度设置不对、Gaussian运行环境有问题、内存不够等原因,自行看报错提示判断。



4 溶剂效应的考虑

为避免一些人犯错,这里特别一提,如果你要在之后的研究中考虑隐式溶剂模型,那么w调控的时候不要直接带着溶剂模型。参考北京科音高级量子化学培训班(http://www.keinsci.com/KAQC)下面这页ppt。在培训里有具体做法的专门的介绍

(, 下载次数 Times of downloads: 106)

作者
Author:
我本是个娃娃    时间: 2016-9-21 08:29
抢沙发,估计日后能用得上
作者
Author:
cgchen    时间: 2016-9-21 08:41
嗯 好物 收藏
作者
Author:
冰释之川    时间: 2016-9-21 09:37
终于可以抛弃自己祖传的定间隔扫描w来确定最佳值了。。。能有optDFTw  确实方便很多,不然 定间距扫描w很耗时……
作者
Author:
LL5312    时间: 2017-5-30 09:10
老师  我想问一下 我调参对两个不同的 但类似的体系  得到了相同的参数 精度是0.01 是巧合 还是 有什么物理背景么   我该如何解释一下 谢谢老师
作者
Author:
LL5312    时间: 2017-5-30 09:18
TTF-TCNQ   Benzene-TCNE 分子间距离都是3.5A
作者
Author:
sobereva    时间: 2017-5-30 10:25
LL5312 发表于 2017-5-30 09:10
老师  我想问一下 我调参对两个不同的 但类似的体系  得到了相同的参数 精度是0.01 是巧合 还是 有什么物理 ...

巧合。
作者
Author:
LL5312    时间: 2017-5-30 10:30
sobereva 发表于 2017-5-30 10:25
巧合。

哦 好的  谢谢
作者
Author:
guoy14iccas    时间: 2018-6-22 18:36
本帖最后由 guoy14iccas 于 2018-6-22 19:17 编辑

有经验师弟给出解决方案是:在PCM下优化w会出现这种情况。改写程序中的优化区间,可以得到解决。
删除不了回复,所以贴出一种解决方案。


作者
Author:
sobereva    时间: 2019-7-12 16:33
添加了for Gaussian 16版
作者
Author:
332544875    时间: 2020-4-18 08:48
本帖最后由 332544875 于 2020-4-18 11:50 编辑

Sob老师您好,我对这个optDFTw工具报告一个BUG,请核实,谢谢。我用Gaussian16优化一个分子后得到输出文件,再分别用Gaussview6和Multiwfn将优化后的结果转换为新的输入文件用于optDFTw工具优化w参数。Gaussview6和Multiwfn转化后的分子坐标相同,但是使用optDFTw工具优化得到的w参数却完全不同。其中使用optDFTw优化Gaussview6转换的输入文件后结果与实验相符,而使用optDFTw优化Multiwfn转化后的转换的输入文件后结果与实验不相符。换成其他分子,所有Multiwfn转换后的gjf文件用于optDFTw优化得到w参数都是在0.05附近,与实验结果不符。附件中包含了转换后的输入文件和optDFTw计算结果。

作者
Author:
sobereva    时间: 2020-4-19 00:37
332544875 发表于 2020-4-18 08:48
Sob老师您好,我对这个optDFTw工具报告一个BUG,请核实,谢谢。我用Gaussian16优化一个分子后得到输出文件 ...

原理上和什么程序做转换无关,gjf里的坐标都是相同的
当前问题应当检查optDFTw中间产生的out看看有什么差异。光从你提供的gjf文件上没法判断原因
作者
Author:
332544875    时间: 2020-4-19 10:04
sobereva 发表于 2020-4-19 00:37
原理上和什么程序做转换无关,gjf里的坐标都是相同的
当前问题应当检查optDFTw中间产生的out看看有什么 ...

谢谢Sob老师指点,我仔细对比了两个程序转换文件的差异,发现问题是出在自选多重度那行。Gaussview6转换后是  0空格1,而Multiwfn转换后是  空格空格0空格空格1。optDFTw工具在调用修改template.gjf文件的自选多重度时没能成功修改Multiwfn转换后的template.gjf文件。我将“空格空格0空格空格1”修改为“0空格1”后就能正常调用修改了。
最后再提一个建议,我在使用optDFTw工具优化参数w过程中,偶尔会遇到SCF不收敛等问题,optDFTw优化中断。Sob老师在后续升级optDFTw工具时能不能加上这个功能:如果SCF不收敛,自动调用template2.gjf模板进行计算(甚至是template3.gjf模板),类似于Molclus程序自动调用不同模板进行优化。
作者
Author:
sobereva    时间: 2020-4-19 20:44
332544875 发表于 2020-4-19 10:04
谢谢Sob老师指点,我仔细对比了两个程序转换文件的差异,发现问题是出在自选多重度那行。Gaussview6转换 ...

你提到的问题我记下了
实际上我是打算之后发布个2.0版,届时有更多功能、支持其它程序,并且写成paper发表。2.0发布时你提到的这些都会实现
作者
Author:
ZHANGZY    时间: 2020-5-4 23:18
求助:
对MN12SX函数使用optdftW时,出现了错误:
Cannot do LC of positive exchange functional values.
Error termination via Lnk1e in /home/zzy/chemsoft/g16/l502.exe at Mon May  4 10:26:58 2020.
求教高手,怎么处理

作者
Author:
sobereva    时间: 2020-5-5 08:03
ZHANGZY 发表于 2020-5-4 23:18
求助:
对MN12SX函数使用optdftW时,出现了错误:
Cannot do LC of positive exchange functional values ...

这是泛函不是函数
此泛函是短程高HF成份而长程HF成分为0,和一般的长程校正泛函明显不同,optDFTw不适合这种情况
作者
Author:
ZHANGZY    时间: 2020-5-7 19:14
多谢,我再试试scanDFTw
作者
Author:
I10140317    时间: 2020-6-28 12:17
sob老师,你好,我看你举的例子是LC-wPBE,我想问下针对一些自定义的LC泛函,如LC-PBE0,程序是否还适用?我推测程序应该也是先把关键词转化成IOP,然后调控IOP里的W值去做优化,如果IOP自定义的话我担心会不会出问题。
作者
Author:
sobereva    时间: 2020-6-29 01:01
I10140317 发表于 2020-6-28 12:17
sob老师,你好,我看你举的例子是LC-wPBE,我想问下针对一些自定义的LC泛函,如LC-PBE0,程序是否还适用?我 ...

程序自动在模板文件的关键词基础上增加对3/107和3/108的设置。
如果对于LC-PBE0你只需要优化w参数,可以用此程序,把下文里提到的LC-PBE0自定义方法中除了IOp(3/107=0300000000,3/108=0300000000) 都写上即可。
在Gaussian中自定义范围分离泛函的方法
http://sobereva.com/550http://bbs.keinsci.com/thread-17256-1-1.html
作者
Author:
wzkchem5    时间: 2020-10-27 10:22
sob老师好,请问可否提供optDFTw的源代码?我最近想把orca和optDFTw结合起来用,如果能有optDFTw的源代码的话,我做一个适用于orca的optDFTw版本应该不麻烦,顺便也能给大家用。如果不方便公开源代码的话,也可以私下给。谢谢!
作者
Author:
让你变成回忆    时间: 2020-10-27 10:56
wzkchem5 发表于 2020-10-27 10:22
sob老师好,请问可否提供optDFTw的源代码?我最近想把orca和optDFTw结合起来用,如果能有optDFTw的源代码的 ...

源代码不是在压缩包里面已经有了吗?
作者
Author:
wzkchem5    时间: 2020-10-27 15:02
让你变成回忆 发表于 2020-10-27 10:56
源代码不是在压缩包里面已经有了吗?

好吧,只在g09的那个压缩包里有,我之前只看了g16的压缩包。
谢谢!
作者
Author:
wzkchem5    时间: 2020-11-1 16:15
我把optDFTw的ORCA版写完了,并且还增加了对于非中性体系,以及自旋多重度不是1的体系的支持(因此除了提供ORCA版的optDFTw以外,我还更新了G09和G16版的optDFTw)。注意对于自旋多重度不是1的体系,程序假定其得到电子和失去电子后自旋多重度都比原来少1,而不是比原来多1。这个假设一般是成立的,但也有例外,如果碰到例外情况可以手动修改源代码解决。

使用方法:
G09和G16:同sob老师主楼的使用方法
ORCA:需要在PATH环境变量中添加ORCA的路径,对于并行计算还需设定LD_LIBRARY_PATH等并行计算需要的环境变量。总之需要设定的参数和单独跑ORCA是相同的。程序调用方式和高斯的情况相同,但模板输入文件应命名为template.inp而非template.gjf。

如果使用中遇到问题,欢迎回复本帖反馈。

(, 下载次数 Times of downloads: 57) (, 下载次数 Times of downloads: 98)


(, 下载次数 Times of downloads: 34)

(, 下载次数 Times of downloads: 54) (, 下载次数 Times of downloads: 26)

(, 下载次数 Times of downloads: 30)




作者
Author:
I10140317    时间: 2021-10-28 19:37
本帖最后由 I10140317 于 2021-10-28 19:39 编辑

sob老师,您好!我想请教一下,如果要考虑溶剂效应(隐式溶剂模型)计算吸收光谱,是应该在基态加上溶剂效应优化w,然后算TDDFT;还是应该在气相优化w,然后在TDDFT增加相应的溶剂关键词?哪种计算方式可能更合理一些?因为优化w的原理是Koopmanns定理,我不太清楚如果用溶剂模型的话,Koopmans定理是不是依旧成立的?
作者
Author:
sobereva    时间: 2021-10-29 06:58
I10140317 发表于 2021-10-28 19:37
sob老师,您好!我想请教一下,如果要考虑溶剂效应(隐式溶剂模型)计算吸收光谱,是应该在基态加上溶剂效 ...

看这篇里的讨论DOI: 10.1002/jcc.26101
作者
Author:
I10140317    时间: 2021-10-29 13:35
本帖最后由 I10140317 于 2021-10-29 16:13 编辑
sobereva 发表于 2021-10-29 06:58
看这篇里的讨论DOI: 10.1002/jcc.26101

sob 老师,谢谢您的回复!我看了您发的这篇文章,文章涉及到优化C0和C1两个参数,但是这两个参数主要是用于调控percent of Long Range HF exchange,也就是beta参数,是仅跟泛函相关的,文章中得到了LC-BLYP等泛函的参数。如果我想用其他该文章中没有拟合的泛函,如LC-wPBE等,那么可不可以通过气相optDFTw优化得到w后,固定w,再用类似optDFTw的算法优化beta参数即可?而不用去管C0和C1的参数。

作者
Author:
张业文    时间: 2021-12-27 20:39
请问老师,调节w值,做单点的初始结构需要什么样的结构,需要基态优化结构去做吗?
作者
Author:
sobereva    时间: 2021-12-28 21:34
张业文 发表于 2021-12-27 20:39
请问老师,调节w值,做单点的初始结构需要什么样的结构,需要基态优化结构去做吗?


作者
Author:
张业文    时间: 2021-12-29 20:26
本帖最后由 张业文 于 2021-12-29 20:27 编辑
sobereva 发表于 2021-12-28 21:34

谢谢sob老师
作者
Author:
函数与激情    时间: 2022-11-11 10:48
卢老师您好,我想计算一个Cl2阴离子自由基的性质(Cl自由基和Cl离子耦合而成),文献中表明这个体系LC-wPBE范围分离泛函描述的挺好,我如果想要优化w参数,对于这种自由基体系调控w的标准还是满足Koopmans定理吗?麻烦卢老师~
作者
Author:
sobereva    时间: 2022-11-11 11:17
函数与激情 发表于 2022-11-11 10:48
卢老师您好,我想计算一个Cl2阴离子自由基的性质(Cl自由基和Cl离子耦合而成),文献中表明这个体系LC-wPBE ...

没必要调控,wB97XD描述得妥妥的。如果只是计算这么小体系,直接用CCSD(T)

作者
Author:
xxxhhh    时间: 2025-5-2 22:19
wzkchem5 发表于 2020-11-1 16:15
我把optDFTw的ORCA版写完了,并且还增加了对于非中性体系,以及自旋多重度不是1的体系的支持(因此除了提供 ...

你好,我已设置路径及环境变量,为何optDFTw-orca依然在试图多核运行时报错?


作者
Author:
wzkchem5    时间: 2025-5-3 20:51
xxxhhh 发表于 2025-5-2 22:19
你好,我已设置路径及环境变量,为何optDFTw-orca依然在试图多核运行时报错?

错误信息已经说得极其清楚了,For parallel runs ORCA has to be "called" with full pathname
注意called这个词。光设对PATH和LD_LIBRARY_PATH是不够的,调用orca时必须用orca的完整路径
作者
Author:
xxxhhh    时间: 2025-5-12 20:30
wzkchem5 发表于 2025-5-3 20:51
错误信息已经说得极其清楚了,For parallel runs ORCA has to be "called" with full pathname
注意call ...

好的好的谢谢,我再研究研究。




欢迎光临 计算化学公社 (http://bbs.keinsci.com/) Powered by Discuz! X3.3