使用GROMACS计算径向分布函数

类别:    标签: gmx bash   阅读次数:   版权: (CC) BY-NC-SA

使用GROMACS计算径向分布函数和配位数比较简单, 只要了解一下gmx rdf程序的使用即可. 但即便是简单的工作, 有时候也需要我们稍微深入思考一下, 如何才能做得更快更好, 如何才能尽量减少重复运行所需的时间. 只有这样, 我们才能不断突破自己的心理舒适区, 提高自己. 否则的话, 每次都暴虎冯河, 自己既做得累, 也易生厌. 毕竟, 人, 是要做些有技巧性, 艺术性的东西的.

就以我最近要处理的问题来说明吧. 我跑了一个纤维素二聚体的模拟, 这个分子共有45个原子,


视图: 投影 正交    速度:
模型: 球棍 范德华球 棍状 线框 线型   名称
左键: 转动   滚轮: 缩放   双击: 自动旋转开关   Alt+左键: 移动

Fig.1

跑完后, 要计算二聚体中每个原子与溶剂水中的氧原子OW, 氢原子HW的RDF, 从而确定每个原子的第一溶剂化层厚度.

要算RDF么, 首先要做出索引文件, 把二聚体的每个原子, OW, HW放到不同的组中. 简单. 执行

gmx make_ndx -f conf.gro

再使用下面的指令

a OW: 选择水中的氧原子

a HW: 选择水中的氢原子

splitat 1: 将组1拆分, 每个原子一个组. 组1是other也就是二聚体.

q: 保存索引文件, 退出程序.

这样我们就得到了索引文件, 可用于指定要计算的RDF了.

打开得到的索引文件看看, 各个组的排列顺序如下

[ System ]
[ Other ]
[ ROH ]
[ 4GB ]
[ 0GB ]
[ Water ]
[ SOL ]
[ non-Water ]
[ OW ]
[ HW ]
[ Other_HO1_1 ]
[ Other_O1_2 ]
[ Other_C1_3 ]
...
[ Other_H2_43 ]
[ Other_O2_44 ]
[ Other_H2O_45 ]

使用这个索引文件指定要计算的RDF不方便, 因为二聚体的原子编号和索引组编号不一致, 我们把它们调整到一致(注意, 索引组编号是从0开始的), 调整后, 索引文件顺序如下

[ System ]
[ Other_HO1_1 ]
[ Other_O1_2 ]
[ Other_C1_3 ]
...
[ Other_H2_43 ]
[ Other_O2_44 ]
[ Other_H2O_45 ]
[ OW ]
[ HW ]
[ Other ]
[ ROH ]
[ 4GB ]
[ 0GB ]
[ Water ]
[ SOL ]
[ non-Water ]

这样在计算RDF时, 只要指定1 46就可以算出二聚体的1号原子与OW的RDF了. 类似操作45次, 就能得到45个原子与OW的RDF了. 再类似操作45次, 得到45个与HW的RDF. 完成.

如果你每次都这么做上90次, 还要保证不输错数字, 那还是很能锻炼耐心的. 我却没有这么大的耐心, 所以就只能偷个懒了, 使用bash脚本

for((i=1; i<=45; i++)); do
	echo -e $i'\n46' | gmx rdf -f -n
done

当然, 如果你bash版本支持, 你也可以使用

for i in {1..45}; do
	echo -e $i'\n46' | gmx rdf -f -n
done

如果你还知道parallel, 那就更好了, 你可以并行执行RDF计算, 这会大大节省时间

paralle echo -e {}'\\n46' | gmx rdf -f -n ::: {1..45}

所以, 学习一点简单的bash语句, 你的生活可以更轻松. 学习的过程并不舒服, 可投入一段时间的不舒服, 换取以后更多时间的舒服, 你要如何选择呢?

补充

上面的示例虽然用了gmx 5.x版本, 但做法还是基于4.x版本的. gmx 5.x版本的rdf实际支持一次计算多个组的RDF, 只要连续给出组编号就可以了. 这样得到的输出文件中, 含有多个RDF.

随意赞赏

微信

支付宝
◆本文地址: , 转载请注明◆
◆评论问题: https://jerkwin.herokuapp.com/category/3/博客, 欢迎留言◆


前一篇: 【转载】[在人间]第48期:一所山村小学的日常
后一篇: GMXTOP:OPLSAA力场的GROMACS拓扑文件生成器

访问人次(2015年7月 9日起): | 最后更新: 2017-08-15 19:57:07 UTC | 版权所有 © 2008 - 2017 Jerkwin