基于amber的MCPB.py构建含金属离子蛋白的力场参数

casea

使用amber、gaussian、MCPB.py对金属蛋白系统参数化(成键模型) 并将其转化为gromacs格式  
前言：此例改编自amber官网教程，如有错误之处，还请批评指教

1. http://ambermd.org/tutorials/advanced/tutorial20/mcpbpy.php

复制代码

       最近在对含有二价铁离子的蛋白体系进行分子动力学模拟，该问题的关键在于铁离子的力场文件如何创建。在很多力场中(如amber forcefield)存在有铁离子的力场参数，但是仅限于血红素铁的使用。在尝试了很多种方法之后，最终选择基于amber和gaussian的MCPB方法。
       在蛋白质（包含 20 个标准氨基酸）和核酸（包含 5 个标准核苷）中，AMBER 力场中有标准库和参数文件。蛋白质中的金属位点变化很大，因此很难为它们参数化标准力场。考虑到这些因素，Merz 研究小组的 Martin Peters 博士在 AmberTools 的 Modeling ToolKit++ (MTK++) 软件包中开发了 Metal Center Parameter Builder (MCPB)它是一种半自动化工作流程，可加快蛋白质系统中金属位点的键合模型参数化过程。MCPB.py是 MCPB 的 Python 版本。它使用优化的工作流程和适当数量的步骤来执行建模。它仅在 AmberTools 15 或更高版本中可用。当前版本支持超过 80 种金属离子和各种 AMBER 力场。它由密歇根州立大学 Merz 研究小组的 Pengfei Li 开发。在此，我们选择MCPB.py进行金属蛋白系统的参数化。
       对于MCPB方法，我们首先推荐获得金属蛋白的解析结构，使用同源建模得出的模型亦可，期间可能需要自己调试金属离子摆放的位置。此教程选取pdb ID: 1mlw的金属蛋白。
1. PDB bank获取1mlw结构
- 首先下载1mlw的pdb格式文件
- 3Dview窗口观看金属与残基的相互作用(需要记住这些残基，后续有用)
2. PDB文件以及mol2文件准备
         所需的PDB文件需要符合PDB 3.0格式并且残基命名格式(e.g. HID, HIP, HIE for HIS)必须与amber一致.同时我们需要为非标准残基建立mol2文件(e.g. 金属离子、水和配体)。
         原始文件即为我们上一步下载的文件: 1mlw.pdb
         这一步十分重要，错误大多出在这。必须小心谨慎的处理pdb文件。
         对于pdb文件，请确保同一残基中的原子不使用相同的原子命名
         MCPB.py支持超过80种离子，对于这些金属离子的残基名以及原子名都必须大写，否则不能被识别。
1）为非标准残基准备pdb和mol2文件
        请注意，为每个 mol2 文件设定正确的总电荷非常重要。因为 MCPB.py会将这些电荷加在一起以确定小模型和大模型的总电荷，并在 RESP 电荷拟合步骤期间对模型应用总电荷的电荷限制。
对于配体
         使用reduce进行配体加氢，或者可以使用Gaussianview手动加氢。
         将1mlw.pdb中的配体HBI提取出来并保存为单一的PDB文件: HBI.pdb。下一步使用reduce进行加氢。

1. reduce HBI.pdb  > HBI_H.pdb

复制代码

如果在第一步查看pdb 3D结构发现离子与配体形成了键，则需要删除配体上相应的氢原子。在此处铁离子与氧原子成键，因此不需要额外的删除氢原子。
然后使用antechamber生成mol2文件，使用AM1-BCC电荷方法生成电荷，并使用GAFF原子类型。
-nc指定配体的电荷，特别重要 。在其他体系，如果需要删除氢原子，就必须使用正确的电荷。

1. antechamber -fi pdb -fo mol2 -i HBI_H.pdb -o HBI_pre.mol2 -c bcc -pf y -nc 0

复制代码

检查HBI_pre.mol2文件后我们注意到HBI_pre.mol2文件中没有“du”原子类型，因此我们可以将其重命名为HBI.mol2文件(提示：非标准残基的mol2文件名应该与其在 PDB 文件中的残基名称相同（字母和大写），后缀为“.mol2”以用于 MCPB.py建模)。
然后使用parmchk2获得配体的力场参数文件

1. parmchk2 -i HBI.mol2 -o HBI.frcmod -f mol2

