求助金颗粒与多肽动力学模拟添加金属原子力场问题

不去重蹈 · 发表于 Post on 2021-3-5 20:43:05

老师您好，请教一个关于力场方面的问题。

背景：实验得纳米金颗粒能和多肽吸附（例如，CCY等肽），我们想在gromacs 2019.3中模拟计算纳米金颗粒和多肽CCY的结合能，比较不同的肽与金颗粒结合能力的强弱。
金原子的力场用 “在Gromacs中模拟金纳米线拉伸过程 http://sobereva.com/153 中 ” 的参数。在gromacs\top\gromos53a6.ff 路径下，将参数添加进ions.itp 文件中 f'f'nonbonded.itp

先是生成了CCY多肽的gro和itp文件。然后生成金颗粒的力场时遇到了报错，
gmx pdb2gmx -f .\Au13fixed.pdb -o .\Au13fixed.gro
报错
Fatal error:
Residue 'AU' not found in residue topology database

问题：适合多肽的力场很好选择，但是在构建Au13的力场时，请问应该怎么操作才能成功得到金的 top 和 itp 文件

附件：Au13的pdb文件。

angervlf · 发表于 Post on 2023-7-11 04:54:28

本帖最后由 angervlf 于 2023-7-11 04:57 编辑

牧生发表于 2023-6-27 09:37
因为文献给出的参数不是标准的参数，所以要转化

是的是的，感谢，sob大佬的那篇博文后面没有看，现在理解了

牧生 · 发表于 Post on 2023-6-27 09:37:43

angervlf 发表于 2023-6-27 06:11
我不太理解，gromacs软件里自带的amber力场不也是已经经过转化的嘛，用的也是C6和C12吧，为什么还要进行 ...

因为文献给出的参数不是标准的参数，所以要转化

angervlf · 发表于 Post on 2023-6-27 06:11:29

lyj714 发表于 2021-3-12 12:22
我怀疑你金的参数不对，你这个top中也看不出你具体非键参数写的什么(ffnonbonded.itp)。如果你是照搬sob的 ...

我不太理解，gromacs软件里自带的amber力场不也是已经经过转化的嘛，用的也是C6和C12吧，为什么还要进行转化？

牧生 · 发表于 Post on 2021-6-20 17:46:03

本帖最后由牧生于 2021-6-20 17:48 编辑

lyj714 发表于 2021-6-20 14:35
就是这个力场的原因啊，非要让我点出问题，我引导你看警告不是没有原因的，如果你做铁，还自己添加了铁的 ...

给大佬递冰可乐。

问题解决了。我输入参数是没有问题的，gromos54a7_atb.ff里面，ffnonbonded.itp只有我输入的参数，先前有的铁参数，我已经删掉的。

真正的原因在NOTE里面提到了，
NOTE 1 [file em.mdp]:
  You have set rlist larger than the interaction cut-off, but you also have
  verlet-buffer-tolerance > 0. Will set rlist using verlet-buffer-tolerance.

于是，将em.mdp里面，将金板适用的这两个0.9，改成1，结果就很正常了。
rcoulomb                   = 0.9 改成1
rvdw                         = 0.9  改成1

lyj714 · 发表于 Post on 2021-6-20 14:35:15

本帖最后由 lyj714 于 2021-6-20 14:39 编辑

牧生发表于 2021-6-20 13:32
对比了一下铝板和铁板与水溶液接触的过程，唯一的区别在于，使用铁板的做em的时候，多了一个note，

NO ...

就是这个力场的原因啊，非要让我点出问题，我引导你看警告不是没有原因的，如果你做铁，还自己添加了铁的原子信息到ffnbonded.itp中中，肯定会有一个atomtype defined previously警告，这就已经体现了错误地方。明显这个力场中已经添加了有一个铁原子类型了，就在ffnobonded.itp中，自己睁大眼睛去看，注释掉所有关于原先这个铁的所有行，包括预先定义的非键参数那些行。
所以我一直都不提倡那些一遇到警告都不仔细看完，就添加什么-maxwarn去忽略的做法，简直害人不浅。有些警告是必须注意的

牧生 · 发表于 Post on 2021-6-20 13:32:43

lyj714 发表于 2021-6-20 11:43
我说了让你看警告信息的，你怎么就不看。这么对比没意义

对比了一下铝板和铁板与水溶液接触的过程，唯一的区别在于，使用铁板的做em的时候，多了一个note，

NOTE 1 [file em.mdp]:
You have set rlist larger than the interaction cut-off, but you also have
verlet-buffer-tolerance > 0. Will set rlist using verlet-buffer-tolerance.

其余的警告都是因为使用gromos54a7_atb.ff力场而出现的警告，这个忽略即可。

那么，可能问题出在这个note上面，这个该怎么解决呢？verlet-buffer-tolerance和rlist是自动设定的吗？

lyj714 · 发表于 Post on 2021-6-20 11:43:47

牧生发表于 2021-6-20 11:30
相同的方法，相同的命令，做金板，铝板和酸性水溶液就完全没有一点问题。唯独做铁板，水分子那些就进入铁 ...

我说了让你看警告信息的，你怎么就不看。这么对比没意义

牧生 · 发表于 Post on 2021-6-20 11:30:12

lyj714 发表于 2021-6-20 11:24
应该是系统某些分子间非键参数有问题，不要忽略警告，看具体的关键性警告

相同的方法，相同的命令，做金板，铝板和酸性水溶液就完全没有一点问题。唯独做铁板，水分子那些就进入铁内部。

