计算化学公社

 找回密码 Forget password
 注册 Register
Views: 7340|回复 Reply: 23
打印 Print 上一主题 Last thread 下一主题 Next thread

[GROMACS] Auto-FEP for Gromacs(Shell 脚本)

[复制链接 Copy URL]

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

本帖最后由 Huschein 于 2023-8-1 19:45 编辑

Auto-FEP for Gromacs
引言:在分子动力学领域的学者都知道,自由能计算是一项非常有用的工具,它在药物设计、生物大分子相互作用分析等领域具有非常重要的地位。而自由能微扰(FEP)则是在经典物理领域(当然我们应该抛开量子化学去谈这件事)目前公认的最准确的计算自由能的方法,而Gromacs很早就实现了FEP的计算。虽然Gromacs为我们降低了很大的工作量,但是FEP所需要的批处理流程依旧十分繁杂,并且我最近的工作又恰好需要基于GMX实现的FEP,所以我在这里写了一些shell脚本来方便大家使用GMX做批处理的FEP计算。其中包括基于服务器的运算脚本(基于slurm调度系统)和基于个人电脑的脚本,但个人电脑脚本我尚未测试,如果出现问题请与我联系(邮箱:Hsuchein0126@outlook.com)。

Github链接:
https://github.com/Hsuchein/AutoFEP_For_GMX

准备工作:
1. Gromacs安装(版本新于5.0.0即可)
2. Alchemical-analysis安装,按照它官网教程即可(需要conda配置python2.7的环境),最后安装完后用我所提供的alchemical-analysis替换原文件的,因为源代码有问题,我已经修改过了。https://github.com/alchemistry/alchemlyb),在做所有工作之前,务必看清楚这里!!!GMX在FEP任务中最后一步跑MD的时候会对每个模拟都输出一个dhdl.xvg这个文件,而alchemical-analysis是一个用于分析dhdl.xvg这个文件的一个python脚本,所以如果你想要处理FEP的输出,你务必要自己去它的官网,去配置好它所需的环境。而我目前的脚本还不支持直接一键自动分析,我目前的脚本是帮助你快速地跑10段FEP模拟,不用你自己手动调mdp文件,一个一个递交任务,所以我这个脚本本质上就是一个遥控器,同时开启10个模拟,没有太多的技术含量。


操作流程:
注意:全部操作都在主文件夹中完成,不需要去子目录中执行脚本命令!!!
1. 主文件夹配置(真的重要,看清楚!!!)
将你配置好的文件(一般建议是Gromacs中完成Genion指令后得到的文件,在此之前你需要自己完成设置盒子大小,添加溶剂,平衡电荷等操作)。脚本的原理就是:你先自己建模,因为模型都是一样的,FEP改变的只是mdp文件,所以你先建模(盒子大小,溶剂,电荷),全部建完之后把下一步(能量最小化)所需的文件全部丢进一个新文件夹。也就是说,需要你把此前最后一步所得的.top .itp .ndx .mdp .gro(注意自己手动修改mdp ndx等文件以适配自己的任务)等文件全部拷贝到新的主目录中,这里切记.gro文件只能有一个,不然脚本无法识别使用哪个.gro文件,然后如果你的力场是自己提供的,请一并复制.ff文件(example中有例子)到主目录中,然后将全部.sh脚本复制到你的主目录中(如果你是个人PC,那就用PC脚本;如果是slurm调度的服务器,就用servers脚本)。然后先把全部的东西都丢到新文件夹(下文也称之为“主目录”)后,再运行脚本,因为每段模拟都需要自己的指令,所以prepare.sh这个脚本,也会对其他脚本做一些修改,所以如果你不把其他脚本都加到主目录,你可能会遇到:
1.提示你找不到相应的.sh脚本
2.你后续任务递交出问题
2. 运行脚本
提示:!!!对于windows用户,bash脚本不能直接在终端运行,必须安装git-bash后用git-bash终端打开命令行,然后再用命令行运行脚本!!!
i)PC的话是一键运行,但是本人没有测试过,不知道是否会报错。并且PC是跑完一段模拟才能跑下一段,效率会低一些。
ii)servers脚本分为多段,我将逐个介绍
Prepare.sh:在你当前主目录下配置10个窗口的子文件夹目录(我设定的窗口是10,当然可以自己改,注意10个连接窗口会有11个模拟要跑),然后依次修改mdp文件的参数,并且复制文件到子文件夹目录。
Minimization.sh: 依次运行10个窗口的能量最小化任务,这个模拟时间极短,所以就不提交到服务器后台运算了,直接在当前终端运算即可。
Npt.sh/npt_sub.shnpt_sub.sh会先检查每个窗口中min的模拟是否成功(是否输出min.gro),如果成功则会执行gmx grompp的任务,随后调用npt.sh来递交任务到服务器后台运算,这里我给每个任务都写了JOB NAME,可以通过查看后台进程来看是否都运行上了(实测的时候有几段模拟会递交不上),如果出现递交不上的情况,那么直接进入相应的子目录运行npt.sh脚本即可。注意:npt_sub.sh中递交任务的部分(第25行)可能需要手动修改以适配自己的任务,同时也注意自己修改npt.shmdrun的命令以发挥最大的效能。
Md.sh/md_sub.sh:与npt同理
Collect.sh:在主目录下运行,它会先检查当前目录是否存在DHDL FILES这个文件夹,如果存在则清空里面的内容(所以有需要注意备份),如果不存在则创建,然后从每个窗口子目录中复制dhdl文件到该目录下,最后按照alchemical-analysis的脚本分析即可。
规划:目前打算规划把analysis也做成一键分析,看看有无时间吧。
备注:
1.其实这个shell脚本很简单,大家都可以写,但是也可以方便很多不会shell脚本又想要批量化处理模拟的同行们。
2.如果模拟有出现错误,我的建议是先把文件(你的.gro .ndx .top)名字改成和我example一样,尤其是ndx,因为我也没打算做一个非常普适的成熟脚本,因为这将需要我考虑各种情况,消耗巨大精力