复制代码

对于金属离子
我们只需要将pdb文件中的FE离子信息单独保存为FE.pdb, 然后使用脚本metalpdb2mol2.py生成mol2文件。使用此脚本的准备工作就是修改第一行的python路径(可使用默认)。该脚本包含在AmbeTools2020的scripts中，之前版本需要点击链接下载。

1. python metalpdb2mol2.py -i FE.pdb -o FE.mol2 -c 3

复制代码

-c指定金属离子所带电荷
对于水
本例中为了简便不保存结晶水
如果想保存结晶水，可以将其复制到一个单独的pdb文件中，然后使用tleap进行加氢操作** 值得注意的是，保存的结晶水mol2文件中的atomtype必须修改为HW和OW(amber格式）

1. tleap -s -f wat_tleap.in > wat_tleap.out
   
2. # wat_tleap.in
   
3. loadoff solvents.lib
   
4. HOH = TP3
   
5. wat = loadpdb WAT.pdb
   
6. savepdb wat WAT_H.pdb
   
7. quit

复制代码

  

使用antechamber生成水分子的mol2文件

1. antechamber -fi pdb -fo mol2 -i WAT_H.pdb -o WAT.mol2 -at amber -c bcc -pf y

复制代码

2) 准备包含标准残基的pdb文件
在这使用WebServer  H++（使用默认设置）将氢原子添加到PDB文件。 H++将在建模过程中删除非标准残基。除了生成含氢pdb文件，还将生成amber格式的拓扑和坐标文件。然后我们可以使用拓扑和坐标文件来生成PDB文件，因为该PDB文件将使用amber命名方式命名残基(例如HID，HIE，HIE，而不是HIS)。在这里，我们执行命令：

1. ambpdb -p 0.15_80_10_pH6.0_1mlw.top -c 0.15_80_10_pH6.0_1mlw.crd > 1mlw_H.pdb

复制代码

H++在加氢过程中忽略金属离子。同时我们需要修改HIP-169为HID-169(删除HE2原子，根据3D结构图删除原子)，修改HIP-174为HID-174(删除HE2原子)，此步骤是为了离子与残基成键，因此使用此方法必须获得金属蛋白的解析结构，同时必须清楚金属离子与周围残基之间的成键关系。修改的pdb文件命名为1mlw_fixed_H.pdb
3) 将pdb文件整合在一个pdb文件中
注意整合顺序，对于之后的拓扑以及模拟很重要。
在这，顺序为：标准残基、金属离子、配体

1. cat 1mlw_fixed_H.pdb FE.pdb HBI_fixed_H.pdb > 1mlw_H.pdb

复制代码

使用antechamber对pdb文件重新标号

1. pdb4amber -i 1mlw_H.pdb -o 1mlw_fixed_H.pdb

复制代码

3. 生成pdb、Gaussian和建模文件

1. # input file
   
2. # 1mlw.in
   
3. original_pdb 1mlw_fixed_H.pdb
   
4. group_name 1mlw
   
5. cut_off 2.8
   
6. ion_ids 4645 #  根据自己pdb文件中的原子编号修改
   
7. ion_mol2files FE.mol2
   
8. naa_mol2files HBI.mol2
   
9. frcmod_files HBI.frcmod
   
10. large_opt 1
    
11. software_version g16

复制代码

我们需要手动创建输入文件。以上是输入文件 1mlw.in。用户可以查看 AMBER 2016 参考手册中的第 287-291 页，了解MCPB.py 使用的输入参数。使用 Gaussian09 或 GAMESS-US 来执行 QM 计算。使用 Gaussian09 时，请将“software_version g09”作为一行添加到1mlw.in文件中，本例我们使用gaussian16进行计算。要使用 GAMESS-US 程序，请将“software_version gms”作为一行添加到1mlw.in文件中。这里我们设置 large_opt 等于 1（默认为 0，表示不优化大模型）来优化大模型的氢位置以进行 RESP 计算，因为配体的氢位置可能不稳定。这里我们使用 ff14SB 力场（默认）来执行建模。我们执行以下命令:

1. MCPB.py -i 1mlw.in -s 1

复制代码