lyj714 · 发表于 Post on 2021-6-20 11:24:41

牧生发表于 2021-6-19 21:57
大佬，请教一下，合理的参数可以使得金属只靠LJ势就维持住簇的状态，且在MD的过程中不跑散。

我目前做 ...

应该是系统某些分子间非键参数有问题，不要忽略警告，看具体的关键性警告

牧生 · 发表于 Post on 2021-6-19 21:57:57

sobereva 发表于 2021-3-13 06:03
不需要任何约束
合理的参数下自然而然就能靠LJ势维持住簇的状态

大佬，请教一下，合理的参数可以使得金属只靠LJ势就维持住簇的状态，且在MD的过程中不跑散。

我目前做了个铁，CTAC，水合氢离子的水溶液接触，尽管铁没有跑散，但是其余的分子进入了铁晶胞内，这个能怎么解决啊。。

具体如下：

以ATB得到CTAB阳离子，水合氢离子的itp文件，

atomtypes.atp里面已经有铁的原子量 FE 55.845
γ铁的参数（该参数来自10.1038/s41524-020-00478-1）

按照培训班的教程，计算得到c6 c12的参数，并写入添加到ffnonbonded.itp里面
FE 26 0 0 A 0.01514608 2.2873928e-06

延展铁板
gmxeditconf -f FE.pdb -o FE_plate_box.gro -box 3.656 3.656 6

往真空部分加入2个CTAC和10个H3O+

gmx insert-molecules -f FE_plate_box.gro -ci CTAB.pdb -oout_box.gro -nmol 2

gmx insert-molecules -f out_box.gro -ci H3O.pdb -box 4 46 -o out_box2.gro -nmol 10

加水填满盒子
gmxsolvate -cp out_box2.gro -p topol.top -o FE_wat.gro
得到tpr
gmx grompp -f em.mdp -c FE_wat.gro -p topol.top -o CTABa.tpr -maxwarn 80
加入反离子，
gmx genion -s CTABa.tpr -p topol.top -o CTABb4em.gro -nname CL -nn 12
铁晶在盒子中间，隐藏水以后得到的图像如下，

能量最小化

gmxgrompp -f em.mdp -c CTABb4em.gro -p topol.top -o CTABb.tpr -maxwarn 80

gmxmdrun -deffnm CTABb -v -ntomp10 -ntmpi 1 -pinoffset 10 -gpu_id 0

铁晶体不再保持在盒子中央了，

隐藏水分子以后的图像如下

铁晶胞不在盒子中央，且非水的分子进入了铁晶体的内部。请教一下怎么解决啊。

使用的文件如下：

npt.mdp (842 Bytes, 下载次数 Times of downloads: 3)

H3O.pdb (660 Bytes, 下载次数 Times of downloads: 3)

FE.pdb (160.54 KB, 下载次数 Times of downloads: 4)

H3O.itp (5.09 KB, 下载次数 Times of downloads: 6)

em.mdp (391 Bytes, 下载次数 Times of downloads: 3)

CTAB.pdb (6.37 KB, 下载次数 Times of downloads: 3)

CTAB.itp (18.44 KB, 下载次数 Times of downloads: 5)

topol.top (502 Bytes, 下载次数 Times of downloads: 9)

牧生 · 发表于 Post on 2021-6-17 11:13:46

lyj714 发表于 2021-6-17 10:55
单位，gromacs用的是nm

好的，谢谢。这下子清晰明白了

lyj714 · 发表于 Post on 2021-6-17 10:55:42

牧生发表于 2021-6-17 09:57
请问一下大佬，我用excel去推算你给这个txt的参数，从文本中抄过来的前三个参数，以Ag为例：

>> DEBUG ...

单位，gromacs用的是nm

牧生 · 发表于 Post on 2021-6-17 09:57:07

本帖最后由牧生于 2021-6-17 10:00 编辑

lyj714 发表于 2021-3-12 12:22
我怀疑你金的参数不对，你这个top中也看不出你具体非键参数写的什么(ffnonbonded.itp)。如果你是照搬sob的 ...

请问一下大佬，我用excel去推算你给这个txt的参数，从文本中抄过来的前三个参数，以Ag为例：

>> DEBUG:
12-6 LJ parameters from [JPCC 2008,112,17281–17290]:
      metal  r0  epsilon  A    B
      Ag 2.955 4.56 2021000  6072,
      Al 2.925 4.02 1577000  5035,
      Au 2.951 5.29 2307000  6987,

For amber/charmm/cvff/oplsaa forcefield of gromacs:
>> DEBUG:
This article 12-6 LJ form:
E = eps*[(r0/r)^12 - 2*(r0/r)^6]
Convert it to below form for gromacs:
E = 4*epsilon*[(sigma/r)^12 - (sigma/r)^6]

That is, sigma = r0 * 2^(-1/6)

#metal    sigma(nm) epsilon(kJ/mol)
Ag    0.263260571 19.079040000
Al    0.260587875 16.819680000
Au    0.262904212 22.133360000

我根据公式计算，大了十倍。

sobereva · 发表于 Post on 2021-3-13 06:03:16

不去重蹈发表于 2021-3-12 14:46
谢谢老师！因为换用了力场，所以没有用sob老师的参数。我用的是朋友给的参数，查看ffnonbonded.itp，AU ...

不需要任何约束
合理的参数下自然而然就能靠LJ势维持住簇的状态

		自动登录 Automatic login	找回密码 Forget password
密码 Password			注册 Register

[GROMACS] 求助金颗粒与多肽动力学模拟添加金属原子力场问题