版本修改纪要:
2023/07/24:AutoFEP4GMX 1.0



图片1.png (16.21 KB, 下载次数 Times of downloads: 77)

图片1.png

图片2.png (31.86 KB, 下载次数 Times of downloads: 68)

图片2.png

评分 Rate

参与人数
Participants 4
威望 +1 eV +11 收起 理由
Reason
ysgwill + 2 赞!
zsu007 + 5 赞!
sobereva + 1
limalvis + 4 赞!

查看全部评分 View all ratings

21

帖子

0

威望

443

eV
积分
464

Level 3 能力者

24#
发表于 Post on 2024-10-24 12:46:22 | 只看该作者 Only view this author
本帖最后由 Jeffrey0521 于 2024-10-24 12:48 编辑

但我查到的文件大多都是fep+做的,这些能作为参考吗?想问你具体参考的是哪些文献呢,无论怎么操作他的首尾窗口的邻接矩阵有的重合度也不太好

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

23#
 楼主 Author| 发表于 Post on 2024-10-22 20:47:33 | 只看该作者 Only view this author
Jeffrey0521 发表于 2024-10-22 12:59
我意思是各个window的范德华、库仑、bonded的变化快慢是怎么确定的?

文献

21

帖子

0

威望

443

eV
积分
464

Level 3 能力者

22#
发表于 Post on 2024-10-22 12:59:25 | 只看该作者 Only view this author
Huschein 发表于 2024-10-21 20:03
我现在就是算绝对结合能 始末就是没有到出现
窗口你可以参考文献,现在一般主流是9个或者11个
梯度是什 ...

我意思是各个window的范德华、库仑、bonded的变化快慢是怎么确定的?

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

21#
 楼主 Author| 发表于 Post on 2024-10-21 20:03:03 | 只看该作者 Only view this author
Jeffrey0521 发表于 2024-10-21 16:47
博主您好,想问下这个fep的始末状态分别是什么呢?就是一个化合物凭空出现吗?另外,选择10个窗口的依据是 ...

我现在就是算绝对结合能 始末就是没有到出现
窗口你可以参考文献,现在一般主流是9个或者11个
梯度是什么?

21

帖子

0

威望

443

eV
积分
464

Level 3 能力者

20#
发表于 Post on 2024-10-21 16:47:51 | 只看该作者 Only view this author
博主您好,想问下这个fep的始末状态分别是什么呢?就是一个化合物凭空出现吗?另外,选择10个窗口的依据是什么呢,梯度又是参考什么来设置的呢?求解,十分感谢!!

19

帖子

0

威望

