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

计算化学公社

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

[GROMACS] 将GROMACS的原子电荷信息读入VMD的方法

[复制链接]

1万

帖子

25

威望

1万

eV
积分
35687

管理员

公社社长

发表于 2017-3-27 07:58:17 | 显示全部楼层 |阅读模式
将GROMACS的原子电荷信息读入VMD的方法

文/Sobereva @北京科音  2017-Mar-27

VMD可以载入GROMACS的gro、trr、xtc文件,但不会载入速度、受力、电荷。以前写过怎么让VMD把GROMACS产生的速度和受力载入并绘制出来的博文:
《使VMD实时显示gromacs轨迹中原子的受力》(http://sobereva.com/36
《使VMD读入Gromacs产生的trr轨迹中速度信息的方法》(http://sobereva.com/117
有很多VMD里面的插件,比如显示偶极矩、计算静电势、绘制红外光谱的插件都依赖于原子电荷,有时候我们还需要编写依赖于原子电荷的脚本、使用依赖于原子电荷的选择语句,因此把GROMACS模拟时用的原子电荷载入VMD也是非常重要的,这里说说怎么做。

首先在这里下载笔者开发的gmxoutchg程序:http://sobereva.com/soft/gmxoutchg_1.0.rar
其中.exe文件是windows版可执行文件,.f90是源文件,可以自行用它编译Linux版。此程序兼容GROMACS 2016.1版,对其它版本兼容性未经测试。

tpr文件里蕴含了模拟所需的一切信息,包括原子电荷。用GROMACS里的dump命令将二进制的tpr文件转化成文本文件dump.txt,执行以下命令:
gmx dump -s test.tpr > dump.txt

然后将dump.txt放到gmxoutchg.exe所在目录下,运行gmxoutchg.exe,程序就会解析其中的数据,在当前目录下产生charge.txt。其中包含体系所有原子的原子电荷,顺序和结构文件里的原子顺序完全一致。其中第一列、第二列、第三列分别是原子电荷、分子序号、原子序号。

将charge.txt放到VMD目录后,在载入对应的结构/轨迹后,可以用以下tcl脚本将原子电荷从charge.txt中读入。

set sel [atomselect top all]
set natom [$sel num]
set rdchg [open "charge.txt" r]
set chglist {}
for {set iatm 0} {$iatm<=[expr $natom-1]} {incr iatm} {
gets $rdchg line
scan $line "%f" chg
lappend chglist $chg
}
$sel set charge $chglist
close $rdchg

如果想验证是否正确载入了,可以用[atomselect top all] get charge命令把所有原子的原子电荷输出出来(当然,也可以把all换成选择语句只输出指定部分的原子电荷)。

也可以将Coloring method设成Charge,直接从颜色上检验是否正确载入了。下面是乙醇盒子,在默认色彩刻度(BWR)下越红代表原子电荷越负,越蓝代表原子电荷越正。
1.png


顺带一提,对于一些体系通过恰当设定显示方式,可以令原子电荷分布一目了然。比如下图的分子,一个显示方式设为Licorice并且把键调细,另一个显示方式是CPK,把圆球调大,用透明材质,用charge来着色,效果挺不错。
2.png


有了原子电荷信息可以直接使用VMD的依赖于原子电荷的插件了,比如Extensions-Visualization-Dipole Moment Watcher可以观看基于原子电荷计算的偶极矩矢量,就是上图的红色箭头。

评分

参与人数 4eV +15 收起 理由
无敌帅超 + 3 牛!
guoy14iccas + 4 牛!
hlmkh + 3
captain + 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!

2363

帖子

9

威望

4027

eV
积分
6570

Level 6 (一方通行)

首席卖萌官

发表于 2017-3-27 08:11:29 | 显示全部楼层
沙发是我的!
She doesn't love me.
Even so,
my heart has been taken away by her.

289

帖子

0

威望

2922

eV
积分
3211

Level 5 (御坂)

发表于 2017-3-27 08:29:57 | 显示全部楼层
非常不错的脚本,显示方式也很nice顶

136

帖子

0

威望

703

eV
积分
839

Level 4 (黑子)

发表于 2017-3-27 08:44:47 | 显示全部楼层
不错,谢谢分享

85

帖子

0

威望

644

eV
积分
729

Level 4 (黑子)

发表于 2017-3-27 14:20:36 | 显示全部楼层
感谢sob老师,亲测GROMACS-4.6.7 可用。

83

帖子

0

威望

1169

eV
积分
1252

Level 4 (黑子)

发表于 2017-3-28 17:11:43 | 显示全部楼层
以前我也用的这种方法

不过前不久在别的群里有人提了个更方便的方法
如果只是想要体系的电荷   gmx可以输出pqr格式文件
gmx editconf   -f a.tpr   -mead   a.pqr  
直接载入vmd就补上了电荷和半径

1万

帖子

25

威望

1万

eV
积分
35687

管理员

公社社长

 楼主| 发表于 2017-3-28 21:41:29 | 显示全部楼层
diaok 发表于 2017-3-28 17:11
以前我也用的这种方法

不过前不久在别的群里有人提了个更方便的方法


这样做的一个缺点是原子电荷保存精度太低,只有两位小数,整个分子所有原子电荷加和会往往能偏离整数不少(特别是进一步通过脚本计算静电作用能的话误差会较大)
北京科音自然科学研究中心: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!

83

帖子

0

威望

1169

eV
积分
1252

Level 4 (黑子)

发表于 2017-3-30 11:09:56 | 显示全部楼层
sobereva 发表于 2017-3-28 21:41
这样做的一个缺点是原子电荷保存精度太低,只有两位小数,整个分子所有原子电荷加和会往往能偏离整数不 ...

的确是这样。。
这种方法只能用来初步可视化电荷的分布了
您需要登录后才可以回帖 登录 | 现在注册!

本版积分规则

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

GMT+8, 2018-11-15 08:57 , Processed in 0.141465 second(s), 29 queries .

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