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

计算化学公社

 找回密码
 现在注册!
查看: 409|回复: 8

[Multiwfn资源与经验] 使用Multiwfn绘制构象权重平均的光谱

[复制链接]

9767

帖子

16

威望

1万

eV
积分
22374

管理员

公社社长

发表于 2017-6-23 03:59:43 | 显示全部楼层 |阅读模式
使用Multiwfn绘制构象权重平均的光谱

文/Sobereva  2017-Jun-23



1 前言

众所周知,柔性分子有大量构象,并且其中许多低能的构象在当前研究的温度下都是有不可忽略的分布比例的,热平衡状态时它们的比例通常使用Boltzmann公式基于高精度计算的各个构象的相对自由能来计算,见《根据Boltzmann分布计算分子不同构象所占比例》(http://sobereva.com/165)。

不同构象的光谱往往是不同的,甚至有定性差异,特别是圆二色谱的构象敏感性颇高。为了得到能够和实际较相符的理论光谱,我们计算柔性体系光谱时应当对分布>10%(更严格时,分布>5%)的所有构象分别算出的光谱曲线,并按照构象分布比例进行加权线性组合。下面管这种方式得到的理论光谱叫做权重平均光谱或构象平均光谱。

Multiwfn(官网http://sobereva.com/multiwfn)虽然主要用处是波函数分析,但经常关注此程序的人一定知道Multiwfn也具有非常强大的绘制光谱的功能,在《使用Multiwfn绘制红外、拉曼、UV-Vis、ECD和VCD光谱图》(http://sobereva.com/224)中做了介绍,没看过的话一定要先仔细看看。从Multiwfn 3.4版开始,Multiwfn支持直接绘制构象权重平均光谱,这大大简化了绘制这种光谱的操作步骤,本文就通过实例介绍怎么实现。以往绘制这种谱通常过程是把每个构象的光谱曲线用Multiwfn分别导出,然后再一起导进Origin等程序,再手动进行权重组合,当要考虑的构象较多时此过程会相当繁琐。

本文的量化计算使用Gaussian09 D.01完成,Multiwfn使用3.4版,可在其官网免费下载。不熟悉Multiwfn的建议参看《Multiwfn入门tips》(http://sobereva.com/167)。本文的计算除了上述博文外还涉及以下博文的知识,缺乏基础者应当先看看:
《谈谈隐式溶剂模型下溶解自由能和体系自由能的计算》(http://sobereva.com/327
《Gaussian中用TDDFT计算激发态和吸收、荧光、磷光光谱的方法》(http://sobereva.com/314
《谈谈谐振频率校正因子》(http://sobereva.com/221


2 实例:绘制plumericin的构象平均电子光谱

这里通过plumericin(鸡蛋花素)作为例子,演示怎么绘制其构象平均的UV-Vis和ECD光谱,为了简单起见不考虑溶剂效应。plumericin有四种值得考虑的构象,分别用a,b,c,d表示,如下所示
1.png
a与b、c与d的差异在于图右端酯基的方向,a与c、b与d的差异在于O11是在内侧还是外侧。

2.1 准备

对这四个构象,依次做以下计算
(1)在B3LYP/6-31G*下进行优化,得到以下计算用的结构
(2)在B3LYP/6-31G*下做振动分析,计算时带着此级别下的ZPE校正因子0.9806,得到自由能校正量
(3)在M062X/def2TZVP下做单点计算
(4)将(2)的自由能校正量与(3)的单点能相加,得到准确度较好的自由能。然后根据玻尔兹曼公式计算常温下的分布比例
(5)用PBE1PBE/TZVP TD(nstates=20)关键词在PBE0/def-TZVP级别下以TDDFT方法计算最低的20个激发态

第(5)步计算的输出文件都已经提供在了Multiwfn压缩包里的examples\spectra\weighted文件夹中。

计算构象分布比例需要把自由能算得非常准。当前计算级别,只能说能给出个定性正确的分布比例,但要求定量结果较好的话,应当用更高级别计算单点能,比如DSD-PBEP86-D3/def2-QZVP。如果考察溶剂下的情况,一定要带着溶剂模型才行,否则分布比例可能严重错误。当前计算得到的a,b,c,d构象分布比率分别为60.46%、19.50%、16.86%、3.17%。因此,构象平均的UV-Vis或ECD光谱应当计算如下:
构象平均光谱=0.6046*a的光谱+0.1950*b的光谱+0.1686*c的光谱+0.0317*d的光谱


2.2 绘制构象平均的UV-Vis谱

我们写一个叫做multiple.txt的文件(必须叫这个名字),内容如下。里面每一行是上面电子激发计算的输出文件路径,之后是相应构象的出现比率。
examples\spectra\weighted\a.out 0.6046
examples\spectra\weighted\b.out 0.1950
examples\spectra\weighted\c.out 0.1686
examples\spectra\weighted\d.out 0.0317

文件路径写相对路径还是绝对路径都行。这里写的是相对路径,是对于“当前目录”而言的,如果不知道什么叫“当前目录”,看《将文件快速载入Multiwfn程序的几个技巧》(http://sobereva.com/237)。

如手册3.13.2节所示,Multiwfn绘制光谱绝不仅仅支持Gaussian输出文件,比如ORCA、sTDA的输出文件也支持,因此multiple.txt里引入的文件也可以是其它类型的,甚至可以不同来源文件混用。

启动Multiwfn,载入multiple.txt,然后依次输入
11  //绘制光谱的主功能
3  //UV-Vis
0  //显示光谱

马上看到下面的图
2.png
图中红色粗曲线是权重平均后的UV-Vis光谱,可以直接和实验对照。如图例所示,绿、蓝、粉、黑色曲线分别是a,b,c,d四种构象的光谱。从这样的图上,非常便于比较各个构象光谱以及权重平均光谱的差异。由于构象a的比例达到60.5%,占了绝对主导,因此权重平均的光谱整体上也和构象a的比较相似。

用Multiwfn作过单个体系的光谱的人都知道图中黑色竖线对应的是量化程序直接算出来的激发能和强度数据,经过展宽就产生了理论光谱。当前图上黑色竖线是所有构象的竖线的集合,但每个构象对应的竖线的高度都已经乘过了相应的权重值。因此,可以认为当前图中红色粗曲线是由图上这些黑色竖线直接展宽产生的。

Multiwfn还可以把各个构象对权重平均的光谱的贡献曲线直接绘制出来,由此便于讨论实际光谱是如何构成的。关闭刚才的图,然后选一次选项18,然后再选0重新绘制光谱,得到下图
3.png

此图中每个构象各自的光谱曲线高度已经乘上了相应的分布比率,因此,图中四个构象光谱曲线高度之和就是红色粗线;换句话说,图中绿、蓝、粉、黑线分别对应于a,b,c,d构象对权重平均光谱的贡献。从此图上明显看出,a构象的贡献是最显著的。

此图中竖线位置和高度和前一张图一样,但是此图里不同构象对应的竖线用了不同的颜色,和曲线颜色相对应。由于此图中无论是曲线高度还是竖线高度都已经是乘过构象权重值的,所以可以认为每个构象的光谱曲线是由相同颜色的竖线展宽产生的。


2.3 绘制构象平均的ECD谱

所有输入文件同上一节。还是启动Multiwfn,载入multiple.txt,然后依次输入
11  //绘制光谱的主功能
4  //ECD
2  //从a.out载入速度表象的转子强度
2  //从b.out载入速度表象的转子强度
2  //从c.out载入速度表象的转子强度
2  //从d.out载入速度表象的转子强度

为了下面便于讨论,我们先选16让程序把权重平均的ECD谱的极大点和极小点位置都找出来,如下所示
Local maximum X:       181.6653      Value:        24.5974
Local maximum X:       223.9452      Value:        44.9789

Local minimum X:       158.1958      Value:        -0.0081
Local minimum X:       200.9999      Value:       -30.7966
Local minimum X:       253.6460      Value:        -5.2244


然后选0显示光谱,看到下图
4.png

曲线和竖线的含义和上一节相同,不再累述。前面说过圆二色谱对构象比较敏感,确实,当前的ECD图上,各个构象光谱的差异比起它们在UV-Vis中明显更大。相对来说,a和b的谱较相似,c和d的较相似。虽然构象a的出现比例占绝对优势,因此其ECD谱(绿曲线)和权重平均的谱(红曲线)相似程度最高,但是,此图也反映,只考虑这一个构象还是不充分的,因为在<190nm的区域,如果不考虑其它构象的影响,则181.7nm处的峰就难以凸显出来。

关闭上图。为了绘制出各个构象对权重平均ECD谱的贡献,选一次18,再选0,看到下图
5.png
从这个图上可以非常直观地看出权重平均ECD谱的构成。181.7nm的峰很明显是粉色曲线(构象c)所贡献的,因为它几乎和红粗线重合,而构象a,b,d在此处的贡献不仅小,而且由于贡献有正有负,彼此也很大程度抵消了。对201.0nm、253.6nm的负峰和223.9nm正峰,构象a(绿线)显然起到了主导性的贡献,因为其曲线很接近红粗线,构象b对201.0nm和223.9nm处的峰的形成也起到了推波助澜的作用,因为其贡献和构象a的符号是相同的。构象c(粉线)在223.9nm处捣蛋,产生显著负贡献,要是没有它,此处的红粗线的峰高度能明显更高。


3 总结

本文通过很简单的例子说明如何结合主流量化程序和Multiwfn绘制构象权重平均的光谱。可见过程极为简单,写个multiple.txt文件列出各个构象的输入文件以及相应比例即可。而且还可以把所有构象的曲线同时绘制出来,或者把各个构象产生的贡献绘制出来,对讨论权重平均光谱的特征极为有益。Multiwfn还提供了大量选项可以调节绘图效果,可以得到能直接用于发表的图。

另外,在绘图界面里用选项2,可以导出spectrum_curve.txt和curveall.txt到当前目录下。前者记录的是权重平均的光谱曲线数据,而curveall.txt记录的是各个构象的曲线数据,每个构象对应一列。这些文件可以导入Origin等程序来作图,届时可以更灵活地调节图像效果。

评分

参与人数 8eV +40 收起 理由
chemhou + 5
librakitty + 5
captain + 5 赞!
zsu007 + 5 牛!
978142355 + 5 GJ!
小范范1989 + 5 好物!
ggdh + 5 赞!
liyuanhe211 + 5 好物!

查看全部评分

思想家公社的门口Blog:http://sobereva.com
北京科音自然科学研究中心:http://www.keinsci.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!

58

帖子

1

威望

497

eV
积分
575

Level 4 (黑子)

发表于 2017-6-23 04:11:00 | 显示全部楼层
楼主辛苦!

295

帖子

6

威望

892

eV
积分
1307

Level 4 (黑子)

发表于 2017-6-23 07:16:49 | 显示全部楼层
构象平均的带振动分辨的光谱应该也能做把?
另外假设我这个体系只有一个local minimal,我跑个动力学,取1000帧,然后每帧权重设为0.001. 应该也能够用这个方法计算平均光谱?

1

帖子

0

威望

17

eV
积分
18

Level 1 能力者

发表于 2017-6-23 07:31:37 | 显示全部楼层
楼主辛苦!

5

帖子

0

威望

58

eV
积分
63

Level 2 能力者

发表于 2017-6-23 07:45:01 | 显示全部楼层
好贴,楼主辛苦!

411

帖子

0

威望

892

eV
积分
1303

Level 4 (黑子)

发表于 2017-6-23 10:41:51 | 显示全部楼层
sob老师,您说的这一句“用PBE1PBE/TZVP TD(nstates=20)关键词在PBE0/def-TZVP级别下以TDDFT方法计算最低的20个激发态”,意思是PBE1PBE/TZVP=PBE0/def-TZVP吗?

1

帖子

0

威望

33

eV
积分
34

Level 2 能力者

发表于 2017-6-23 14:35:36 | 显示全部楼层
对于ECD计算的数据处理,的确是省事了不少,谢谢sob老师的讲解!

9767

帖子

16

威望

1万

eV
积分
22374

管理员

公社社长

 楼主| 发表于 2017-6-23 15:56:41 | 显示全部楼层
xujian 发表于 2017-6-23 10:41
sob老师,您说的这一句“用PBE1PBE/TZVP TD(nstates=20)关键词在PBE0/def-TZVP级别下以TDDFT方法计算最低的 ...

y
思想家公社的门口Blog:http://sobereva.com
北京科音自然科学研究中心:http://www.keinsci.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!

9767

帖子

16

威望

1万

eV
积分
22374

管理员

公社社长

 楼主| 发表于 2017-6-23 15:59:23 | 显示全部楼层
ggdh 发表于 2017-6-23 07:16
构象平均的带振动分辨的光谱应该也能做把?
另外假设我这个体系只有一个local minimal,我跑个动力学,取1 ...

Multiwfn目前不支持绘制振动分辨的光谱,除非自行把相应的能量、强度数据写成Multiwfn可以认的文本文件,然后再照常绘制。这种情况,也是可以如文中这样多个构象直接混合的。

可以
思想家公社的门口Blog:http://sobereva.com
北京科音自然科学研究中心:http://www.keinsci.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!
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

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

GMT+8, 2017-7-25 18:49 , Processed in 0.101647 second(s), 26 queries .

Powered by Discuz! X3.2

© 2001-2013 Comsenz Inc.

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