计算化学公社

 找回密码 Forget password
 注册 Register
Views: 46131|回复 Reply: 32
打印 Print 上一主题 Last thread 下一主题 Next thread

[综合交流] 计算适用于OPLS-AA力场做模拟的1.2*CM5原子电荷的懒人脚本

[复制链接 Copy URL]

6万

帖子

99

威望

6万

eV
积分
125141

管理员

公社社长

跳转到指定楼层 Go to specific reply
#
计算适用于OPLS-AA力场做模拟的1.2*CM5原子电荷的懒人脚本
A lazy script to calculate 1.2*CM5 atomic charges suitable for OPLS-AA force field simulation

文/Sobereva@北京科音  2021-Jan-28


根据J. Phys. Chem. B., 121, 3864 (2017)中的测试,基于OPLS-AA力场的模拟很适合结合1.2*CM5原子电荷。虽然用1.2*CM5的时候水合自由能的计算精度稍逊于文中测试的另一种原子电荷1.14*CM1A-LBCC(因为这种原子电荷本来就是针对水合自由能训练的校正参数),但在计算其它属性上,如蒸发焓、密度,用1.2*CM5时的误差都更小,而且经验性更低、更有普适性。所以当大家用OPLS-AA力场研究新的小分子的时候,我比较推荐用1.2*CM5电荷。1.2*CM5电荷就是在Truhlar等人提出的原始的CM5电荷基础上乘上1.2,这相当于增大了原子电荷的数量级,等效地体现出了溶剂环境对溶质的极化作用。

