脚本分享：sobtop联合Gromacs对天然产物进行分子动力学构象搜索

欢乐多

脚本分享：sobtop联合Gromacs对天然产物进行分子动力学构象搜索

前言

上一篇博文在《Windows系统Gromacs对有机小分子进行分子动力学构象搜索》略有不便，因此本次脚本是在Linux系统Centos7.6 和9中使用的，在 1sobtop.zip (8.71 MB, 下载次数 Times of downloads: 115)

2023-12-8 05:29 上传 Uploaded

点击下载
Click to download

压缩包中名为00_sob_autoall.sh。脚本主要进行了三步计算：1，用Gaussian计算电荷；2，用Sobtop产生top文件；3，Gromacs进行动力学构象搜索。如何使用脚本，B站视频演示《sobtop联合Gromacs对天然产物进行分子动力学构象搜索》。非常欢迎大家使用交流。

1. 首次使用所需要进行的基本设置

Linux系统所需要的设置

Sobtop的安装使用：http://sobereva.com/soft/Sobtop/

Gaussain的安装：《Gaussian的安装方法及运行时的相关问题》

Gromacs程序的安装： 参考社长大人Sobereva老师博文http://sobereva.com/458

脚本使用Gaussain16，Gromacs 2022.6，Sobtop 1.0(dev3.1)

以上软件加入到Linux系统的环境变量。

Window系统所需要的设置:

安装VMD，并将其加入到系统环境变量，以便用命令行调用VMD

VMD下载地址https://www.ks.uiuc.edu/Development/Download/download.cgi?PackageName=VMD

下载Windows64-bit的文件，安装即用。

2. 使用方法

解压压缩包，进入目录，赋予可执行权限，得到以下目录：

修改RESP2_noopt.sh文件中12行为自己系统的Gaussian版本。脚本使用的是g16。

将自己结构的*.mol2格式的文件放入该目录中，需要Gaussian View 保存的.mol2的文件，比如abc.mol2，文件名可以任意，后缀名需要是.mol2，可以放多个, 比如1.mol2，2.mol2，3.mol2。脚本是一个for循环，会遍历当前目录下所有的mol2格式的文件。

运行脚本

./ 00_sob_autoall.sh

将会自动进行以下三步计算：

1，用Gaussian计算电荷；

2，用Sobtop产生top文件；

3，Gromacs进行动力学构象搜索。

md01.mdp

md02.mdp

md03.mdp

以上3个文件，是动力学的输入文件，三个文件中的温度不同，将会进行三次动力学搜索。

脚本还对运行的时长进行了记录，如果您的系统可以发邮件，修改脚本中最后一行为自己的邮箱，模拟结束，还会将运行的时长发送到您的邮箱，进行提醒，以便及时处理计算结果，非常方便。

运行结束后，将相应的*_md01.gro, *_md01.xtc文件拖到Window系统，用VMD打开即可。

3. Gromacs构象搜索结果处理

打开VMD，载入abc_md01.gro，在VMD main界面—file—loaddata into molecule—Browse—选择abc_md01.xtc。

载入构象搜索的运动轨迹之后。在VMD main界面—file—save coordinates，选择.xyz格式文件。将轨迹保存成*_traj.xyz ，比如2R4S_traj.xyz。后缀名“_traj.xyz”是必须的以便博文《chatGPT辅助生成molclus懒人脚本：一键完成对复杂天然产物NMR和ECD计算》使用，由此我们得到2R4S_traj.xyz文件。

只保存其中的一个构象，last选0，保存为*_traj2.xyz，比如2R4S_traj2.xyz，后缀名“_traj2.xyz”是必须的，以供博文《chatGPT辅助生成molclus懒人脚本：一键完成对复杂天然产物NMR和ECD计算》进行甲基的自动识别。

压缩包中还提供一个Window系统04showgro.bat文件，省事的话，将VMD加入到系统环境变量，将*_md01.gro, *_md01.xtc与之放入同一个文件夹，双击，即可依次打开相应的*_md01.gro, *_md01.xtc，非常方便，

4. Gromacs构象搜索细节

对md01.mdp, md02.mdp, md03.mdp中动力学模拟时间和温度调节，使得模拟更充分，尽可能获得较为全面的构象。

