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

计算化学公社

 找回密码
 现在注册!
查看: 4519|回复: 7

[Gaussian/gview] Gaussian中非内置的理论方法和泛函的用法

[复制链接]

1万

帖子

25

威望

1万

eV
积分
35687

管理员

公社社长

发表于 2016-9-4 06:18:48 | 显示全部楼层 |阅读模式
Gaussian中非内置的理论方法和泛函的用法

文/Sobereva @北京科音  2016-Sep-4

本文把Gaussian中没有内置关键词,但能利用IOp或其它手段等效实现的一些理论方法的用法进行介绍,本文对应于Gaussian 09。文中提到的DSD-PBEP86-D3(BJ)在Gaussian 16中已经可以直接使用。


1 微扰相关的方法

2.1 SCS-MP2

这是Grimme在J. Chem. Phys., 118, 9095 (2003)提出的对MP2的改进。MP2相关能可以写为电子对的加和,SCS-MP2把自旋相同电子对的贡献乘上1/3以削弱之,而把自旋相反的贡献乘上6/5以增强之。SCS-MP2算反应能比MP2有了很大提高,达到QCISD级别甚至个别时达到QCISD(T)级别。几何结构也比MP2略有改进。但势垒计算精度比MP2反倒略低一点,弱相互作用计算精度恶化不少。其它方面和MP2差不多。整体性能明显不如目前最好的双杂化泛函,所以如今用处也不太大。

SCS-MP2的使用方法是在MP2计算时额外加上IOp(3/125=0333312000)。3/125=MMMMMNNNNN代表将自旋平行项与自旋反平行项的贡献以MMMMM/10000比NNNNN/10000的比例组合。

MP2计算时从输出中看到
     alpha-alpha T2 =       0.1336200196D-01 E2=     -0.4731575741D-01
     alpha-beta  T2 =       0.8298753242D-01 E2=     -0.3004708246D+00
     beta-beta   T2 =       0.1336200196D-01 E2=     -0.4731575741D-01
SCS-MP2计算时为
     alpha-alpha T2 =       0.1336200196D-01 E2=     -0.1577034194D-01
     alpha-beta  T2 =       0.8298753242D-01 E2=     -0.3605649895D+00
     beta-beta   T2 =       0.1336200196D-01 E2=     -0.1577034194D-01
可见alpha-alpha、beta-beta的能量贡献都被乘上了1/3,而alpha-beta乘上了6/5。


2.2 SOS-MP2

于J. Chem. Phys., 121, 9793(2004)中提出,忽略了自旋平行贡献,反平行部分乘以1.3,结果比MP2烂,特别是弱相互作用。不过号称通过Laplace变换方法可以令标度降到O(N^4)从而使花费比MP2的O(N^5)低。不过再便宜也没普通DFT便宜,MP2精度尚且不及选得最合适的DFT泛函,所以SOS-MP2这个破玩意儿并没什么卵用。大家知道有这么个垃圾就够了。

用法是MP2 IOp(3/125=0000013000)。由于Gaussian并未对之进行优化,所以耗时不会比MP2低。


2.3 SCSN-MP2

在J. Chem. Theory Comput., 3, 80 (2007)中提出,忽略自旋反平行部分,自旋平行部分的系数为1.76,拟合自核酸碱基对相互作用能。弱相互作用计算精度挺不错,误差比MP2低两三倍,有使用价值。另外还有Mol. Phys., 105, 1073 (2007)中提出的SCS(MI)-MP2,MI代表molecular interaction,它与SCSN-MP2的提出目的和精度都很类似,并对不同基组拟合了不同的系数,这里就不提了。从SCSN-MP2和SCS(MI)-MP2的系数都能看出,对于计算弱相互作用的目的,应该削减MP2自旋反平行而增大自旋平行的贡献,这和一般目的的SCS-MP2恰好相反。

用法是MP2 IOp(3/125=1760000000)。


2.4 MP2.5

是在ChemPhysChem, 10, 282-289 (2009)中提出的,将MP2总能量加上乘以0.5的MP3三阶微扰能。计算量和MP3一样。系数0.5来自于分析计算精度、基组依赖性和理论意义。各方面性能比MP2强,特别是计算弱相互作用准确度很高,号称在中等基组下(不加弥散亦可)就能接近CCSD(T)/CBS。