227

eV
积分
246

Level 3 能力者

19#
发表于 Post on 2024-9-14 11:02:19 | 只看该作者 Only view this author
Huschein 发表于 2024-9-14 02:25
md是正常动力学 sd是随机动力学,sd采样更广,这样可以覆盖窗口之间的空隙,你做完FEP是要算窗口之间重叠 ...

好的,谢谢您!

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

18#
 楼主 Author| 发表于 Post on 2024-9-14 02:25:09 | 只看该作者 Only view this author
yyf123 发表于 2024-9-9 20:45
我还想请教一下,integrator = md和sd的区别在哪里,对精度的影响大吗

md是正常动力学 sd是随机动力学,sd采样更广,这样可以覆盖窗口之间的空隙,你做完FEP是要算窗口之间重叠的,如果重叠的不够那就说明会有跳跃,最好是每个窗口能重叠上

19

帖子

0

威望

227

eV
积分
246

Level 3 能力者

17#
发表于 Post on 2024-9-9 20:45:01 | 只看该作者 Only view this author
我还想请教一下,integrator = md和sd的区别在哪里,对精度的影响大吗

19

帖子

0

威望

227

eV
积分
246

Level 3 能力者

16#
发表于 Post on 2024-9-9 19:33:21 | 只看该作者 Only view this author
博主您好,最近想用fep计算配体和蛋白的结合自由能,请问gromacs的fep可以使用gpu加速吗,速度会损失多少, 谢谢您

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

15#
 楼主 Author| 发表于 Post on 2023-7-30 18:30:52 | 只看该作者 Only view this author
xptracy 发表于 2023-7-30 09:28
我是在主目录下运行bash AutoFEP4PC.sh,子目录下没有min.tpr

留一个联系方式吧,qq或微信,我具体看看

49

帖子

0

威望

1094

eV
积分
1143

Level 4 (黑子)

14#
发表于 Post on 2023-7-30 09:28:47 | 只看该作者 Only view this author
Huschein 发表于 2023-7-30 09:10
你要去各个Windows窗口里看有没有min.tpr,不是在主目录里看

我是在主目录下运行bash AutoFEP4PC.sh,子目录下没有min.tpr

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

13#
 楼主 Author| 发表于 Post on 2023-7-30 09:10:44 | 只看该作者 Only view this author
xptracy 发表于 2023-7-29 22:49
嗯嗯,我把example的都放进来了,运行也是这两个错误
AutoFEP_For ...

你要去各个Windows窗口里看有没有min.tpr,不是在主目录里看

49

帖子

0

威望

1094

eV
积分
1143

Level 4 (黑子)

12#
发表于 Post on 2023-7-29 22:49:51 | 只看该作者 Only view this author
Huschein 发表于 2023-7-29 22:11
你去对应的子目录查有没有min.tpr这个文件,你这压根就没有执行minimization,2021版本我觉得gmx mdrun - ...

嗯嗯,我把example的都放进来了,运行也是这两个错误
AutoFEP_For_GMX-main里没有min.tpr,21版之前数据是有生成min.tpr的

214

帖子

4

威望

1462

eV
积分
1756

Level 5 (御坂)

11#
 楼主 Author| 发表于 Post on 2023-7-29 22:11:51 | 只看该作者 Only view this author
xptracy 发表于 2023-7-29 21:22
请教下您是这些文件吗,还有缺什么吗,有两个错误?

alchemical-analysis替换的原文件是在哪个位置 ...

你去对应的子目录查有没有min.tpr这个文件,你这压根就没有执行minimization,2021版本我觉得gmx mdrun -deffnm是可以的,难道说非要-s指定tpr? 你先看看有没有min.tpr吧,有的话可能我要改一下脚本适配更低版本的gmx,如果没有的话那根本就是你在prepare的时候出了问题,我的建议是,你自己先拿我的example去跑,确定是版本或者代码的问题再找我

手机版 Mobile version|北京科音自然科学研究中心 Beijing Kein Research Center for Natural Sciences|京公网安备 11010502035419号|计算化学公社 — 北京科音旗下高水平计算化学交流论坛 ( 京ICP备14038949号-1 )|网站地图

GMT+8, 2026-2-17 17:31 , Processed in 0.191788 second(s), 25 queries , Gzip On.

快速回复 返回顶部 返回列表 Return to list