此步会生成小模型、标准模型和大模型的PDB和fingerprint文件。还会生成小模型和大模型的高斯输入文件。大模型的高斯输入文件将首先对氢原子进行优化，以纠正任何放置不当的氢原子。标准模型fingerprint文件（1mlw_standard.fingerprint）包含标准模型中原子的原子类型信息（第3列是原始原子类型，最后一列是最终分配的原子类型）。金属离子的原子类型（名称以 M 开头）和连接原子（名称以 Y/Z 开头）会自动重命名，以将它们与 AMBER 力场中的其他原子类型区分开来。您可以将 step_number 设为 1n 来强制不重命名任何原子类型并生成fingerprint文件，然后手动修改（在此之前，请检查 AMBER force field parm*.dat 文件以确保您的金属离子和连接原子存在amber已有的原子类型）。fingerprint文件在文件末尾还包含金属离子与周围原子之间的链接信息，开头字母为“LINK”。用户还可以根据自己的需要手动修改它们（例如，如果 QM 几何优化后有任何连接更改）

4. 执行Gaussian计算
         根据第三步获得的gaussian输入文件执行gaussian计算，对输入文件的cpu核数、内存用量、方法、基组进行修改以获得高效的计算。此例使用B3LYP/6-31G*方法，用户也可自行选定计算方法。但小模型的几何优化和力常数计算应当使用同一水平的方法，以确保找到局部最小值。为确保参数化合理，建议用户在优化后检查小模型的结构。这是因为几何优化后金属的一些配位键可能会断裂，在这种情况下，最终的参数化结果可能有问题，可以尝试不同水平的理论或其他策略（例如执行多步优化）以获得有意义的结果。在大模型的高斯计算中，我们只对氢原子进行优化，以节省计算任务。可以通过修改 MCPB.py 输入文件中的“large_opt”变量来更改此设置。此步十分耗时，建议使用并行计算。
           

1. g16 < 1mlw_small_opt.com > 1mlw_small_opt.log
   
2. g16 < 1mlw_small_fc.com > 1mlw_small_fc.log
   
3. formchk 1mlw_small_opt.chk > 1mlw_small_opt.fchk
   
4. # 对大模型执行RESP电荷计算
   
5. g16 < 1mlw_large_mk.com > 1mlw_large_mk.log

复制代码

5. 执行最后的建模
使用 半经验方法来生成力场参数。也可使用其他选项：Z 矩阵（设置 step_number 2z）和经验（设置 step_number 2e）方法。经验方法不需要任何高斯计算来获得力常数（但仍然需要高斯计算来获得 RESP 电荷），但它在当前版本中仅支持锌离子建模。

1. MCPB.py -i 1mlw.in -s 2
   
2. # 此步骤我们将获得1mlw_mcpbpy.frcmod文件

复制代码

使用ChgModB方法执行RESP电荷计算，并为金属离子残基生成mol2文件

1. MCPB.py -i 1mlw.in -s 3
   
2. # 此步骤获得金属位点的五个mol2文件：HD1.mol2、HD2.mol2、HE1.mol2、FE1.mol2、GU1.mol2

复制代码

生成tleap输入文件

1. MCPB.py -i 1mlw.in -s 4

复制代码

我们将获得一个新的pdb文件1mlw_mcpbpy.pdb以及一个leap输入文件1mlw_tleap.in。可以修改in文件中使用的力场等参数。

然后生成拓扑和坐标文件

1. # 1mlw_tleap.in
   
2. source oldff/leaprc.ff14SB
   
3. source leaprc.gaff
   
4. addAtomTypes {
   
5.         { "M1"  "Fe" "sp3" }
   
6.         { "Y1"  "N" "sp3" }
   
7.         { "Y2"  "N" "sp3" }
   
8.         { "Y3"  "O" "sp3" }
   
9. }
   
10. HD1 = loadmol2 HD1.mol2
    
11. HD2 = loadmol2 HD2.mol2
    
12. GU1 = loadmol2 GU1.mol2
    
13. FE1 = loadmol2 FE1.mol2
    
14. loadoff HBI.lib
    
15. loadamberparams HBI.frcmod
    
16. loadamberparams frcmod.ions1lm_126_tip3p
    
17. loadamberparams 1mlw_mcpbpy.frcmod
    
18. mol = loadpdb 1mlw_mcpbpy.pdb
    
19. bond mol.169.NE2 mol.291.FE
    
20. bond mol.174.NE2 mol.291.FE
    
21. bond mol.214.OE2 mol.291.FE
    