在Gaussian中用MP3计算完后,会看到比如
E2 =    -0.3951023394D+00 EUMP2 =    -0.11430702788798D+03
Keep R2 and R3 ints in memory in symmetry-blocked form, NReq=9641770.
DD1Dir will call FoFMem   1 times, MxPair=        42
NAB=    21 NAA=     0 NBB=     0.
E3=       -0.48480932D-02        EUMP3=      -0.11431187598D+03
MP2.5能量就是-114.30702788798-0.0048480932*0.5=-114.3094519

后来在Phys. Chem. Chem. Phys., 13, 21121 (2011)中还提出了MP2.X,X是给三阶微扰能乘上的系数。MP2.5用在小基组上结果不如在大基组好,为了使得较小基组下就能得到不错的弱相互作用计算结果,MP2.X对从小到大的基组都通过S66测试集拟合了MP3校正能的权重系数,这使得不同基组下(乃至低至6-31G*)得到的弱相互作用能精度都相仿佛,和MP2.5/aug-cc-pVTZ下差不多。虽看似很好,但对更多的体系的可靠性还有待广泛验证。

由于MP2.5、MP2.X能量只能在MP3单点任务结束后根据输出信息手动计算,没法像SCS-MP2那样直接用IOp,所以这两种方法只能算能量,没法用来优化、计算频率等,除非自己去改Gaussian代码实现。


2.5 SCS-MP3

Grimme在J. Comput. Chem., 24, 1529 (2003)中提出的,是SCS-MP2能量加上乘以0.25的三阶微扰校正能。热化学性能比SCS-MP2好很多,往往接近QCISD(T)。但弱相互作用能精度不算理想。

用法是MP3 IOp(3/125=0333312000),然后把输出文件中EUMP2能量加上乘以了0.25的E3。


顺带一提,SCS的思路也被用在CI和耦合簇方法上。SCS-CCSD的正反自旋参数来自于拟合几十个反应能。算反应能,特别是算弱相互作用能比CCSD好很多。Phys. Chem. Chem. Phys., 12, 9611 (2010)中提出的SCS-MI-CCSD则令参数向S22弱相互作用测试集拟合,算弱相互作用极为接近金标准CCSD(T),比MP2.5、SCS-CCSD又明显更好。但可惜Gaussian09中没法实现这些SCS的耦合簇方法,因为程序也不把耦合簇的自旋平行和反平行的贡献独立输出。



2 非内置的普通DFT泛函的用法

2.1 规则

在Gaussian中可以通过IOp在一定程度上自定义泛函。简单来说,若定义IOp(3/76=ab)、IOp(3/77=cd)、IOp(3/78=ef),这里abcdef都已经除以了10000,则所用泛函的实际形式是(下文说的GGA也包含meta-GGA)
E_XC = a*[ d*E_X_LDA + c*ΔE_X_GGA ] + b*E_X_HF + f*E_C_LDA + e*ΔE_C_GGA
每一项的含义如下
E_X_LDA:LDA交换项
ΔE_X_GGA:GGA交换项对LDA交换项的梯度修正部分
E_X_HF:HF交换项。前头的系数b就是常说的HF杂化成份
E_C_LDA:LDA相关项
ΔE_C_GGA:GGA相关项对LDA相关项的梯度修正部分

无论对于交换还是相关部分,都是E_GGA=ΔE_GGA+E_LDA。比如E_X_B88=E_X_LDA+ΔE_X_B88,E_C_LYP=E_C_LDA+ΔE_C_LYP。

对于DFT,a总为1(其实这个参数纯属多余,即便有用也可以直接融合到d、c里嘛)。对于纯泛函,满足d=c=1、f=e=1、b=0;对于单参数杂化泛函,比如PBE0、BHandHLYP,满足d=c、d+b=1、f=e。对于B3LYP这样的三参数杂化泛函才可能d!=c、f!=e。而对于HF,显然b=1,其它皆为0。

用#P在计算时会有诸如这样的输出:
IExCor= 402 DFT=T Ex=B+HF Corr=LYP ExCW=0 ScaHFX= 0.200000
ScaDFX= 0.800000 0.720000 1.000000 0.810000
其中ScaHFX就是b,ScaDFX后面的值是d c f e。Ex就是所用的交换泛函,当前是B88和HF杂化。Corr是所用的相关泛函,当前是LYP。


2.1 B3LYP及变体

老不死的B3LYP定义为
E_XC_B3LYP = (1-a0)*E_X_LDA + a0*E_X_HF + aX*ΔE_X_B88 + E_C_VWN + aC*ΔE_C_LYP
其中a0=0.2、aX=0.72、aC=0.81。

