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

sobereva

高性价比热力学组合方法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)

2018-5-27 00:04 上传 Uploaded

点击下载
Click to download

。需要将此文件中的坐标、电荷和自旋多重度改成自己分子的情况
%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)

2018-5-27 00:03 上传 Uploaded

点击下载
Click to download

。把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

谢谢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

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

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

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

wuzhiyi

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

wuzhiyi

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

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

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

在gaussian输出中对应Atomic Type

在我的输出中永远是0.

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

sobereva

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

...

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

是ONIOM不是ONION

wuzhiyi

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

是ONIOM不是ONION

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

zmm0418

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

sobereva

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

精度很理想

1355447622

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

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

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

sobereva

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

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

一条君

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

sobereva

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

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

Melvin

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

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