TRAPPE_small力场包

喵星大佬

==========2023-04-08更新===========
由于有用户反映该力场常与水模型合用以模拟水溶液的情况，因此更新了力场中的设置并添加了水模型的相关支持

对于默认使用SETTLE的硫化氢和环氧乙烷分子，添加了使用SHAKE约束的参数

如果需要模拟这两种分子在水溶液中的情况并对水使用SETTLE约束，可以在mdp文件中设置

1. define                                        = -DNOSETTLES

复制代码

选项激活相关设置以保证可以正确对水使用SETTLE

另一项更新是在二氧化碳中使用了虚拟质量点。
由于二氧化碳为直线分子，如果直接使用SHAKE/RATTLE算法对3个键长进行约束将，180度的键角会导致约束步的Jacobian矩阵奇异而无法施加约束。因此，为了保证直线约束的稳定性，重新调整了约束方式。采用新参数后，可以正确的对二氧化碳分子施加LINCS约束
该参数已通过二氧化碳水溶液NPT模拟中的RDF进行了验证

此外，由于SHAKE约束方法无法在能量极小化过程中使用，因此建议如果需要做能量极小化，可以使用布朗动力学进行近似能量极小化


============分割线==============

TRAPPE力场是明尼苏达大学Siepmann课题组开发的一系列力场，主要用于包括相平衡在内的一系列模拟
其中的Small版本包含了7种简单分子，为硫化氢，氮气，氧气，环氧乙烷，氨气，二氧化碳和氦气
其官网为： http://trappe.oit.umn.edu/#small

本人根据官网和原文中提供的数据(不包含氦气，不管原文还是官网都没见到参数，但是页面上又说有)，将其转换为Gromacs的格式力场包，可进行相关模拟


力场包中包含了所有6种分子的itp文件和pdb文件

6种气体单独存在的情况下已经过测试(测试条件为常温常压，盒子大小约30 nm，内含1000个各种分子)无误，注意其建议范德华截断半径为1.4 nm，且须结合PME计算长程静电作用

由于其均为刚性模型，考虑到约束的运行效率，对硫化氢和环氧乙烷设置了SETTLE约束，因此不能同时加入这两种分子，当然由于本来这两种分子之间就会发生反应，不应该共存，所以实际并没有影响。但是注意，此时不能和任意刚性水模型合用，除非修改水模型的约束，因为Gromacs中的水模型默认用SETTLE的参数，但是一个模拟中只能有一个分子使用SETTLE

对其他的默认为LINCS约束，对氨气和二氧化碳由于涉及非对角元，LINCS算法的近似矩阵求逆不成立，只能使用SHAKE约束迭代求解验证中CO2较同样原子数的H2S慢了相当多(约20-30%)

在forcefield.doc中对各个参数的来源作了说明

欢迎使用

ByTheWay

谢谢老师，帮大忙了

KRAS

谢谢老师

喵星大佬

更新了版本，欢迎大家测试和使用

喵星大佬

之前CO2的键长设置错误，现已修正

xiaoruoqiu

请问这些力场参数能否移植到诸如oplsaa、amber、charmm力场里？谢谢

喵星大佬

应该不可以吧

如果需要其他分子可以先看看trappe的网站有没有

喵星大佬

有同学反应在刚结束的gromacs培训班文件包中trappe_small的力场文件有误，经检查为4月5日的初版，各位同学如需使用请在本帖中下载最新版，在新版中已对错误进行修正并更改约束形式以使动力学稳定性增加，修改内容见本帖一楼

感谢卢老师分享并同时提醒更新一下文件，谢谢

KRAS

喵星大佬 发表于 2023-11-18 16:45
有同学反应在刚结束的gromacs培训班文件包中trappe_small的力场文件有误，经检查为4月5日的初版，各位同学 ...

感谢楼主大大在力场文件中留了邮箱

lmch

请问
1. CO2分子能否使用LINCS约束?
2. 拓扑文件中的虚拟位点参数是如何确定的？
3. 能够推广使用到其他形式类似的co2力场参数吗？

喵星大佬

lmch 发表于 2024-1-10 21:54
请问
1. CO2分子能否使用LINCS约束?
2. 拓扑文件中的虚拟位点参数是如何确定的？

可以使用lincs约束
仔细看文件，实际上这个CO2是由两个质量22的实际原子构成的，距离通过使其转动惯量与原力场定义的CO2相等来设置的，也就是说实际上有运动方程的只是一个双原子分子，而且这两个假原子不直接受到任何力
3个在力场设置的原子的位置是虚拟作用位点，其通过两个假原子的坐标线性组合来确定位置，其受力会通过线性组合给到两个假原子上(计算上就是链式求导)。
这种方式可以扩展到任何刚体线性分子上，因为实际上这个刚体3原子分子只有5个自由度(有3个约束条件和两个简并自由度)，而这和刚体双原子分子是没有区别的(本质上就是3个质心位置坐标和两个欧拉角)。因此这个方法完全可以扩展，保证所有原始原子(相互作用位点)的力矩，分子总质量，分子转动惯量与原始分子一样即可，所以原则上任线性刚体分子都可以转换成无相互作用双原子分子+虚拟相互作用位点的形式

lmch

喵星大佬 发表于 2024-1-11 03:00
可以使用lincs约束
仔细看文件，实际上这个CO2是由两个质量22的实际原子构成的，距离通过使其转动惯量与 ...

感谢您的解答
还有个问题哈，您在文件包中的这些设置是否是正确的？

ffnonbonded.itp中
CC           6       0.0000  0.0000  A   0.280000000  0.224478000 ;C as virtual site in CO2

元素质量也是0吗？

对其他的默认为LINCS约束，对氨气和二氧化碳由于涉及非对角元，LINCS算法的近似矩阵求逆不成立，只能使用SHAKE约束迭代求解验证中CO2较同样原子数的H2S慢了相当多(约20-30%)

这里的意思是指？

喵星大佬

lmch 发表于 2024-1-11 15:45
感谢您的解答
还有个问题哈，您在文件包中的这些设置是否是正确的？

当然是0，真正在原子位置的是虚拟点，只是用来提供真实相互作用的力矩。真正用来描述质心运动和分子转动的是那两个22质量啥相互作用都没有的原子。lincs约束的也是这个玩意的运动

对于CO2那句已经作废了，当时用SHAKE算法约束三点共线，后来发现在求解约束方程的时候矩阵不是绝对对角占优且不保证正定，此时SHAKE方法采用的高斯-赛德尔迭代不成立，产生的结果都是异常的，结构会变得很奇怪，因此改用了后来的方法

lmch

喵星大佬 发表于 2024-1-11 21:59
当然是0，真正在原子位置的是虚拟点，只是用来提供真实相互作用的力矩。真正用来描述质心运动和分子转动 ...

OC           8      16.00    0.0000  A   0.305000000  0.656806000 ;O as virtual site in CO2

那请问为何这里不为0呢？

喵星大佬

lmch 发表于 2024-1-12 17:29
那请问为何这里不为0呢？

可能当时写错了。。。。。我查了下我这里的是对的，我重新传一下

不过无所谓其实，因为在CO2的itp里面重新写了一遍质量，是0，所以在[ atomtypes ]里面写的质量会被覆盖掉