显然,换算成上一节的系数表达方式,就是a=1.0,b=0.2,c=0.72,d=0.8,e=0.81,f=1.0。相应地,计算的时候关键词就写成BLYP IOp(3/76=1000002000,3/77=0720008000,3/78=0810010000)。

B3LYP在提出的时候没有明确说清楚LYP用的VWN到底是VWN3还是VWN5,二者结果有一点差异。Gaussian里默认用的是前者,比后者更好一点点,见Chem. Phys. Lett., 268, 345 (1997)的对比。但有的程序默认用的则是VWN5,如GAMESS-US。如果想用VWN5的B3LYP的话,就可以写BV5LYP IOp(3/76=1000002000,3/77=0720008000,3/78=0810010000),这里V5LYP相关泛函是指把VWN5作为LYP里的LDA部分。

在J. Chem. Phys., 117, 4729 (2002)中提出了B3LYP*,它把B3LYP的HF成份从20%降低到15%,这大大解决了B3LYP高估过渡金属配合物高自旋态稳定性的问题(即低估了高、低自旋态间的能量差)。用B3LYP*就写BLYP IOp(3/76=1000001500,3/77=0720008500,3/78=0810010000),即把X_HF的系数b降到0.15,相应地把X_LDA的系数d提升到0.85,使二者加和仍保持为1。

HF成份对结果有很多系统性的影响,比如HF成份越大TDDFT算的激发能越高。因此可以容易地调节b并相应地修改d,使得算出的光谱和实际尽量接近(虽然洒家鄙夷这种刻意往实验结果上凑的勾当)。


2.2 PBE0及变体

PBE0在Gaussian里写作PBE1PBE,是把PBE泛函中引入25% HF成份,即a=1.0,b=0.25,c=d=0.75,e=f=1.0。相当于PBEPBE IOp(3/76=1000002500,3/77=0750007500)。由于IOp(3/78=1000010000)是默认的,所以可以不用写。

J. Chem. Phys., 138, 021104 (2013)定义了一个无聊的泛函PBE0-1/3,就是把HF成份增加到1/3,故对应于PBEPBE IOp(3/76=1000003333) IOp(3/77=0666706667)。

DFT-D3原文里顺带定义了PBE38,是把HF成份改为3/8=37.5%,测试表明算含频极化率是最好的,对应于PBEPBE IOp(3/76=1000003750) IOp(3/77=0625006250)。


2.3 明尼苏达系列老泛函

Truhlar组之前搞过一大堆乱七八糟的泛函,什么MPW1K、PBE1KCIS,也没多少人用,就他们自己玩得happy。从M05起我看才算真正步入正轨。那些烂七八糟的老明尼苏阿达系列泛函我强烈不推荐使用,早被更好的泛函完爆了,初学者也甭老效仿文献试图去用。如果非要用,去看他们网页http://comp.chem.umn.edu/info/DFT.htm,能通过Gaussian实现的泛函都给了用法。而在Gaussian中即便通过IOp也使不了的,死活非要用就必须找他们要修改版的Gaussian了。


2.4 对长程校正泛函调节ω

Gaussian中对纯泛函可以用LC-前缀做长程校正,比如LC-BLYP。此方法将1/r12算符划分为短程和长程部分,短程不用HF交换项,而长程用100% HF交换项。此做法有一个参数ω很重要,此值越大长程部分涵盖范围就越广,会明显影响结果。LC-加到各种纯泛函上默认是ω=0.47,ωB97X默认是ω=0.3,LC-wPBE默认是0.4。如果要调节这个参数,就把3/107和3/108都设为MMMMM00000,代表ω值为MMMMM/10000。比如要把LC-wPBE的ω设0.3就写LC-wPBE IOp(3/107=0300000000,3/108=0300000000)。



3 非内置的双杂化泛函的用法

2.1 老式双杂化泛函

比较老的双杂化泛函,如B2PLYP、mPW2PLYP、B2GP-PLYP没有考虑色散校正也没有对掺入的二阶微扰能E2以SCS方式考虑。这类老双杂化泛函的交换相关能通式为
E_XC = (1-c1)*E_X_GGA + c1*E_X_HF + (1-c2)*E_C_GGA + c2*E2