如果发现搜索到的构象不够理想，分子运动范围不充分，

三个温度可以进行设定，可以更改其中的温度和步长。

nsteps = 500000

annealing_temp = 0 300 300 0

运行

./ 04chg.sh

这个脚本是比./ 00_sob_autoall.sh少一步计算原子电荷，因此运行该脚本，目录中需要包含相应的*.chg，比如需要同时含有abc.mol2和abc.chg。

模拟结束，要观察结构是否变形，一般认为温度高一些获得构象更为充分，但是如果温度设置过高，结构将会裂解。请自己掌握好温度和步长的设定。

请各位多多指教，谢谢大家！

Duoo

TJ老师您好,我在尝试使用你分享的脚本的时候,发现第一步计算迟迟没有结束,即用Gaussian计算电荷.目前已经过去了接近两个小时了,请问我需要修改配置么?

sob top.ini如下:
  nthreads= 4  // Number of threads used for parallel calculation
  iskipgendih= 0  // 1: Skip generating dihedral terms, 0: Do not skip
  ioutatminfo= 0  // 1: Output atomic coordinates and connectivities when to atminfo.txt in current folder when loading input file, 0: Do not output to file but shown on screen
  ichggeom= 1  // 1: When loading .chg file, replace current geometry with that in .chg file, 0: Do not replace
  k_method= 2  // Default method of determining k. 1: Seminario, 2: mSeminario, 3: m2Seminario, 4: DRIH
  bondcrit= 1.15  // When pdb/pqr is used as input, two atoms are considered as bonded if their distance is smaller than sum of their covalent radii multiplied by this factor. Priority is lower than the criteria defined in bondcrit.dat
  Multiwfn_cmd= "/home/yyp/opt/Multiwfn_3.8_dev_bin_Linux_noGUI/Multiwfn_noGUI" // Path of executable file of Multiwfn
  OpenBabel_cmd= "/usr/bin/obabel"  // Path of executable file of OpenBabel

欢乐多

Duoo 发表于 2024-10-3 03:55
TJ老师您好,我在尝试使用你分享的脚本的时候,发现第一步计算迟迟没有结束,即用Gaussian计算电荷.目前已经过 ...

可以用这个脚本 RESP2_noopt.sh (2.98 KB, 下载次数 Times of downloads: 6)

2024-10-5 03:54 上传 Uploaded

点击下载
Click to download

替换一下，原始脚本高斯要在气相和溶剂中计算两次电荷，如果想要快速，只计算气相条件的电荷也可以。

sobereva

Duoo 发表于 2024-10-3 19:55
TJ老师您好,我在尝试使用你分享的脚本的时候,发现第一步计算迟迟没有结束,即用Gaussian计算电荷.目前已经过 ...

RESP2_noopt.sh都还没跑完，用sobtop是后话。另外，不叫“用Gaussian计算电荷”，而是“基于Gaussian产生的信息通过Multiwfn计算电荷”，仔细看下文搞清楚原子电荷到底是什么程序算的、计算过程是什么
RESP2原子电荷的思想以及在Multiwfn中的计算
http://sobereva.com/531（http://bbs.keinsci.com/thread-16190-1-1.html）

没具体体系细节、多少核并行，别人没法判断为什么计算时间很长

Duoo

欢乐多 发表于 2024-10-5 03:57
可以用这个脚本替换一下，原始脚本高斯要在气相和溶剂中计算两次电荷，如果想要快速，只计算气相条件的电 ...

好的,谢谢老师.我用的是学院的超算,可能没有并行运算直接在用户目录下运算下了.我再研究研究.

Duoo

sobereva 发表于 2024-10-5 04:04
RESP2_noopt.sh都还没跑完，用sobtop是后话。另外，不叫“用Gaussian计算电荷”，而是“基于Gaussian产生 ...

好的,谢谢社长!!! 我再多学习学习

Duoo

似乎是因为我用的是超算平台,所以得用slurm脚本提交任务,而不是./ 00_sob_autoall.sh提交脚本.我也将超算平台版的00_sob_autoall.sh脚本上传供各位大佬参考吧