高性价比热力学组合方法G4(MP2)-6X的计算方法

sobereva · 发表于 Post on 2018-5-27 00:04:18

高性价比热力学组合方法G4(MP2)-6X的计算方法

Cost-effective thermodynamics composite method G4(MP2)-6X

文/Sobereva @北京科音

First release: 2018-May-26 Last update: 2019-Jul-30

在2011年，J. Chem. Theory Comput., 7, 112-120中，几个非Gn系列方法的官方人员提出了G4(MP2)-6X热力学组合方法，号称耗时与G4(MP2)相仿佛，而精度很接近昂贵得多的G4。其相对于G4(MP2)最主要变化是把MP2换成了SCS-MP2，把CCSD(T)的CCSD和(T)部分的相关能乘了系数，把优化和振动分析的泛函从B3LYP改为BMK。此方法竟然一直到当下最新的G16 B.01里都仍然没有被加入，好在文章的补充材料里给出了基于Gaussian做此方法计算的Perl脚本，使用很简单，这里简单说一下用法。原作者给的脚本只能输出H(0)和H(T)，没法输出指定温度下的自由能和内能，因此我对脚本进行了一些修改使之能够输出。

首先需要编辑模板.gjf文件，内容如下，也可以直接从这里下：

G4MP2_6x.gjf (686 Bytes, 下载次数 Times of downloads: 93) 。需要将此文件中的坐标、电荷和自旋多重度改成自己分子的情况
%chk=mol.chk
# BMK/6-31+G(2df,p) Opt

A molecule G4(MP2)-6X calculation

0 1
C                0.00000000 0.00000000 -0.56221066
H                0.00000000 -0.92444767 -1.10110537
H                -0.00000000 0.92444767 -1.10110537
O                0.00000000 0.00000000 0.69618930

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read BMK/6-31+G(2df,p) Freq

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read CCSD(T,FrzG4)/GTBas1

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read MP2(FrzG4)/GTMP2LargeXP

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read HF/GFHFB3

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read HF/GFHFB4

用Gaussian运行此脚本，产生比如G4MP2_6x.out。G09和G16经测试都可以用。

下载此Perl脚本：

G4MP2_6x.pl (9.65 KB, 下载次数 Times of downloads: 81) 。把Gaussian输出文件文件和这个.pl文件都拷到Linux下，运行G4MP2_6x.pl G4MP2_6x.out，这个Perl脚本就会自动把相关数据从Gaussian输出文件中提取出来并进行处理，默认是在标况下算的。输出信息例子如下
Temperature (K) 298.15
Pressure (atm) 1
NImag 0
E_ele -114.40450436
H(0K) -114.37806386
H(T) -114.37424742
U(T) -114.37519161
G(T) -114.39906286

Nimag就是虚频数目，E_ele就是电子能量，其它的都不言自明，单位是Hartree。

如果需要计算别的温度和大气压的情况，用文本编辑器打开.pl脚本，修改开头的$temp和$pres即可。

注：后来此方法的作者又提出了G4(MP2)-XK，把G4(MP2)-6X用的Pople基组改为了def2系列，使得此方法可以用于H~Rn的主族体系，对前四周期精度和G4(MP2)-6X相仿佛。在其原文DOI: 10.1021/acs.jctc.9b00449的补充材料里提供了相应的结合Gaussian使用的Perl脚本。

wuzhiyi · 发表于 Post on 2019-7-25 15:35:39

谢谢sob老师的分享有一问题请教我的小分子液相和气相优化结构有较大的差异想问能否用热力学的方式来算？
我的想法是这样

%chk=mol.chk
# BMK/6-31+G(2df,p) Opt SCRF=SMD

A molecule G4(MP2)-6X calculation

0 1
C                0.00000000 0.00000000 -0.56221066
H                0.00000000 -0.92444767 -1.10110537
H                -0.00000000 0.92444767 -1.10110537
O                0.00000000 0.00000000 0.69618930

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read BMK/6-31+G(2df,p) Freq SCRF=SMD

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read CCSD(T,FrzG4)/GTBas1

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read MP2(FrzG4)/GTMP2LargeXP

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read HF/GFHFB3

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read HF/GFHFB4

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read M052X/6-31G*

--Link1--
%chk=mol.chk
# Geom=AllCheck Guess=Read M052X/6-31G* SCRF=SMD

最后加上溶解自由能

sobereva · 发表于 Post on 2019-7-25 21:18:58

wuzhiyi 发表于 2019-7-25 15:35
谢谢sob老师的分享有一问题请教我的小分子液相和气相优化结构有较大的差异想问能否用热力学的方式来算？ ...

不是“热力学的方式”而是“热力学组合方法的方式”

这么做从实际角度上可行，但原理上不十分严格，毕竟热力学组合方法里的参数是对于气相优化的情况拟合的