最老、最知名,但也是最差的双杂化泛函B2PLYP是J. Chem. Phys., 124, 034108 (2006)中提出的,也是Gaussian09中内置的,定义为
E_XC=0.47*E_X_B88+0.53*E_X_HF+0.73*E_C_LYP+0.27*E2
计算时会看到相应的泛函定义信息,只要看懂了2.1节就自然能理解:
IExCor=  419 DFT=T Ex+Corr=B2PLYP ExCW=0 ScaHFX=  0.530000
ScaDFX=  0.470000  0.470000  0.730000  0.730000 ScalE2=  0.270000  0.270000
这里比2.1节看到的多了ScalE2,后面两个值是给E2的自旋平行和反平行部分分别乘的系数,这里都是0.27,所以说没有以SCS方式考虑E2。

B2GP-PLYP是在B2PLYP之后提出的,原文见J. Phys. Chem. A, 112, 12868 (2008),精度比B2PLYP好不少,但比起后来以SCS方式考虑E2的双杂化泛函还是差多了,其系数c1=0.65,c2=0.36,纯泛函部分还是BLYP,可明确写为
E_XC=0.35*E_X_B88+0.65*E_X_HF+0.64*E_C_LYP+0.36*E2

Gaussian09支持的很可惜只有最老最烂的B2PLYP和后来提出的单精度与之半斤八两的mPW2PLYP。B2GP-PLYP在Gaussian09中没有内置,好在若想用的话可以写
B2PLYP/cc-pVTZ IOp(3/125=0360003600,3/76=1000006500,3/77=0350003500,3/78=0640006400)
这里还写B2PLYP是为了让Gaussian做双杂化计算,但是参数已经被我们改了。3/125=MMMMMNNNNN代表给E2的自旋平行和反平行部分分别乘上MMMMM/10000和NNNNN/10000,其它参数前面已经讲过了。

注:也可以把3/76=1000006500,3/77=0350003500改为3/76=0350006500,3/77=1000010000,根据2.1节的表达式就知道这两种写法是等价的。


2.2 较新的双杂化泛函

J. Phys. Chem. C, 114, 20801 (2010)中提出了DSD-BLYP,它考虑了色散校正,并且以SCS方式考虑了E2的贡献,精度比B2GP-PLYP有一定提升。从此之后各种通式类似的双杂化泛函冒出来不少,比如DSD-PBEP86、DSD-PBEhB95等,此二者都比DSD-BLYP精度又有进一步提升。这些泛函通式写为这样:
E_XC = (1-cX)*E_X_GGA + cX*E_X_HF + cC*E_C_GGA + cO*E2_OS + cS*E2_SS + E_D
这里E2_OS和E2_SS分别是E2中自旋相反和相同部分的贡献,E_D是DFT-D3色散校正能,一般都用BJ阻尼。

J. Comput. Chem., 34, 2327 (2013)是一篇好文章,在其中把各种常见纯泛函组合代入到上面的公式来拟合了系数,包括色散校正参数。经过比较,发现DSD-PBEP86-D3(BJ)精度名列前茅,而且最难得的是能直接在许多主流量化程序里用,包括Gaussian09。在此文的补充材料里可以找到各种纯泛函对应的双杂化泛函的参数,在“DSD-DFT-D3BJ”表格中可见对于DSD-PBEP86-D3(BJ),cX=0.69,cC=0.44,cO=0.52,cS=0.22,DFT-D3(BJ)的参数为s6=0.48,a2=5.6。

下面说怎么在G09中用DSD-PBEP86-D3(BJ)。DFT-D3(BJ)的参数在Gaussian09输入文件里没法直接设,自定义的话必须得设定环境变量。对于Linux版,bash环境,在计算之前先把以下内容输入到控制台中
export GAUSS_DFTD3_S6=480000
export GAUSS_DFTD3_SR6=0
export GAUSS_DFTD3_S8=0
export GAUSS_DFTD3_ABJ1=0
export GAUSS_DFTD3_ABJ2=5600000

之后输入文件里的关键词就写
B2PLYP/cc-pVTZ IOp(3/125=0220005200,3/76=1000006900,3/77=0310003100,3/78=0440004400,3/74=1004) empiricaldispersion=GD3BJ
这里3/74=1004代表把X_GGA设为PBE,C_GGA设为P86。

注1:也可以像此泛函提出者的网页http://www.compchem.me/dsd-pbep86-functional中那样写为IOp(3/125=0220005200,3/78=0440004400,3/76=0310006900,3/74=1004),是等价的。

