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

不去重蹈 · 发表于 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文件。

sobereva · 发表于 Post on 2021-3-5 22:11:06

看我博文里的例子top里是什么样的。当前体系的金团簇部分明显不是直接靠pdb2gmx来产生拓扑信息的。
先拿pdb2gmx产生多肽的拓扑文件，然后自行把金的信息加入其中。

不去重蹈 · 发表于 Post on 2021-3-11 15:50:03

sobereva 发表于 2021-3-5 22:11
看我博文里的例子top里是什么样的。当前体系的金团簇部分明显不是直接靠pdb2gmx来产生拓扑信息的。
先拿pd ...

谢谢sobereva老师的指导，力场问题已经解决，由于gromos53a6力场不能识别部分氢原子（如HA等），改用amber99sb-lidn力场，做动力学后有两个问题。
1: 分子在穿过边界时，键被拉长了，请问周边应该怎么设置？
2: 金颗粒跑散了。请问金颗粒该怎么限制或固定，金颗粒才不散开？

附图，
第一帧，第n帧，最后一帧。动力学使用mdp文件

sobereva · 发表于 Post on 2021-3-12 02:17:45

不去重蹈发表于 2021-3-11 15:50
谢谢sobereva老师的指导，力场问题已经解决，由于gromos53a6力场不能识别部分氢原子（如HA等），改用ambe ...

1 用dynamic bonds方式显示，或者让每帧重新判断连接关系，看下文
VMD初始化文件(vmd.rc)我的推荐设置
http://sobereva.com/545（http://bbs.keinsci.com/thread-16834-1-1.html）

2 并不需要任何固定。跑散了肯定是拓扑文件有问题

不去重蹈 · 发表于 Post on 2021-3-12 11:16:21

sobereva 发表于 2021-3-12 02:17
1 用dynamic bonds方式显示，或者让每帧重新判断连接关系，看下文
VMD初始化文件(vmd.rc)我的推荐设置
...

谢谢sobereva老师，
请问拓扑文件改怎么修改？

lyj714 · 发表于 Post on 2021-3-12 12:22:55

本帖最后由 lyj714 于 2021-3-12 12:27 编辑

我怀疑你金的参数不对，你这个top中也看不出你具体非键参数写的什么(ffnonbonded.itp)。如果你是照搬sob的那篇博文的参数那么你就错的离谱，因为sob写的那个参数只能用于gromos力场(用的是C6和C12)，你这里用的amber系力场，理应看那篇文献中的参数，进行合适的转换得到ε和σ
给你个本人根据那篇原文献转的各种参数，供参考:
https://github.com/liuyujie714/temp/blob/main/param.txt

不去重蹈 · 发表于 Post on 2021-3-12 14:46:09

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

谢谢老师！因为换用了力场，所以没有用sob老师的参数。我用的是朋友给的参数，查看ffnonbonded.itp，AU 79 11967 0.0000 A 2.93400e-01 1.63300e-01，应该不对。改用您给的再试一试。
请问金颗粒（Au13.pdb）文件中需要给Au原子之间添加键或者做什么约束吗？建模时13个金原子之间没有任何关系，是不是动力学时13个金原不是一个整体，是13个独立的单体？

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

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

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

牧生 · 发表于 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

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

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 11:13:46

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

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

牧生 · 发表于 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)

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

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

我目前做 ...

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

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

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

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

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

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

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

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

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