笔者之前在《计算RESP原子电荷的超级懒人脚本(一行命令就算出结果)》(http://sobereva.com/476)和《RESP2原子电荷的思想以及在Multiwfn中的计算》(http://sobereva.com/531)中分别给出了超级懒人计算RESP和RESP2电荷的Linux脚本,完全不会用Gaussian的人都能轻松计算,已经有不少人在用。最近碰到完全不会用Gaussian的人试图计算1.2*CM5电荷向我求助,我遂又写了个类似的计算1.2*CM5电荷的懒人Linux脚本,在这里说一下。

这个脚本是Multiwfn文件包里的examples\scripts\1.2CM5.sh,在2021-Jan-28及之后更新的Multiwfn中才有,Multiwfn可以在其主页http://sobereva.com/multiwfn免费下载。

此脚本使用很简单:先确保Gaussian在当前机子里已经装了,见《Gaussian的安装方法及运行时的相关问题》(http://sobereva.com/439),也确保将Multiwfn按照手册2.1.2节的说明在机子里装了。之后假设1.2CM5.sh和一个结构文件phenol.xyz都在当前目录下,只需要运行./1.2CM5.sh phenol.xyz,之后在屏幕上就会看到脚本的运行过程:
Net charge was not defined. Default to 0
Spin multiplicity was not defined. Default to 1
Running optimization task via Gaussian...
Done!
Running formchk...
Running Multiwfn...
Finished! The optimized atomic coordinates with 1.2*CM5 charges (the last column) have been exported to phenol.chg in current folder


最终得到的phenol.chg用文本编辑器打开后可见,其中2、3、4列是优化后的XYZ坐标(埃),最后一列就是1.2*CM5电荷了。

脚本原理:这个脚本会将输入文件里的结构作为初猜结构用Gaussian在B3LYP-D3(BJ)/def2-SVP下做几何优化,然后调用formchk将得到的chk转化为fchk文件,然后调用Multiwfn计算CM5电荷,最后给出1.2*CM5电荷。

几点相关事项:
• 用的输入文件可以是任意含有结构信息而且Multiwfn支持的文件格式,比如pdb/mol/mol2/xyz/fch/gjf/wfn等等等等,见《详谈Multiwfn支持的输入文件类型、产生方法以及相互转换》(http://sobereva.com/379)。
• 默认情况下电荷当成0,自旋多重度当成1。如果是比如不带电的二重态体系,应当写比如./1.2CM5.sh phenol.xyz 0 2。
• 如果运行脚本时提示没可执行权限,先运行chmod +x ./1.2CM5.sh
• 此脚本调用的是Gaussian 16,如果当前机子里装的是诸如Gaussian 09,需要将文件的第22行改为Gaussian=g09。
• 从脚本里的关键词可见,计算是在真空下进行的,这是刻意而为之,不要自行改成在溶剂模型下进行。
• 如果脚本运行失败,有这么几种可能:(1)Gaussian根本没恰当安装,应当确保能手动调用Gaussian运行一个最简单的计算任务 (2)Multiwfn没按照手册2.1.2节恰当安装。其中有些图形库不方便装的话,可以用Multiwfn的noGUI版,此时不需要装图形库 (3)Gaussian计算失败,可能是版本太老而不支持em=GD3BJ关键词,可将脚本中此关键词去掉。也可能几何优化或SCF不收敛,注意看计算中途产生的gau.out文件内容判断原因。
• 如果你手头已经有Gaussian或其它程序产生的Multiwfn支持的波函数文件了(fch/wfn/wfx/molden/mwfn等),就没必要再用这个脚本了。直接载入Multiwfn,依次输入7、16、1就能得到CM5电荷,再手动乘上1.2即可。
• 前述J. Phys. Chem. B.文章中只是考虑了中性体系,没有考虑带电体系。带电体系用什么原子电荷结合OPLS-AA没有统一说法,一般情况下不能用1.2*CM5电荷,因为此时总电荷都不是整数了。笔者建议此时用Multiwfn算RESP2电荷,虽然与OPLS-AA兼容性没有充分测试,但原理上问题不大。
• 如果使用此脚本计算1.2*CM5电荷,请按照Multiwfn程序包里How to cite Multiwfn.pdf文档的说明恰当引用Multiwfn。

如果你没买Gaussian,也可以用免费的ORCA量子化学程序结合Multiwfn算1.2*CM5原子电荷,笔者也提供了相应的傻瓜式脚本,见《ORCA结合Multiwfn计算RESP、RESP2和1.2*CM5原子电荷的懒人脚本》(http://sobereva.com/637)。

评分 Rate

参与人数
Participants 3
eV +15 收起 理由
Reason
laoman + 5 精品内容
zsu007 + 5
朙天儿 + 5

查看全部评分 View all ratings

北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

63

帖子

0

威望

303

eV
积分
366

Level 3 能力者

32#
发表于 Post on 2024-11-15 09:56:14 | 只看该作者 Only view this author
sobereva 发表于 2024-11-15 04:53
1 没有。只能告诉你从原理上这是合适的

2 Na+的原子电荷只能为1.0,这是显而易见的,没什么可“计算” ...

谢谢sob老师

6万

帖子

99

威望

6万

eV
积分
125141

管理员

公社社长

31#
 楼主 Author| 发表于 Post on 2024-11-15 04:53:23 | 只看该作者 Only view this author
xishaofan 发表于 2024-11-14 17:34
sob老师,请问
1)针对带电的十二烷基磺酸根,RESP2电荷计算是更合理,有没有文献支撑呀
2)若十二烷基 ...

1 没有。只能告诉你从原理上这是合适的

2 Na+的原子电荷只能为1.0,这是显而易见的,没什么可“计算”的
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

63

帖子

0

威望

303

eV
积分
366

Level 3 能力者

30#
发表于 Post on 2024-11-14 17:34:25 | 只看该作者 Only view this author
sobereva 发表于 2021-1-29 19:34
笔误改了

Na用+1.2明显说不通。用RESP2明显有意义得多

sob老师,请问
1)针对带电的十二烷基磺酸根,RESP2电荷计算是更合理,有没有文献支撑呀
2)若十二烷基磺酸根用RESP2电荷计算,其抗衡离子钠离子,可以直接使用opls力场文件里面的ions.itp中的钠离子参数吗?如果不可以,也要用RESP2电荷去计算吗?