22. bond mol.168.C mol.169.N
    
23. bond mol.169.C mol.170.N
    
24. bond mol.173.C mol.174.N
    
25. bond mol.174.C mol.175.N
    
26. bond mol.213.C mol.214.N
    
27. bond mol.214.C mol.215.N
    
28. savepdb mol 1mlw_dry.pdb
    
29. saveamberparm mol 1mlw_dry.prmtop 1mlw_dry.inpcrd
    
30. solvatebox mol TIP3PBOX 10.0
    
31. addions mol Na+ 0
    
32. addions mol Cl- 0
    
33. savepdb mol 1mlw_solv.pdb
    
34. saveamberparm mol 1mlw_solv.prmtop 1mlw_solv.inpcrd
    
35. quit

复制代码

1. tleap -s -f 1mlw_tleap.in > 1mlw_tleap.out

复制代码

我们将获得坐标和拓扑文件(1mlwsolv.prmtop and 1mlw_solv.inpcrd)

注：需要注意的是，在本例中会报错FATAL:  Atom .R<HBI 292>.A<N1 1> does not have a type.

解决方式建立HBI的lib文件，在leap.in输入文件中加载此lib文件即可。

1. # 使用如下命令生成lib文件
   
2. tleap
   
3. HBI =loadmol2 HBI.mol2
   
4. saveoff HBI HBI.lib

复制代码

6. 检查模型
在进行最小化和MD之前，需要使用VMD以及cpptraj进行快速检查以确保模拟的顺利/

1）使用VMD检查拓扑文件中是否存在金属离子的成键信息

1. vmd -parm7 1mlwsolv.prmtop -rst7 1mlwsolv.inpcrd

复制代码

2) 使用cpptraj检查原子编号

1. cpptraj -p 1mlw_solv.prmtop

复制代码

3) 使用parmed检查金属离子参数

1. # mcpbpy_parmed.in
   
2. printBonds :FE1
   
3. printAngles :FE1
   
4. printDihedrals :FE1
   
5. printDetails :FE1

复制代码

1. parmed -i mcpbpy_parmed.in -p 1mlw_solv.prmtop

复制代码

金属离子相关参数遵循以下标准：

A) 金属离子与其连接原子之间的bond force constants小于200 kcal/(mol*Angstrom^2)，且平衡键距小于2.8 A；

B) 与金属离子相关的angle force constants通常小于100 kcal/(mol*Rad^2)，而平衡角值大于100度；

C) 对于涉及金属的二面体，所有或大部分二面体势垒为零；

D）金属离子的RESP电荷小于其氧化态，通常甚至小于+1；

E) 一种金属离子的 LJ 半径通常大于 1.0 埃

7. 将amber参数文件转换为gromacs格式

1. import parmed as pmd
   
2. #convert prmtop and inpcrd into top and gro
   
3. amber = pmd.load_file('1mlw_solv.prmtop','1mlw_solv.inpcrd')
   
4. amber.save('1mlw.top')
   
5. amber.save('1mlw.gro')

复制代码

Joysharing

Joysharing 发表于 2024-11-18 02:12
想请问一下各位老师，最后一步使用parmed将amber文件转换为gromacs文件出现这种情况怎么办呢，下面是报错文 ...

继续求问，感谢各位老师

Joysharing

想请问一下各位老师，最后一步使用parmed将amber文件转换为gromacs文件出现这种情况怎么办呢，下面是报错文件，326,SC1是一个与FE配位的配体，5004是FE的位置，所以这里Could not find <Atom C1 [5004]; In SC1 326>我不太明白
Traceback (most recent call last):
  File "E:\bioalign\amber-gro.py", line 8, in <module>
    amber.save('CAS2-PC-1112_solv.gro')
  File "C:\Users\97490\anaconda3\envs\align\lib\site-packages\parmed\structure.py", line 1540, in save
    gromacs.GromacsGroFile.write(self, fname, **kwargs)
  File "C:\Users\97490\anaconda3\envs\align\lib\site-packages\parmed\gromacs\gromacsgro.py", line 289, in write
    raise RuntimeError(f"Could not find {atom}")
RuntimeError: Could not find <Atom C1 [5004]; In SC1 326>

light_123

M107 发表于 2024-9-26 16:00
老师您好，想问一下在使用pdb4amber重新标号时出现以下报错：
========================================= ...

