- 2017-06-29 21:46:49
- 2017-11-29 21:42:29 修订扩充
- 2017-12-09 09:03:00 集成MKTOP, 增加生成输入文件功能
使用GROMACS进行分子动力学模拟时, 获得体系的拓扑文件是关键, 也是难点. 虽然有些辅助工具可以直接获得一些分子的拓扑文件, 但你仍然需要对得到的拓扑文件进行仔细检查, 否则的话, MD过程中出现问题很难排查.
tppmktop是获得有机分子OPLSAA力场拓扑文件的好工具, 但也存在一定的不足. 一是不支持周期性分子, 二是提供的网络服务同时只能运行一个任务. 为此, 我觉得还是有必要自己做一个简单OPLSAA力场拓扑生成器GMXTOP, 一则可用于理解力场拓扑的生成方法, 二则可用于检查其他方法得到的拓扑是否合适.
GMXTOP是个在线工具, 运行在浏览器中, 不依赖于任何其他环境. 它并不支持自动生成拓扑文件, 而是需要你先指定每个原子的原子类型, 然后它会根据原子类型和力场文件中的参数信息生成GROMACS的拓扑文件. 虽然GMXTOP不是完全自动的, 使用起来有点不方便, 但使用它你可以明确地控制如何指定原子类型, 遇到其他工具无法自动匹配的原子时, 你就可以根据情况选用相近的原子类型, 保证最终能得到合适的拓扑文件.
使用说明
打开GMXTOP网页 https://Jerkwin.github.io/prog/gmxtop.html. 建议使用Chrome, 因为我没有测试其他浏览器. 界面尚未美化, 看起来有点简陋.
分子构型格式
暂时支持读入.pdb
和.mol
格式的分子构型文件, 而且文件中必须包含原子间的连接信息, 因为程序要使用这些连接信息来确定原子间的成键信息.
这两种格式的分子构型文件都可以使用GaussView获得, 程序只测试过GaussView的输出格式, 尚未测试其他分子编辑软件. 鉴于各种软件给出的.pdb
格式不尽相同, 建议优先使用.mol
格式.
指定原子类型
有三种方式指定原子类型:
- 右上方分子结构窗口中点击原子后, 原子会高亮, 同时弹出所有可能匹配的原子类型及其相应的示意图. 图片中对应的原子类型以红色表示. 双击图片即可完成指定. 也可选点击图片前面的按钮,
OK
确认后完成指定. - 如果需要同时指定多个原子的原子类型, 可在左上方文本框内选中多个原子, 然后点击
Assign Atom Type
, 然后指定. 当然, 这种方法也可以用于单个原子. - 使用MKTOP的判定方法自动指认原子类型. 点击
Try AUTO Assign
即可. 由于运行时间稍长, 请耐心. 有关说明见GMXTOP:集成MKTOP的原子类型判定代码
注意, 自动指认原子类型时可能存在无法判定的情况, 必要时附加手动指定. 如果判断原子类型存在困难, 你可以试着先用tppmktop来处理一下, 参考它给出的原子类型. 在大多数情况下, tmmktop给出的原子类型都是正确的.
生成拓扑文件
指定好每个原子的原子类型后, OPLSAA Atom Types
文本框内会列出每个原子的原子类型, 这些原子类型也可以直接进行修改. 确认正确后, 点击Create Topology File
即可生成GROMACS拓扑文件. 点击Save .top File
即可下载生成的拓扑文件.
如果你要使用拓扑文件中的电荷, 一个判断拓扑文件好坏的简单标准就是看整个分子的净电荷qtot
与实际是否相符. 如果基本相符, 那说明原子类型的指认还是比较正确的. 否则的话, 那就需要你仔细检查每个原子的原子类型了.
对一些特殊的原子类型, 它们之间的成键相互作用项可能缺失, 得到的拓扑文中会标识!!! NOT DEFINED !!!
. 在使用拓扑文件前, 你需要将这些缺失项补充完整. 补充的方法一个是拟合, 一个是采用相近原子类型的成键参数.
测试拓扑文件
如果对生成的拓扑文件满意, 可以依次点击Save .gro File
和Save .mdp File
保存运行GROMACS所需的结构文件和参数文件, 然后gmx grompp; gmx mdrun
就可以快速地运行模拟来测试得到的拓扑是否合适了.
其他辅助功能
- 鼠标滚轮: 缩放
Ctrl
+鼠标左键: 平移Labels
: 可显示每个原子的编号, 用于区分Hide H
: 隐藏氢原子Reset View
: 重新将分子居中显示Rot X/Y/Z
: 自动旋转开关
简单示例
MeNO2
MeNO2.pdb | |
---|---|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 | REMARK 1 File created by GaussView 5.0.9
HETATM 1 C 0 -2.245 0.945 -0.011 C
HETATM 2 H 0 -1.888 -0.063 -0.011 H
HETATM 3 H 0 -1.888 1.450 -0.885 H
HETATM 4 H 0 -3.315 0.946 -0.012 H
HETATM 5 N 0 -1.755 1.638 1.189 N
HETATM 6 O 0 -1.548 2.858 1.158 O
HETATM 7 O 0 -1.549 1.002 2.230 O
END
CONECT 1 2 3 4 5
CONECT 2 1
CONECT 3 1
CONECT 4 1
CONECT 5 6 7 1
CONECT 6 5
CONECT 7 5
|
自动指定原子类型即可, 程序判断无误.
PhCNHNH2
PhCNHNH2.pdb | |
---|---|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 | REMARK 1 File created by GaussView 5.0.9
HETATM 1 C 0 -1.372 1.261 1.088 C
HETATM 2 N 0 -0.527 2.236 0.993 N
HETATM 3 H 0 0.072 2.456 1.763 H
HETATM 4 N 0 -1.452 0.477 2.329 N
HETATM 5 H 0 -1.206 1.057 3.106 H
HETATM 6 H 0 -2.386 0.138 2.449 H
HETATM 7 C 0 -2.293 0.922 -0.098 C
HETATM 8 C 0 -3.208 -0.134 0.005 C
HETATM 9 C 0 -2.216 1.669 -1.281 C
HETATM 10 C 0 -4.046 -0.443 -1.075 C
HETATM 11 H 0 -3.267 -0.705 0.908 H
HETATM 12 C 0 -3.054 1.360 -2.361 C
HETATM 13 H 0 -1.517 2.475 -1.360 H
HETATM 14 C 0 -3.969 0.304 -2.258 C
HETATM 15 H 0 -4.745 -1.250 -0.996 H
HETATM 16 H 0 -2.995 1.931 -3.264 H
HETATM 17 H 0 -4.610 0.068 -3.082 H
END
CONECT 1 2 4 7
CONECT 2 3 1
CONECT 3 2
CONECT 4 5 6 1
CONECT 5 4
CONECT 6 4
CONECT 7 8 9 1
CONECT 8 7 10 11
CONECT 9 7 12 13
CONECT 10 8 14 15
CONECT 11 8
CONECT 12 14 9 16
CONECT 13 9
CONECT 14 10 12 17
CONECT 15 10
CONECT 16 12
CONECT 17 14
|
自动指定的原子类型有误, 需要手动重新指定. 原子类型位于opls_736
附近.
石墨烯
gra.pdb | |
---|---|
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 | TITLE tip
REMARK 1 File created by GaussView 5.0.9
HETATM 1 C 0 0.500 0.500 0.500 C
HETATM 2 C 0 1.710 1.200 0.500 C
HETATM 3 C 0 1.710 2.600 0.500 C
HETATM 4 C 0 0.500 3.300 0.500 C
HETATM 5 C 0 2.920 0.500 0.500 C
HETATM 6 C 0 4.130 1.200 0.500 C
HETATM 7 C 0 4.130 2.600 0.500 C
HETATM 8 C 0 2.920 3.300 0.500 C
HETATM 9 C 0 0.500 4.700 0.500 C
HETATM 10 C 0 1.710 5.400 0.500 C
HETATM 11 C 0 1.710 6.800 0.500 C
HETATM 12 C 0 0.500 7.500 0.500 C
HETATM 13 C 0 2.920 4.700 0.500 C
HETATM 14 C 0 4.130 5.400 0.500 C
HETATM 15 C 0 4.130 6.800 0.500 C
HETATM 16 C 0 2.920 7.500 0.500 C
HETATM 17 C 0 5.350 0.500 0.500 C
HETATM 18 C 0 6.560 1.200 0.500 C
HETATM 19 C 0 6.560 2.600 0.500 C
HETATM 20 C 0 5.350 3.300 0.500 C
HETATM 21 C 0 7.770 0.500 0.500 C
HETATM 22 C 0 8.980 1.200 0.500 C
HETATM 23 C 0 8.980 2.600 0.500 C
HETATM 24 C 0 7.770 3.300 0.500 C
HETATM 25 C 0 5.350 4.700 0.500 C
HETATM 26 C 0 6.560 5.400 0.500 C
HETATM 27 C 0 6.560 6.800 0.500 C
HETATM 28 C 0 5.350 7.500 0.500 C
HETATM 29 C 0 7.770 4.700 0.500 C
HETATM 30 C 0 8.980 5.400 0.500 C
HETATM 31 C 0 8.980 6.800 0.500 C
HETATM 32 C 0 7.770 7.500 0.500 C
END
CONECT 1 2 22 12
CONECT 2 1 3 5
CONECT 3 2 4 8
CONECT 4 3 9 23
CONECT 5 2 6 16
CONECT 6 5 7 17
CONECT 7 6 8 20
CONECT 8 3 7 13
CONECT 9 4 10 30
CONECT 10 9 11 13
CONECT 11 10 12 16
CONECT 12 1 11 31
CONECT 13 8 10 14
CONECT 14 13 15 25
CONECT 15 14 16 28
CONECT 16 5 11 15
CONECT 17 6 18 28
CONECT 18 17 19 21
CONECT 19 18 20 24
CONECT 20 7 19 25
CONECT 21 18 22 32
CONECT 22 1 21 23
CONECT 23 4 22 24
CONECT 24 19 23 29
CONECT 25 14 20 26
CONECT 26 25 27 29
CONECT 27 26 28 32
CONECT 28 15 17 27
CONECT 29 24 26 30
CONECT 30 9 29 31
CONECT 31 12 30 32
CONECT 32 21 27 31
|
周期性的体系, 无法使用自动的方法指定原子类型(否则程序会陷入死循环), 手动指定所有原子类型为opls_147
即可.
待完成
- 反常二面角原子顺序可能与其他程序所给的不一致, 待考.
- OPLS原子类型示例图已完成, 但可能存在错误之处, 待查.
- 扩展到支持其他力场的原子类型, 如amber, charmm等.
- 对分子结构进行分析, 显示原子类型时按匹配度高低排序, 尽可能自动化, 减少选择时的纠结.
- 美化界面
- 分子结构显示可换用Three.js, CH5M3D效果和效率不佳.
致谢
没有下面这些人的热心付出, GMXTOP工具是很难完成的. 如果你使用这个工具, 请感谢他们的努力与付出.
- 张楠@北京 : 58-96
- 蒲中机@大连 : 101-130
- 黄建湘@杭州 : 141-173
- 叶盛@合肥 : 178-192, 197-212
- 梅龙灿@武汉 : 217, 222-241
- 刘恒江@上海 : 247-282
- 郝阳@上海 : 285-318
- 马郑@天津 : 320-348
- 康文渊@成都 : 349-380
- 刘清南@西安 : 398-424
- 杜春保@西安 : 425-459
- 李正@西安 : 460-473, 490-496
- 吴思晋@大连 : 497-537
- 刘凤海@成都 : 569-598
- 刘胜堂@苏州 : 603-645
- 张国成@成都 : 645-681
- 席昆@武汉 : 542-565, 677-713, 913-940
- 吴念@武汉 : 714-748
- 李乐乐@成都 : 749-758, 771-779, 785
- 王新宇@天津 : 941-MW