注2:虽然JCC文中测试表明DSD-PBEhB95-D3(BJ)综合性能最好,但在G09中通过上述方法使用会报错,需要对代码进行特殊修改。

类似地,还可以根据补充材料里的信息使用各种其它的双杂化泛函,但大多不如DSD-PBEP86-D3(BJ)。另外,如果不想用色散校正的版本(比如只有G09 D.01以前版本或者嫌设定环境变量麻烦),即只用DSD-PBEP86,参数是为cX=0.72,cC=0.44,cO=0.51,cS=0.36,故应当写
B2PLYP/cc-pVTZ IOp(3/125=0360005100,3/76=1000007200,3/77=0280002800,3/78=0440004400,3/74=1004)
如文中的数据所示,不用DFT-D3(BJ)的话整体误差会有所增加,但个别项目可能反倒误差会降低(不排除巧合可能)。

评分

参与人数 9eV +46 收起 理由
abdoman + 1 好物!
aqhuangry + 5
chrischen1128 + 5 赞!
cgchen + 5 好物!
卡开发发 + 10 好物!
librakitty + 5 好物!
zsu007 + 5 GJ!
captain + 5 赞!
greatzdk + 5

查看全部评分

北京科音自然科学研究中心: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!

1557

帖子

2

威望

2934

eV
积分
4531

Level 6 (一方通行)

给dalao们倒茶

发表于 2016-9-4 09:03:57 | 显示全部楼层
sob出品,必是精品,顶顶顶,受教了。
淡泊以明志,宁静以致远。

2363

帖子

9

威望

4027

eV
积分
6570

Level 6 (一方通行)

首席卖萌官

发表于 2016-9-4 09:54:50 | 显示全部楼层
板凳,学习受教了
She doesn't love me.
Even so,
my heart has been taken away by her.

28

帖子

0

威望

158

eV
积分
186

Level 3 能力者

发表于 2016-9-4 14:01:32 | 显示全部楼层
学习一个

230

帖子

0

威望

2854

eV
积分
3084

Level 5 (御坂)

发表于 2016-9-4 21:10:56 | 显示全部楼层
用双杂化泛函的时候,需要注意PT2部分到底是用frozen core还是full。

比如Prof. Grimme的双杂化泛函的参数都是采用full拟合的。

有些特定的问题,可能对forzen core或者full比较敏感。

评分

参与人数 1eV +2 收起 理由
卡开发发 + 2

查看全部评分

1万

帖子

25

威望

1万

eV
积分
35687

管理员

公社社长

 楼主| 发表于 2017-6-28 06:17:24 | 显示全部楼层
刚测试了下,G16的DSDPBEP86关键词的结果和此文的方法精确相同(都用int=fine时)
北京科音自然科学研究中心: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!

8

帖子

0

威望

110

eV
积分
118

Level 2 能力者

发表于 2018-8-23 10:26:22 | 显示全部楼层
“在G09 D.01版之前不支持DFT-D3,但是支持了一些标配DFT-D2校正的泛函,即B97D、ωB97XD、B2PLYPD,另外还支持了2008年提出的对弱相互作用很好的M06-2X。靠ωB97XD和M06-2X,基本上就足矣对付各种弱相互作用体系了,二者的性能在伯仲之间,后者略微占优。这两个泛函在G09里分别写为wB97XD和M062X。不过这两个泛函计算速度比起B3LYP这样传统的泛函要慢很多。若想精度更高可以用双杂化结合色散校正的泛函B2PLYPD。”
老师您好,我用的是09A版本,我想考虑色散校正,想用M062X泛函,请问该如何调用?
我是在windows下提交任务到学校的超算中心,请问具体该如何实现色散校正?谢谢

1万

帖子

25

威望

1万

eV
积分
35687

管理员

公社社长

 楼主| 发表于 2018-8-24 17:08:24 | 显示全部楼层
木子之玉ccc 发表于 2018-8-23 10:26
“在G09 D.01版之前不支持DFT-D3,但是支持了一些标配DFT-D2校正的泛函,即B97D、ωB97XD、B2PLYPD,另外还 ...

手册里写得明明白白
M062X不带色散校正,只不过能考虑色散作用

写对应的关键词就完了
DFT-D色散校正的使用
http://sobereva.com/210
北京科音自然科学研究中心: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!
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

手机版|北京科音自然科学研究中心|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949-1号 )

GMT+8, 2018-11-15 09:00 , Processed in 0.349742 second(s), 24 queries .

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