请选择 进入手机版 | 继续访问电脑版

计算化学公社

 找回密码
 现在注册!
查看: 575|回复: 20

[综合交流] 自由能计算对优化与单点级别的敏感性测试

[复制链接]

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

发表于 2017-8-27 21:29:15 | 显示全部楼层 |阅读模式
本帖最后由 liyuanhe211 于 2017-9-2 10:26 编辑

之前在群里讨论通常所做 低级别优化 + 高级别单点 的合理性问题,一时不知如何解释,遂当时提交了一系列计算验证此事。这期间 Sobereva 老师已发布讨论相同问题的帖子《浅谈为什么优化和振动分析不需要用大基组》(http://bbs.keinsci.com/forum.php?mod=viewthread&tid=6600),以不同级别所得键长和频率的差异讨论其影响;本文则增加了从能量角度进行的讨论——即相同单点级别下、使用不同级别优化得到的结构的能量有多大差异,这样两个方面比较的物理量相同、更加直观一些,算是对前文的补充。下面讨论其结果


计算以乙醛烯醇化过程的自由能变为测试对象。
snap0000577.png

测试的几何优化级别自低到高依次包括:
  • PBE0/6-31G
  • HF-3c
  • PBEh-3c
  • RI-BLYP-D3/def2-SV(P)-gCP
  • B3LYP-D3/def2-SV(P)
  • PBE0-D3/def2-SV(P)
  • PBE0-D3/def2-SVP
  • PBE0-D3/def-TZVP
  • PBE0-D3/def2-TZVP
  • APFD/def-TZVP
  • M06-2X-D3/def-TZVP
  • M06-2X-D3/def2-TZVP
  • DSD-PBEP86-D3/def2-TZVP
  • DSD-PBEP86-D3/def2-QZVPP

测试的单点计算级别自低到高包括:
  • B3LYP-D3/def2-SV(P)
  • PBE0-D3/def2-SV(P)
  • PBE0-D3/def2-SVP
  • PBE0-D3/def2-TZVP
  • APFD/def2-TZVP
  • ωB97X-D3/def2-TZVP
  • M06-2X-D3/def-TZVP
  • M06-2X-D3/def2-TZVP
  • RI-B2PLYP-D3/def2-TZVP
  • RI-DSD-PBEP86-D3/def2-TZVP
  • RI-PWPB95-D3/def2-TZVP
  • RI-B2PLYP-D3/def2-QZVPP
  • RI-DSD-PBEP86-D3/def2-QZVPP
  • RI-PWPB95-D3/def2-QZVPP
  • DLPNO-CCSD(T)/cc-pVTZ
  • DLPNO-CCSD(T)/aug-cc-pVTZ
  • DLPNO-CCSD(T)/cc-pVQZ
  • DLPNO-CCSD(T)/aug-cc-pVQZ
  • CCSD(T)/cc-pVTZ
  • CCSD(T)/aug-cc-pVTZ
  • CCSD(T)/cc-pVQZ
  • CCSD(T)/jul-cc-pVQZ


测试结果


  • 几何参数:

原子编号:
dsaffdaf.png


下表中列出了挑选的几个几何参量,包括各个键长,C-C-O键角和一个O-C-H键角:

表:挑选的几个几何参量在各级别下的数值(单位分别为 pm 和 °)

表:与DSD-PBEP86-D3/def2-QZVPP级别的差值

以 DSD-PBEP86-D3/def2-QZVPP 的结果为参照比较可发现:
  • 在一般杂化泛函/def-TZVP 级别以上,几何结构对级别和基组是基本收敛的。
  • 杂化泛函在DZ级别以上的结果也基本可以接受,下文可见其误差对自由能影响不大
  • 纯泛函 BLYP即使加上D3校正,结果也不太好
  • 基于 HF 的 HF-3c 方法、虽然经过"3c"的校正,结果也不能与杂化泛函并论
  • 半经验方法结果凑合,考虑耗时可以说很令人满意了,但其没法明确的算自由能,所以没有参与下文的比较


  • 几何与单点对自由能的影响:

表:各级别的组合算出的Gibbs自由能变原始值(kJ/mol)(可单击放大观看)

经 @ggdh 提议,下面的两个表格列出了各级别组合计算出的自由能与参比级别 CCSD(T)/jul-cc-pVQZ//DSD-PBEP86-D3/def2-QZVPP ,从单点和优化级别两个维度上的差异:

表:单点能级别的影响
(参比第一行:CCSD(T)/jul-cc-pVQZ级别,单击放大)

从单点级别的比较看:
  •     DLPNO-CCSD(T) 与 真·CCSD(T) 的差别很小,小于 0.1 kcal/mol,在 0.3 kJ/mol 左右
  •     双杂化泛函的结果意外的差,双杂化/QZ 级别的结果与CCSD(T)接近,但 双杂化/TZ 级别竟然普遍比普通DFT/TZ级别还差一些!本案例只有一个反应,可能有偶然因素,但这一问题确实在我之前的计算中也发现过,有文章(10.1021/acs.jctc.7b00335)对 DSD-PBEP86-D3 的测试中亦发现这个问题。不知道这与双杂化往往在较高基组下优化参数、而其又对参数比较敏感是否有关系。这一结果更深层次的原因还有待讨论。
  •     M06-2X/TZ 级别做单点的自由能结果甚好,如此之好确实可能有偶然因素,但结合其他benchmark来看,对计算自由能来说其确实是比较好用的普通泛函
  •     仅从此例讲,相对来说,瓜斯吹捧的 APFD 泛函、以及 wB97X-D3 泛函都没有什么亮点,是TZ基组级别结果较差的两个,与 PBE0 配 DZ 基组级别的结果相当、


表:优化级别的影响
(参比第一列:DSD-PBEP86-D3/def2-QZVPP,单击放大)

从优化级别的表格看:
  •     HF-3c 方法虽然加上了校正,但误差仍较大;不给重原子加极化的PBE0/6-31G级别也有明显的误差,且在其他案例中会给出定性都不对的几何结构。均不宜使用。
  •     APFD 泛函相比 M06-2X 和 PBE0 较差

从两个表格的比较可以看出
  •     优化级别带来的影响远小于单点级别的影响(第二个表格的偏差远小于第一个表格)。从标准差看(不严格),单点的敏感程度大致是几何优化敏感程度的 5~15 倍。
  •     对单点来说,H原子不带极化函数的 def2-SV(P) 级别完全不能用,误差约 30kJ/mol(因特意选择了带有氢转移的反应,所以对H原子极化函数敏感);但对优化来说,即使基组不带H的极化函数,影响也不大,在1~2 kJ/mol 的量级。


结论

做二维比较工作量较大,所以本文只考察了一个反映,有较大的偶然性,但下面的结论与日常实践中的感受也比较接近,与大家交流,欢迎讨论:

  • 几何优化在 普通泛函/TZ 级别(如PBE0/def-TZVP)以上是基本收敛的,再提升优化级别对结果的改善有限。
  • 从计算反应自由能的角度,其结果对优化级别的敏感程度远远低于单点级别;故 低级别优化 + 高级别单点 从这个角度讲是很划算的
  • 对涉氢的反应,优化可以不给氢加极化,但单点则必须加
  • DLPNO-CCSD(T)是个好东西
  • 用纯泛函优化有机反应的几何结构不太好、用来算自由能凑合;HF-3c 不太好
  • 似乎存在一些案例表明 DSD-PBEP86-D3 的可靠性、普适性需要进一步探讨(或许在使用时辅以其拟合参数时的基组更安全一些?)
  • 别用瓜斯吹到天上去的 APFD,没看到好处



计算细节


双杂化单点用ORCA完成,其他计算使用 G09 E.01 完成;
几何初猜经微调破坏所有对称性;除 M06-2X 和 APFD 使用 int=ultrafine 之外使用默认积分精度 (即 int=(acc2E=10,fine));
DFT 方法中,有 D3-BJ 参数的用了 empiricaldispersion=GD3BJ,否则用 empiricaldispersion=GD3;
ORCA 的双杂化单点开 RI 近似,用 def/JK,def/C 的组合作为辅助基组,积分格点 Grid6 FinalGrid7,SCF 收敛限设为 tightSCF;
ORCA 的 DLPNO 计算用了 TightPNO,jul-cc-pVQZ 的辅助基组为 aug-cc-pVQZ/C

在有限的精力下,测试仅取一个案例,可能存在一定的巧合因素,应看大致趋势;级别测试的不很全面,主要挑选我平日关心、用得到的级别测试,对其他级别感兴趣的可依据这些信息重复或进行其他测试。








评分

参与人数 12威望 +1 eV +60 收起 理由
啃骨魔 + 5 赞美行动党
杜黎小松 + 5 赞!
zsu007 + 5 赞!
ChrisZheng + 5 GJ!
冰释之川 + 5 好物!
978142355 + 5 好物!
airyang + 5 可视化甚好~
ggdh + 5 GJ!
winterzen + 5 赞!
byymem + 10 谢谢分享
sobereva + 1
ABetaCarw + 5 谢谢分享

查看全部评分

51

帖子

1

威望

553

eV
积分
624

Level 4 (黑子)

那个考了三年没考上的

发表于 2017-8-27 22:26:02 | 显示全部楼层
李老师,我想问个量化以外的问题。那个表格里面的条形图是怎么绘制的?很新奇。
我比你多一点搜主意

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

 楼主| 发表于 2017-8-27 22:43:08 | 显示全部楼层
ABetaCarw 发表于 2017-8-27 22:26
李老师,我想问个量化以外的问题。那个表格里面的条形图是怎么绘制的?很新奇。

就Excel

1万

帖子

16

威望

1万

eV
积分
24039

管理员

公社社长

发表于 2017-8-27 22:59:45 | 显示全部楼层
好文。略有遗憾没把wB97X-D3顺带测了
北京科音自然科学研究中心:http://www.keinsci.com  不定期开办各层次量子化学、分子动力学、Multiwfn程序培训
思想家公社的门口Blog:http://sobereva.com
Multiwfn量子化学波函数分析程序主页:http://sobereva.com/multiwfn
计算化学公社论坛:http://bbs.keinsci.com
思想家公社QQ群1号:18616395,2号:466017436。用于讨论理论、计算化学,两个群讨论内容相同,可加入任意其一但不可都加入,申请信息必须注明研究方向,否则一概不批。

此账号为诸Sobereva共用
Money and papers are rubbish, get a real life!

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

 楼主| 发表于 2017-8-27 23:09:39 | 显示全部楼层
sobereva 发表于 2017-8-27 22:59
好文。略有遗憾没把wB97X-D3顺带测了

机组空出来补一个

336

帖子

6

威望

1034

eV
积分
1490

Level 4 (黑子)

发表于 2017-8-27 23:22:25 | 显示全部楼层
本帖最后由 ggdh 于 2017-8-27 23:28 编辑

。。。之前看错,看到这个可以放心的用b3lyp了,话说为啥没有6-31g(d),可以认为6-31g(d)就是def2-sv(p)等价把?

336

帖子

6

威望

1034

eV
积分
1490

Level 4 (黑子)

发表于 2017-8-27 23:26:51 | 显示全部楼层
liyuanhe211 发表于 2017-8-27 23:25
第二组表格就是为了说明这个问题

恩。。看到这个矩阵了,,这工作量大啊

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

 楼主| 发表于 2017-8-27 23:30:22 | 显示全部楼层
ggdh 发表于 2017-8-27 23:26
恩。。看到这个矩阵了,,这工作量大啊

之前干别的事有现成脚本了,一行一行的做,还好

51

帖子

1

威望

553

eV
积分
624

Level 4 (黑子)

那个考了三年没考上的

发表于 2017-8-28 00:07:02 | 显示全部楼层

条件格式,找到了。
我比你多一点搜主意

7

帖子

0

威望

138

eV
积分
145

Level 2 能力者

发表于 2017-8-28 00:49:20 | 显示全部楼层
以后我也试试!这种计算比较好玩
刚毅坚卓

46

帖子

0

威望

398

eV
积分
444

Level 3 能力者

发表于 2017-8-28 01:54:02 | 显示全部楼层
好的数据可视化之后结论一目了然啊~

话说,例子下的M06-2X-D3/def-TZVP 算撞大运了;

也能看出来3-zeta下,双杂化很可能不如M06-2X-D3+3-zeta稳健,有MP2部分还是乖乖RI-4-zeta走起的好。

44

帖子

0

威望

840

eV
积分
884

Level 4 (黑子)

发表于 2017-8-28 19:59:05 | 显示全部楼层
能不能请李老师把表格赐给我学习一下

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

 楼主| 发表于 2017-8-28 20:06:17 | 显示全部楼层
ChrisZheng 发表于 2017-8-28 19:59
能不能请李老师把表格赐给我学习一下

我正在加其他的级别,现在调乱了,先看图吧。我加完其他级别之后传上来

1904

帖子

20

威望

6145

eV
积分
8449

Level 6 (一方通行)

 楼主| 发表于 2017-8-29 12:26:43 | 显示全部楼层
ggdh 发表于 2017-8-27 23:22
。。。之前看错,看到这个可以放心的用b3lyp了,话说为啥没有6-31g(d),可以认为6-31g(d)就是def2-sv(p)等 ...

BLYP更新了,还是有那么点问题。
平时我完全用def2-SV(P)替代6-31G(d)了

336

帖子

6

威望

1034

eV
积分
1490

Level 4 (黑子)

发表于 2017-8-29 16:23:32 | 显示全部楼层
ΔΔG那个表,都是以CCSD(T)/jul-cc-pVQZ//DSD-PBEP86-D3/def2-QZVPP 为参比。
但是有两个变量(单点方法和优化方法),因此是否可能存在两个方法的误差抵消导致结果虚假好的情况
(但这种误差抵消不能依靠,毕竟只有一个体系)
能不能把ΔΔG表拆开成两个,一个第一列全部弄成0 (只看结构优化的影响)
一个第一行全部弄成0(只看单点的影响)


您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

Archiver|手机版|小黑屋|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949-1号 )

GMT+8, 2017-9-25 03:07 , Processed in 0.115818 second(s), 26 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

快速回复 返回顶部 返回列表