求助：蛋白配体复合物动力学能量最小化过程不收敛且配体结构发生改变（部分解决）

chenwei

各位老师好，我是利用gromacs做蛋白质-钙离子-配体复合物动力学模拟，复合物的结构是通过分子对接得到的，在Amber F14SB力场和TIP3P水模型下生成蛋白-离子的拓扑文件，
利用sobtop生成GAFF力场下配体的拓扑文件，并在配体.itp文件中补充了用gaussian和Multiwfn计算的resp电荷。在能量最小化过程中，出现以下信息：Energy minimization has stopped, but the forces have not converged to the requested precision Fmax < 100 (which may not be possible for your system).
It stopped because the algorithm tried to make a new step whose size was too small, or there was no change in the energy since last step. Either way, we
regard the minimization as converged to within the available machine precision, given your starting configuration and EM parameters.
Double precision normally gives you higher accuracy, but this is often not needed for preparing to run molecular dynamics.
Polak-Ribiere Conjugate Gradients converged to machine precision in 51 steps,
but did not reach the requested Fmax < 100.
Potential Energy  = -1.2216996e+06
Maximum force     =  1.1368452e+05 on atom 6139
Norm of force     =  5.7332058e+02
我用pymol打开gro文件后发现配体的结构发生了变化。因为我已经做了50ns的单独蛋白、蛋白-离子、蛋白-离子-配体（未给配体添加resp电荷）的动力学，并未出现报错。
所以我觉得应该是配体文件有不合理的地方，但我个人水平有限，不懂是那部分内容不合理，请教各位老师，谢谢。



（附图为配体能量最小化em前/后的结构，标注为结构变化的地方，相关文件已上传附件，由于蛋白的top文件太大，未上传，只上传了跟配体相关的，配体命名为pa）

dzdhp

    遇到这种情况可以先用steep，把emtol 调大点，比如150，然后再cg，cg不成继续steep，总之就是多做最小化，慢慢减少emtol 的数值，交替steep和cg。等差不多了，直接上npt，把步长调小，比如0.001等平衡好了再回归正常做md。
   或者安装双精度版本做最小化，sob老师有帖子教怎么装。

chenwei

谢谢指导，我是试过先steep，emtol=1000，然后cg，再把emtol调成100，试了几次还是不行，我再按照您说的方法试试。

暮晨

楼主大大，这个问题请问一下您解决了吗？

chenwei

暮晨 发表于 2024-10-24 14:44
楼主大大，这个问题请问一下您解决了吗？

目前还没有，找别人做了，等待结果出来，我再分析是哪一步错了，再分享经验

littlefive

后面MD能正常跑，应该没有什么问题。以前社长回答过这类问题http://bbs.keinsci.com/forum.php ... F%D7%EE%D0%A1%BB%AF

chenwei

littlefive 发表于 2024-10-28 15:08
后面MD能正常跑，应该没有什么问题。以前社长回答过这类问题http://bbs.keinsci.com/forum.php?mod=viewthr ...

我试过接着跑MD，刚开始nvt平衡，就出现大量的错误，过程中生成的pdb文件打开后发现配体的结构已经跑散了，也按照社长的检查方法尝试找原因，结果都是跑散，我觉得问题就是出在小分子上，加了resp电荷就出错，可能是我处理的不对，不是很懂这方面

暮晨

chenwei 发表于 2024-10-28 10:58
目前还没有，找别人做了，等待结果出来，我再分析是哪一步错了，再分享经验

好的楼主大大，谢谢你

student0618

配体的polar hydrogen是怎么决定的？也检查一下配体有没有跟蛋白的原子重叠。

chenwei

student0618 发表于 2024-10-28 22:43
配体的polar hydrogen是怎么决定的？也检查一下配体有没有跟蛋白的原子重叠。

好的，谢谢

chenwei

暮晨 发表于 2024-10-28 17:11
好的楼主大大，谢谢你

你好，我目前关于本帖子的问题已经解决。主要问题出现在小分子的拓扑文件上。我先用sobtop生成GAFF力场下的拓扑文件，在itp文件中没有improper项的内容（显示No improper needs to generate），能量最小化出现上述问题。 后来，找人代做，在利用高斯软件优化配体结构计算resp电荷过程中，未加小分子六元环上的氢，在41个原子基础上计算的，后面用的acpype生成GAFF力场的配体拓扑文件，经过对比，各原子电荷与之前不同，itp文件中出现二面角的improper项，其他我没仔细看。然后接着跑em，未出现能量不收敛的情况，nvt、npt、md能跑下去，但是出现新的问题（box不断变大），我也发帖问各位老师了，就不在赘述。我后面再试试小分子计算去掉6个氢试试。以上是我的一点经验和看法，不确定是否正确，仅供参考。

Seyilaxa

chenwei 发表于 2024-12-2 22:53
你好，我目前关于本帖子的问题已经解决。主要问题出现在小分子的拓扑文件上。我先用sobtop生成GAFF力场下 ...

量化计算没有任何理由不加本应该有的氢，否则只会优化出不合理的结构
原来的问题多半是生成itp时部分键连关系或原子类型指认有误，导致的拓扑不合理

chenwei

Seyilaxa 发表于 2024-12-2 23:53
量化计算没有任何理由不加本应该有的氢，否则只会优化出不合理的结构
原来的问题多半是生成itp时部分键 ...

非常感谢您的指导，我再试试

暮晨

chenwei 发表于 2024-12-2 22:53
你好，我目前关于本帖子的问题已经解决。主要问题出现在小分子的拓扑文件上。我先用sobtop生成GAFF力场下 ...

谢谢楼主大大