计算化学公社

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

[GROMACS] 使用Gromacs模拟碳纳米管的一个简单例子

  [复制链接 Copy URL]

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

跳转到指定楼层 Go to specific reply
楼主
:本文的做法如今已经完全过时,如今强烈建议用Sobtop(http://sobereva.com/soft/Sobtop)建立纳米球/管/板的拓扑文件,比x2top灵活得多,参考网页里的例子。另外,北京科音分子动力学与GROMACS培训班(http://www.keinsci.com/workshop/KGMX_content.html)里有模拟纳米管、纳米板的完整讲解和详细的例子。


使用Gromacs模拟碳纳米管的一个简单例子
A simple example of using Gromacs to simulate carbon nanotubes

文/Sobereva @北京科音   2014-Dec-15


很早以前写过一个帖子《amber与gromacs读入碳纳米管的方法》(http://sobereva.com/39)介绍了怎么用gromacs处理碳纳米管,不过貌似还是有很多人不会。这里就提供一个在真空中和在水中模拟碳纳米管的简单且完整的例子。这个小教程本来是给一个国际友人写的,所以都是英文。本文用的是Gromacs 4.6.5,文中的CNT10指的是一个很普通的碳纳米管,cnt10.pdb是其结构文件,这可以用Nanotube Modeler很容易地产生。

文中涉及到的文件在此:
md.mdp (898 Bytes, 下载次数 Times of downloads: 404)

md_vacuum.mdp (444 Bytes, 下载次数 Times of downloads: 328)

cnt10.pdb (22.53 KB, 下载次数 Times of downloads: 413)

CNT10.top (150.83 KB, 下载次数 Times of downloads: 372)


Create a new file named atomname2type.n2t in "share/gromacs/top/gromos54a7.ff" folder in gromacs directory, and then copy below content into it
C    CR1    0.0    12.011         3    C 0.142   C 0.142  C 0.142
C    CR1    0.0    12.011         2    C 0.142   C 0.142
C    CR1    0.0    12.011         1    C 0.142

That means when we use x2top to generate topology file under G54A7 forcefield, if a carbon atom has one or two or three neighbours with a distance about 0.142nm (normal C-C bond length in CNT), then this carbon will be recognized as CR1 atom, which is the atomtype corresponding to the carbon in aromatic CH-group of G54A7. This atomtype is suitable for representing vdW interaction between CNT carbons and environmental atoms.

Run below command to generate the CNT10.top:
g_x2top -f cnt10.pdb -o CNT10.top -ff select -nopbc -name CNT10 -kb 400000 -kt 600 -kp 150
Select G54A7 from the forcefield list, then CNT10.top will be yielded in current folder, and the molecular name is CNT10 (specified by "-name"). -nopbc have to be specified, otherwise x2top can't work normally. -kb, -kt and -kp are used to define the force constant of bond, angle and dihedral terms.

In the CNT10.top, the [ bonds ] field looks like below
    1     3     1  1.420000e-01  4.000000e+05  1.420000e-01  4.000000e+05
    2     3     1  1.420000e-01  4.000000e+05  1.420000e-01  4.000000e+05
    2     4     1  1.420000e-01  4.000000e+05  1.420000e-01  4.000000e+05
    2    27     1  1.420000e-01  4.000000e+05  1.420000e-01  4.000000e+05
    3   218     1  1.420000e-01  4.000000e+05  1.420000e-01  4.000000e+05
...
The 4th column is the equilibrium length determined based on the input geometry (cnt10.pdb), the 5th column is the force constant set by -kb. The last two columns are redundant, you can simply ignore or even directly delete them. The content of [ angles ] and [ dihedrals ] fields are similar to [ bonds ].

The bond, angle and dihedral forcefield parameters we set above are not important, they are mainly used to keep the nearly rigid structure of CNT during simulation, so the force constant can be simply set to a large value; however, too large values will cause high-frequency vibrations and make the simulation unstable. If you want to make the CNT more flexible (rigid), you can decrease (increase) the dihedral force contant. Currently, AFAIK, no forcefield contains bond, angle and dihedral parameters specifically optimized for CNT modelling.


Then there are two cases, you can follow either one

1 MD simulation in vaccum

grompp -f md_vacuum.mdp -c cnt10.pdb -p CNT10.top -o md_vacuum.tpr
mdrun -v -deffnm md_vacuum


2 MD simulation in solvated box

editconf -bt triclinic -f cnt10.pdb -o cnt10_box.gro -d 2
genbox -cp cnt10_box.gro -cs spc216.gro -o cnt10_wat.gro -p CNT10.top

Add below sentence at the head of CNT10.top
#include "gromos54a7.ff/spce.itp"

Start simulation:
grompp -f md.mdp -c cnt10_wat.gro -p CNT10.top -o md.tpr
mdrun -v -deffnm md


模拟和石墨烯、碳球的过程和此文没有任何区别,把文中的cnt10.pdb替换为相应的结构文件即可。石墨烯结构也可以由Nanotube Modeler产生。Nanotube Modeler自带的富勒烯库里有大量碳球结构,各种原子数的各种富勒烯异构体也可以用CaGe产生,CaGe的基本使用方法见《生成富勒烯的螺旋算法简介以及使CaGe中的编号与Fowler-Manolopoulos编号相符的方法》(http://sobereva.com/104

评分 Rate

参与人数
Participants 5
eV +23 收起 理由
Reason
一条君 + 5 好物!
heqiandai + 4 谢谢
Cladoss + 4
hlmkh + 5 <font><font>با &
ter20 + 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

45

帖子

0

威望

365

eV
积分
410

Level 3 能力者

61#
发表于 Post on 2022-10-16 21:03:20 | 只看该作者 Only view this author
谢谢sob老师

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

60#
 楼主 Author| 发表于 Post on 2022-10-16 04:06:43 | 只看该作者 Only view this author
RES 发表于 2022-10-15 15:53
sob老师,我是想要进行一个简单的模拟,是在真空中,对于无氢的开口碳管,想要用分子动力学模拟查看氢原子 ...

至少不适合GROMACS,GROMACS不支持反应力场,没法描述动力学动力学过程中的成键
北京科音自然科学研究中心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

45

帖子

0

威望

365

eV
积分
410

Level 3 能力者

59#
发表于 Post on 2022-10-15 15:53:59 | 只看该作者 Only view this author
sobereva 发表于 2022-10-14 20:07
g_x2top是5.0以前的版本才用的,早已过时
这种体系目前都建议用sobtop建立拓扑文件,看http://sobereva. ...

sob老师,我是想要进行一个简单的模拟,是在真空中,对于无氢的开口碳管,想要用分子动力学模拟查看氢原子或者离子更易于跟那个位置的C成键,老师您看我这种思路适合用分子动力学来模拟吗?

6

帖子

0

威望

57

eV
积分
63

Level 2 能力者

58#
发表于 Post on 2022-10-15 13:14:22 | 只看该作者 Only view this author
sobereva 发表于 2022-10-15 00:31
gmx editconf结合-translate可以平移体系

好的,谢谢sob老师

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

57#
 楼主 Author| 发表于 Post on 2022-10-15 00:31:22 | 只看该作者 Only view this author
YPL 发表于 2022-10-14 22:41
再请教一下,使用Nanotube Modeler软件建立的碳纳米管默认坐标位置就带负值,请问这个怎么更改呢?

gmx editconf结合-translate可以平移体系
北京科音自然科学研究中心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

6

帖子

0

威望

57

eV
积分
63

Level 2 能力者

56#
发表于 Post on 2022-10-14 22:41:19 | 只看该作者 Only view this author
sobereva 发表于 2022-10-14 22:16
搞清楚PBC是什么自然就明白了
不想分成四部分的话,一开始的坐标不要出现负值,否则会被挪到盒子另一头

再请教一下,使用Nanotube Modeler软件建立的碳纳米管默认坐标位置就带负值,请问这个怎么更改呢?

6

帖子

0

威望

57

eV
积分
63

Level 2 能力者

55#
发表于 Post on 2022-10-14 22:18:06 | 只看该作者 Only view this author
好的,谢谢sob老师

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

54#
 楼主 Author| 发表于 Post on 2022-10-14 22:16:16 | 只看该作者 Only view this author
YPL 发表于 2022-10-14 22:11
我想请问一下,建立好碳纳米管后,里面塞入水分子。做能量最小化的时候,碳纳米管分裂成四块,这个是什么问 ...

搞清楚PBC是什么自然就明白了
不想分成四部分的话,一开始的坐标不要出现负值,否则会被挪到盒子另一头
北京科音自然科学研究中心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

6

帖子

0

威望

57

eV
积分
63

Level 2 能力者

53#
发表于 Post on 2022-10-14 22:11:20 | 只看该作者 Only view this author
我想请问一下,建立好碳纳米管后,里面塞入水分子。做能量最小化的时候,碳纳米管分裂成四块,这个是什么问题呢?

202210142210533973..png (17.2 KB, 下载次数 Times of downloads: 96)

202210142210533973..png

202210142211115011..png (33.63 KB, 下载次数 Times of downloads: 98)

202210142211115011..png

45

帖子

0

威望

365

eV
积分
410

Level 3 能力者

52#
发表于 Post on 2022-10-14 21:12:33 | 只看该作者 Only view this author
sobereva 发表于 2022-10-14 20:07
g_x2top是5.0以前的版本才用的,早已过时
这种体系目前都建议用sobtop建立拓扑文件,看http://sobereva. ...

谢谢sob老师

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

51#
 楼主 Author| 发表于 Post on 2022-10-14 20:07:55 | 只看该作者 Only view this author
RES 发表于 2022-10-14 15:33
sob老师,我第一次尝试做分子模拟,对于整个模拟过程还是很模糊,想跟您请教,我输入第一个命令,为什么显 ...

g_x2top是5.0以前的版本才用的,早已过时
这种体系目前都建议用sobtop建立拓扑文件,看http://sobereva.com/soft/Sobtop里面的例子
北京科音自然科学研究中心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

45

帖子

0

威望

365

eV
积分
410

Level 3 能力者

50#
发表于 Post on 2022-10-14 15:33:16 | 只看该作者 Only view this author
sob老师,我第一次尝试做分子模拟,对于整个模拟过程还是很模糊,想跟您请教,我输入第一个命令,为什么显示的是没有g_x2top这个命令,请问这是为什么,是我Gromacs装错了吗?但是我输入gxm可以正常显示版本。

6万

帖子

99

威望

5万

eV
积分
120134

管理员

公社社长

49#
 楼主 Author| 发表于 Post on 2022-9-23 05:18:27 | 只看该作者 Only view this author
吧唧爱吃糖 发表于 2022-9-22 22:06
老师,这个是我的结构文件,我在min.mdp文件中按照您的知道把EnerPres项注释掉了,但是问题还是没有解决 ...

除非文件搞错了,否则不可能提示47091号原子
注意拓扑文件和结构文件的对应关系
而且用这么长的管很莫名其妙
北京科音自然科学研究中心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

54

帖子

0

威望

131

eV
积分
185

Level 3 能力者

48#
发表于 Post on 2022-9-22 22:06:16 | 只看该作者 Only view this author
本帖最后由 吧唧爱吃糖 于 2022-9-22 22:11 编辑
sobereva 发表于 2022-9-21 22:42
没结构文件不好说
肉眼看上去,垂直于管的方向盒子尺寸太小了。当前不该用EnerPres

老师,这个是我的结构文件,我在min.mdp文件中按照您的知道把EnerPres项注释掉了,但是问题还是没有解决。而且输出力原子编号出现了47091,但是我Gro文件里没有这么多原子啊。。。
在SWNT.mdp文件中注释掉EnerPres后问题也没有解决

CNT.zip

1.07 MB, 下载次数 Times of downloads: 3

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

GMT+8, 2025-8-15 09:28 , Processed in 0.194660 second(s), 25 queries , Gzip On.

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