wuzhiyi · 发表于 Post on 2019-7-30 16:02:02

本帖最后由 wuzhiyi 于 2019-7-30 17:05 编辑

https://pubs.acs.org/doi/10.1021/acs.jctc.9b00449
这里好像有人提出了一个将G4(MP2)-6X拓展到第五行以内的所有的主族元素
作者说在有些测试集表现比原本的G4(MP2)-6X好

wuzhiyi · 发表于 Post on 2019-9-18 19:12:22

想请教社长一个比较细节的问题，对于HLC的计算，作者分成了三种情况考虑，闭壳，开壳，原子和单电子对。

在作者的perl脚本中，为了区分这几种情况，作者考虑了四个变量$na，$mult，$nh和$ex

$nh和$ex是根据$xyz[2]这个变量生产的

在gaussian输出中对应Atomic Type

在我的输出中永远是0.

我想问一下这Atomic Type到底是啥，我翻manual好像是onion的东西，但肯定这里做的是纯QM的，和onion无关。

sobereva · 发表于 Post on 2019-9-19 16:25:48

wuzhiyi 发表于 2019-9-18 19:12
想请教社长一个比较细节的问题，对于HLC的计算，作者分成了三种情况考虑，闭壳，开壳，原子和单电子对。

...

量化计算这个始终都是0，不用管（对于Bq原子是1000）

是ONIOM不是ONION

wuzhiyi · 发表于 Post on 2019-9-19 16:45:07

sobereva 发表于 2019-9-19 16:25
量化计算这个始终都是0，不用管（对于Bq原子是1000）

是ONIOM不是ONION

谢谢社长我发现在perl中$xyz[2]对应的是第二个值不是第三个值
所以这应该对应原子序号

zmm0418 · 发表于 Post on 2019-10-8 17:08:47

请问社长，想用这个方法计算单点能，不知道精度怎么样呀？

sobereva · 发表于 Post on 2019-10-9 08:37:21

zmm0418 发表于 2019-10-8 17:08
请问社长，想用这个方法计算单点能，不知道精度怎么样呀？

精度很理想

1355447622 · 发表于 Post on 2020-5-13 09:26:04

本帖最后由 1355447622 于 2020-5-13 09:29 编辑

sobereva 老师您好，我想用这个方法计算一个单独硫原子的热力学数据，但是算完后我在Liunx里运行脚本时他没有数据，我算的其他分子倒是挺正常的，请问怎么办呢？

C:\Users\admin\Desktop\233[scu@localhost ~]$ cd wzk/[scu@localhost wzk]$ perl G4MP2_6x.pl S.out
Can't take log of 0 at G4MP2_6x.pl line 298, <LOG> line 1749.
S.out [scu@localhost wzk]$

linq5203 · 发表于 Post on 2020-5-13 22:27:17

本帖最后由 linq5203 于 2020-5-13 23:26 编辑

1355447622 发表于 2020-5-13 09:26
sobereva 老师您好，我想用这个方法计算一个单独硫原子的热力学数据，但是算完后我在Liunx里运行脚本时他没 ...

单原子连频率都没有，当然报错。详情可看Shermo手册或自行复习物化（不是阉割版的话）

sobereva · 发表于 Post on 2020-5-14 15:58:54

1355447622 发表于 2020-5-13 09:26
sobereva 老师您好，我想用这个方法计算一个单独硫原子的热力学数据，但是算完后我在Liunx里运行脚本时他没 ...

单个原子虽然能做振动分析（此时只有平动贡献），但没法优化，显然脚本调用Gaussian时会报错。你可以修改脚本，把优化的步骤去掉。

一条君 · 发表于 Post on 2020-7-4 11:36:33

请问老师此方法如何在windows下输出能量呢？谢谢

sobereva · 发表于 Post on 2020-7-5 02:03:49

一条君发表于 2020-7-4 11:36
请问老师此方法如何在windows下输出能量呢？谢谢

何必非要用windows，本来这种耗时高的任务就严重不适合windows版Gaussian
不想装实体linux就用虚拟机就完了

Melvin · 发表于 Post on 2022-5-25 11:47:36

sobereva 发表于 2020-5-14 15:58
单个原子虽然能做振动分析（此时只有平动贡献），但没法优化，显然脚本调用Gaussian时会报错。你可以修改 ...

sobereva老师请问这个“可以修改脚本，把优化的步骤去掉”这个怎么实现呢（比如用的CBS-QB3），是在gjf文件里加什么关键字就可以实现吗？

		自动登录 Automatic login	找回密码 Forget password
密码 Password			注册 Register

[量化理论] 高性价比热力学组合方法G4(MP2)-6X的计算方法

评分 Rate

评分 Rate