您好，我也遇到了这样的问题，想问下您解决了吗？

M107

老师您好，想问一下在使用pdb4amber重新标号时出现以下报错：
==================================================
Summary of pdb4amber for: 1mlw_H.pdb
===================================================
Traceback (most recent call last):
  File "/opt/anaconda3/envs/AmberTools23/bin/pdb4amber", line 33, in <module>
    sys.exit(load_entry_point('pdb4amber==22.0', 'console_scripts', 'pdb4amber')())
             ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/pdb4amber/pdb4amber.py", line 819, in main
    run(
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/pdb4amber/pdb4amber.py", line 516, in run
    parm = parmed.read_PDB(pdbin)
           ^^^^^^^^^^^^^^^^^^^^^^
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/parmed/formats/pdb.py", line 447, in parse
    inst.struct.assign_bonds()
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/parmed/structure.py", line 886, in assign_bonds
    templ = _res_in_templlib(res, all_residues)
            ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/parmed/structure.py", line 4302, in _res_in_templlib
    if residue.DNAResidue.has(res.name):
       ^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/parmed/residue.py", line 51, in has
    cls.get(thing)
  File "/opt/anaconda3/envs/AmberTools23/lib/python3.12/site-packages/parmed/residue.py", line 216, in get
    if key[-1] in '35':
       ~~~^^^^
IndexError: string index out of range

想问一下老师这个报错是怎么回事呢，需要怎样解决？

不甚感激！！

YuhangYao

1009836241 发表于 2024-1-1 20:43
请问一下 我的蛋白结晶里面有两个镁和一个钠 在用MCPB.py时 执行第一步出现报错
***Creating the small mo ...

请问你有解决多个金属的问题吗？

t11b41

感谢您分享的教程，非常有帮助。
现在已经按照教程成功生成了top和crd文件，自动输出的pdb看着是正常的，使用ambpdb重新生成的pdb文件结构也是正常的。但是使用vmd读取top和crd文件检查成键参数时结构是乱的，想请问问题出在哪里了？
top和crd文件随附件上传。

1009836241

请问一下 我的蛋白结晶里面有两个镁和一个钠 在用MCPB.py时 执行第一步出现报错
***Creating the small model...
It contains the residue 41-ASP as sidechain coordinated.
It contains the residue 147-GLU as sidechain coordinated.
Creating the residue 193-SER into ACE...
Creating the residue 194-ALA into NME...
It contains the residue 202-ASP as sidechain coordinated.
It contains the residue 250-NA as normal.
It contains the residue 251-MG1 as normal.
It contains the residue 252-MG2 as normal.
It contains the residue 253-ATP as normal.
It contains the residue 254-nam as normal.
It contains the residue 255-W1 as normal.
It contains the residue 256-W2 as normal.
It contains the residue 257-W3 as normal.
It contains the residue 258-W4 as normal.
Totally there are 138 atoms in the small model.
Traceback (most recent call last):
  File "/root/miniconda3/envs/AmberTools23/bin/MCPB.py", line 651, in <module>
    gene_model_files(orpdbf, ionids, addres, addbpairs, gname, ff_choice,
  File "/root/miniconda3/envs/AmberTools23/lib/python3.12/site-packages/pymsmt/mcpb/gene_model_files.py", line 1902, in gene_model_files
    build_small_model(mol, reslist, ionids, cutoff, smresids, smresace, smresnme,
  File "/root/miniconda3/envs/AmberTools23/lib/python3.12/site-packages/pymsmt/mcpb/gene_model_files.py", line 1380, in build_small_model
    AtNum = Atnum[gatm.element]
            ~~~~~^^^^^^^^^^^^^^
KeyError: 'M'
但是我将其中一个镁删去MCPB.py就可以执行了 老师们知道怎么回事吗？

liu19950830

nianbin 发表于 2022-3-22 16:26
现在有了卢老师的sobtop了，这个基本可以被替代了

就是我也是金属蛋白与配体模拟。然后是3个金属镁离子与5个氨基酸的原子和配体的原子形成配位，请问用sobtop怎么构建拓扑进行GROMACS模拟。万分感谢，期待您的回复

dayu8278

casea 发表于 2023-4-28 07:33
是的，你的水分子的分子名以及原子类型都一样，因此只需制作单个水分子的top就行了。之后无论体系内有多 ...

谢谢老师

casea

dayu8278 发表于 2023-4-27 18:33
谢谢老师，保留任意一个都可以吗？

是的，你的水分子的分子名以及原子类型都一样，因此只需制作单个水分子的top就行了。之后无论体系内有多少水，都会使用该top。可以看一下top的使用原理

dayu8278

casea 发表于 2023-4-26 11:40
amtechamber只支持单个分子，你想在模型中保留多个水，只需要创建一个水分子的top就行了。

谢谢老师，保留任意一个都可以吗？

casea

dayu8278 发表于 2023-4-25 19:28
您好！

我在运行antechamber -fi pdb -fo mol2 -i WAT_H.pdb -o WAT.mol2 -at amber -c bcc -pf y这一步 ...

amtechamber只支持单个分子，你想在模型中保留多个水，只需要创建一个水分子的top就行了。

dayu8278

您好！

我在运行antechamber -fi pdb -fo mol2 -i WAT_H.pdb -o WAT.mol2 -at amber -c bcc -pf y这一步时，amber报错了。

报错信息为：

Welcome to antechamber 22.0: molecular input file processor.

Info: acdoctor mode is on: check and diagnose problems in the input file.
Info: The atom type is set to amber; the options available to the -at flag are
      gaff, gaff2, amber, bcc, and sybyl.

-- Check Format for pdb File --
   Status: pass
Warning: Detected more than 10 Residue sequence numbers;
         this may be a large multiple residue PDB file;
         large multiple residue PDB files are not supported.
         This warning usually indicates a conceptual misunderstanding.
         We recommend reviewing the Information flow in Amber documentation
         and the antechamber tutorials.
         Continuing, but problems may be encountered.
-- Check Unusual Elements --
   Status: pass
-- Check Open Valences --
   Status: pass
-- Check Geometry --
      for those bonded   
      for those not bonded   
   Status: pass
-- Check Weird Bonds --
   Status: pass
-- Check Number of Units --
/opt/mamba/envs/AmberTools22/bin/wrapped_progs/antechamber: Fatal Error!
This molecule may have more than one unit.
       antechamber can only handle one unit.  If the input is a single unit
       then the connectivity is wrong and the geometry may be bad.
       Please convert your molecule to a mol2 file via:
       antechamber -j 5 -at sybyl -dr no
       And then check your molecule with a visualization program;
       manually add missing bonds or delete unwanted bonds as appropriate.

我使用wat_tleap.in脚本对WAT.pdb加氢后的水文件WAT_H.pdb中水分子格式大致为：
ATOM      1  O   WAT     1      44.426  25.269  23.516  1.00  0.00
ATOM      2  H1  WAT     1      45.383  25.269  23.516  1.00  0.00
ATOM      3  H2  WAT     1      44.186  26.196  23.516  1.00  0.00
TER   
ATOM      4  O   WAT     2      34.238  31.812  45.122  1.00  0.00
ATOM      5  H1  WAT     2      35.195  31.812  45.122  1.00  0.00
ATOM      6  H2  WAT     2      33.998  32.739  45.122  1.00  0.00
TER   
ATOM      7  O   WAT     3      35.805  24.387  48.950  1.00  0.00
ATOM      8  H1  WAT     3      36.762  24.387  48.950  1.00  0.00
ATOM      9  H2  WAT     3      35.565  25.314  48.950  1.00  0.00
TER   
ATOM     10  O   WAT     4      43.121  54.497  38.143  1.00  0.00
ATOM     11  H1  WAT     4      44.078  54.497  38.143  1.00  0.00
ATOM     12  H2  WAT     4      42.881  55.424  38.143  1.00  0.00
TER   
。。。。。

请问怎么解决，谢谢

casea

sunmeng 发表于 2023-3-15 10:32
老师，请问 如果要指定金属离子和哪些原子结合是要怎么做呢？（还是不可指定呢）

见第5步是如何生成topol文件的
使用bond指定成键

1. bond mol.169.NE2 mol.291.FE
   
2. bond mol.174.NE2 mol.291.FE
   
3. bond mol.214.OE2 mol.291.FE
   
4. bond mol.168.C mol.169.N
   
5. bond mol.169.C mol.170.N
   
6. bond mol.173.C mol.174.N
   
7. bond mol.174.C mol.175.N
   
8. bond mol.213.C mol.214.N
   
9. bond mol.214.C mol.215.N
   

复制代码