6万

帖子

99

威望

6万

eV
积分
125141

管理员

公社社长

29#
 楼主 Author| 发表于 Post on 2024-1-22 02:48:13 | 只看该作者 Only view this author
一条君 发表于 2024-1-21 23:41
老师,(1)针对【Na用+1.2明显说不通】延展的一个问题,比如NaPF6,阴离子用了0.8*OPLS-2009IL 的电荷, ...

对一个体系,用诸如1.2、-0.8这样的净电荷,除非有文章证明过对此类体系这种情况的电荷能得到合理结果,否则不要随便用。
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

191

帖子

0

威望

1098

eV
积分
1289

Level 4 (黑子)

28#
发表于 Post on 2024-1-21 23:41:35 | 只看该作者 Only view this author
sobereva 发表于 2021-1-29 19:34
笔误改了

Na用+1.2明显说不通。用RESP2明显有意义得多

老师,(1)针对【Na用+1.2明显说不通】延展的一个问题,比如NaPF6,阴离子用了0.8*OPLS-2009IL 的电荷,也就是阴离子电荷为-0.8,那么对应的钠离子电荷用+0.8是没有问题的吧?
(2)针对【聚丙烯酸钠,要用opls力场——最后建议用1.2*CM5】延展的一个问题,此时聚丙烯酸阴离子是属于带电体系吧,仍然是建议1.2*CM5吗 {3# 3楼的回答1(Na用+1.2...)和回答3区别在哪里呢),谢谢老师
科音学员

327

帖子

2

威望

2811

eV
积分
3178

Level 5 (御坂)

27#
发表于 Post on 2023-11-20 13:06:32 | 只看该作者 Only view this author
sobereva 发表于 2023-11-19 01:21
有机体系基本都可以用OPLS-AA,但极少有非用OPLS-AA不可的时候

谢谢社长

6万

帖子

99

威望

6万

eV
积分
125141

管理员

公社社长

26#
 楼主 Author| 发表于 Post on 2023-11-19 01:21:25 | 只看该作者 Only view this author
七尺贱 发表于 2023-11-18 10:46
**** 本内容被作者隐藏 ****

有机体系基本都可以用OPLS-AA,但极少有非用OPLS-AA不可的时候
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

327

帖子

2

威望

2811

eV
积分
3178

Level 5 (御坂)

25#
发表于 Post on 2023-11-18 10:46:53 | 只看该作者 Only view this author
游客,本帖隐藏的内容需要积分高于 25 才可浏览,您当前积分为 0

513

帖子

0

威望

4386

eV
积分
4899

Level 6 (一方通行)

truffle

24#
发表于 Post on 2021-12-11 13:48:34 | 只看该作者 Only view this author
sobereva 发表于 2021-12-11 04:39
只有拟合静电势电荷在原理上才能实现等价性约束,其它原子电荷都不能,只能之后求平均。
原理上人为平均 ...

谢谢Sob老师解惑~
No problem is insoluble in all conceivable circumstances.

513

帖子

0

威望

4386

eV
积分
4899

Level 6 (一方通行)

truffle

23#
发表于 Post on 2021-12-11 13:48:00 | 只看该作者 Only view this author
喵星大佬 发表于 2021-12-11 00:23
我觉得这两个方法算出来的电荷看上去应该区别不大

不过聚合物的话你要怎么算这两种电荷

我想按Sob老师在17楼的回复先试试,但是我的体系是2-3官能度教练聚合物,可能搞起来比较麻烦
No problem is insoluble in all conceivable circumstances.

6万

帖子

99

威望

6万

eV
积分
125141

管理员

公社社长

22#
 楼主 Author| 发表于 Post on 2021-12-11 04:39:37 | 只看该作者 Only view this author
naoki 发表于 2021-12-10 20:04
Sob老师好,我想请教您一下,我想对聚合物结构单元计算1.2*CM5电荷,但是发现Multiwfn主功能7下的子功能- ...

只有拟合静电势电荷在原理上才能实现等价性约束,其它原子电荷都不能,只能之后求平均。
原理上人为平均不如等价性约束来得理想,但差距不大,可以接受
北京科音自然科学研究中心http://www.keinsci.com)致力于计算化学的发展和传播,长期开办极高质量的各种计算化学类培训:初级量子化学培训班中级量子化学培训班高级量子化学培训班量子化学波函数分析与Multiwfn程序培训班分子动力学与GROMACS培训班CP2K第一性原理计算培训班,内容介绍以及往届资料购买请点击相应链接查看。这些培训是计算化学从零快速入门以及进一步全面系统性提升研究水平的高速路!培训各种常见问题见《北京科音办的培训班FAQ》
欢迎加入北京科音微信公众号获取北京科音培训的最新消息,并避免错过网上有价值的计算化学文章!
欢迎加入人气极高、专业性特别强的理论与计算化学综合交流群思想家公社QQ群(群号见此链接),合计达一万多人。北京科音培训班的学员在群中可申请VIP头衔,提问将得到群主Sobereva的最优先解答。
思想家公社的门口Blog:http://sobereva.com(发布大量原创计算化学相关博文)
Multiwfn主页:http://sobereva.com/multiwfn(十分强大、极为流行的量子化学波函数分析程序)
Google Scholar:https://scholar.google.com/citations?user=tiKE0qkAAAAJ
ResearchGate:https://www.researchgate.net/profile/Tian_Lu

1665

帖子

5

威望

4788

eV
积分
6553

Level 6 (一方通行)

喵星人

21#
发表于 Post on 2021-12-11 00:23:45 | 只看该作者 Only view this author
naoki 发表于 2021-12-10 23:42
大佬好,您觉得对于聚合物OPLS-AA算1.2*CM5电荷手动平均等价原子更合理,还是等价约束算RESP2更合理呢, ...

我觉得这两个方法算出来的电荷看上去应该区别不大

不过聚合物的话你要怎么算这两种电荷

513

帖子

0

威望

4386

eV
积分
4899

Level 6 (一方通行)

truffle

20#
发表于 Post on 2021-12-10 23:42:57 | 只看该作者 Only view this author

大佬好,您觉得对于聚合物OPLS-AA算1.2*CM5电荷手动平均等价原子更合理,还是等价约束算RESP2更合理呢,谢谢~
No problem is insoluble in all conceivable circumstances.

1665

帖子

5

威望

4788

eV
积分
6553

Level 6 (一方通行)

喵星人

19#
发表于 Post on 2021-12-10 22:41:49 | 只看该作者 Only view this author
naoki 发表于 2021-12-10 20:04
Sob老师好,我想请教您一下,我想对聚合物结构单元计算1.2*CM5电荷,但是发现Multiwfn主功能7下的子功能- ...

手动平均等价原子

513

帖子

0

威望

4386

eV
积分
4899

Level 6 (一方通行)

truffle

18#
发表于 Post on 2021-12-10 20:04:21 | 只看该作者 Only view this author
sobereva 发表于 2021-1-29 19:34
笔误改了

Na用+1.2明显说不通。用RESP2明显有意义得多

Sob老师好,我想请教您一下,我想对聚合物结构单元计算1.2*CM5电荷,但是发现Multiwfn主功能7下的子功能-16计算1.2*CM5电荷不像RESP模块可以设置等价性约束,如果想约束的话该如何实现呢,谢谢您!
No problem is insoluble in all conceivable circumstances.

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2026-2-21 05:28 , Processed in 0.186865 second(s), 23 queries , Gzip On.

快速回复 返回顶部 返回列表 Return to list