计算化学公社

 找回密码 Forget password
 注册 Register
Views: 9752|回复 Reply: 9

[综合交流] 批量转化二维分子结构到量化输入文件

[复制链接 Copy URL]

312

帖子

1

威望

3106

eV
积分
3438

Level 5 (御坂)

发表于 Post on 2020-1-18 11:18:13 | 显示全部楼层 Show all |阅读模式 Reading model
本帖最后由 exity 于 2021-2-2 11:53 编辑

有的同学可能没有正版的chemoffice套件,又因为洁癖等原因不愿意安装DAO版,那么Chemsketch软件有免费版,可以生成cdx文件,理论上应该可以在这里替代chemoffice。--24/Feb/2020


关于open babel的使用请参考社长的帖子:http://sobereva.com/440
==================================================================

最近做了一个课题,考察不同官能团对蒽环结构上的电子云影响(其实是我瞎编的,并没有做这个课题)。
利用chemoffice里的套件,随手就画了若干种可能性,如下图所示(随手画的,不要抬杠):
1.png

再想到还有其他各种官能团要考察(硝基,胺基,甲基,磺酸基等等等),排列组合,要考察的对象数量相当可观,脑补了一下,心态爆炸。
这TMD的一个一个存为单独的文件,再用chem3D打开,另存为gjf,再写关键词,岂不是要搞到明年今日?
冷静了一下,摸索出一条路径,相当方便,在此分享给大家,一起进步。


1. 打开你的chemdraw,尽情的绘制你的分子,全部画在一起,就和上图一样,然后另存为1.cdx。
2. 登录你的服务器,安装epel源,然后yum install openbabel.x86_64,安装openbabel,视情况使用sudo权限。
3. 输入命令 obabel,看看你安装好了没,if安装好了,就继续,没安装好,自己想办法去。
4. 把1.cdx上传到你的服务器某目录,这里我放在/home/username/temp/obabeltest/目录下,放哪里随你,能找到就行。
5. 在此目录下输入,obabel 1.cdx -O new.gjf -m --gen3D --minimize --ff UFF,然后稍等,视你cdx文件中的分子数量多少耗时不尽相同。
6. 在此目录下再次输入命令 touch keyword.sh,此命令的意思是创建一个名为keyword.sh的脚本文件。
7. 用vim等命令打开keyword.sh,写入如下内容后保存:(复制粗体的部分,上下空行不要复制)

#!/bin/bash
    for i in *.gjf
    do
    echo Decent G09 inputfile of ${i} is being prepared ...
    sed -i "s/#Put Keywords Here, check Charge and Multiplicity./%nproc=48\n%MEM=185GB\n%chk=${i//gjf/chk}\n#p freq opt scrf=(smd,solvent=methanol) EmpiricalDispersion=GD3 int=ultrafine M062X\/6-311G\**/g" ${i}
    echo Your refined ${i} was prepared!
    echo
    done


8. 用chmod +x keyword.sh命令给此脚本运行权限后执行。
9. ls看看目录下的文件,检查一下关键词和电荷等内容,就可以啦。
10. 谢谢收看,祝新年快乐!

========================================
刚刚发现,直接用chemdraw存的cdx文件进行转化的话,容易翻车,见下图
2.png
一般有chemdraw的人都会有chem3D,在chemdraw绘制好后用chem3D打开一下后再另存为cdx文件就不会有这样的问题了,见下图:
3.png

当然了,你在chem3D里跑下能量最小化也行,但是有时候反而会把有序的分子排序搞乱搞错,注意一下这一点就好了。



===============================把两个脚本文件整理了一下============================
使用obabel脚本如下粗体所示:

time obabel *.cdx -O new.gjf -m --gen3D --minimize --ff UFF

  • 解释如下:
  • time         对obabel程序的转化过程进行计时
  • obabel       调用obabel程序
  • *.cdx         表示使用当前目录下的cdx文件作为obabel的输入文件
  • -O             obabel程序的输出参数,后面接具体的文件名,例如 a.gjf等
  • new.gjf      上述“-O”参数的具体输出的文件名,在我们这里的实践中,则会输出new1.gjf new2.gjf new3.gjf...等等
  • -m             obabel程序参数,表示把一个输入文件,分割成多个输出文件,每个独立的分子构型都写入到一个独立的文件中
  • --gen3D     obabel程序参数,表示生成三维构型的分子坐标,如果不加这个参数,用可视化软件打开转化后的gjf文件,会发现所有的原子都在一个平面上,而且没有H原子
  • --minimize   表示进行能量最小化
  • --ff UFF      表示用UFF分子力场进行能量最小化,我的分子比较特殊,在UFF下产生的构型比较好,其他同学可以使用MMFF94力场,可能效果更好,即--ff MMFF94就可以了

========================================


写入keyword的脚本如下粗体所示(太麻烦了,就不解释了,红色部分是你自己能定义的部分,具体意义看看就明白了,都是高斯输入文件的东西)


#!/bin/bash
    for i in *.gjf
    do
    echo Decent G09 inputfile of ${i} is being prepared ...
    sed -i "s/#Put Keywords Here, check Charge and Multiplicity./%nproc=48\n%MEM=185GB\n%chk=${i//gjf/chk}\n#p freq opt scrf=(smd,solvent=methanol) EmpiricalDispersion=GD3 int=ultrafine M062X\/6-311G\**/g" ${i}
    echo Your refined ${i} was prepared!
    echo
    done


========================================



评分 Rate

参与人数
Participants 7
威望 +1 eV +26 收起 理由
Reason
harrisonmo + 3 谢谢
ggdh + 5 赞!
Satoru + 5 谢谢
pika02 + 5 谢谢茄子
Novice + 3 好物!
sobereva + 1
zsu007 + 5 赞!

查看全部评分 View all ratings

1020

帖子

0

威望

3508

eV
积分
4528

Level 6 (一方通行)

发表于 Post on 2020-1-18 12:02:14 来自手机 | 显示全部楼层 Show all
要是机器学习手写结构式转为gjf文件就牛p了

312

帖子

1

威望

3106

eV
积分
3438

Level 5 (御坂)

 楼主 Author| 发表于 Post on 2020-1-18 13:10:45 | 显示全部楼层 Show all
granvia 发表于 2020-1-18 12:02
要是机器学习手写结构式转为gjf文件就牛p了

这牵涉到图形识别了,期待有高人出现造福大家

57

帖子

1

威望

2227

eV
积分
2304

Level 5 (御坂)

发表于 Post on 2020-1-19 00:04:01 | 显示全部楼层 Show all
本帖最后由 阿甘 于 2020-1-19 00:07 编辑

这个问题还可以利用内坐标来造这个结构,其中蒽环就相当于一个framework,让位置对应到1-10,取代基对应到不同的整数,编个程序写内坐标一劳永逸~~

顺便推荐一下杨伟涛老师的Best-First-Search算法~~~

261

帖子

0

威望

2073

eV
积分
2334

Level 5 (御坂)

计算化学路人甲

发表于 Post on 2020-1-19 10:43:34 | 显示全部楼层 Show all
友情提示,即使使用Chem3D转存一下可能也会出现分子结构键连出问题的现象,建议存完后确认一下。

评分 Rate

参与人数
Participants 1
eV +1 收起 理由
Reason
exity + 1 确实!还是应该把生成的gjf手动检查一遍!

查看全部评分 View all ratings

1020

帖子

0

威望

3508

eV
积分
4528

Level 6 (一方通行)

发表于 Post on 2020-1-19 12:35:59 来自手机 | 显示全部楼层 Show all
阿甘 发表于 2020-1-19 00:04
这个问题还可以利用内坐标来造这个结构,其中蒽环就相当于一个framework,让位置对应到1-10,取代基对应到 ...

Best 还是 Breadth ?

57

帖子

1

威望

2227

eV
积分
2304

Level 5 (御坂)

发表于 Post on 2020-1-19 20:37:27 | 显示全部楼层 Show all

261

帖子

0

威望

3756

eV
积分
4017

Level 6 (一方通行)

发表于 Post on 2020-1-30 21:14:17 | 显示全部楼层 Show all
请教一下楼主脚本中${i//gjf/chk}这里是什么含义呢?

183

帖子

4

威望

1457

eV
积分
1720

Level 5 (御坂)

发表于 Post on 2020-1-30 21:17:22 | 显示全部楼层 Show all
函数与激情 发表于 2020-1-30 21:14
请教一下楼主脚本中${i//gjf/chk}这里是什么含义呢?

改后缀。比如1.gjf变为1.chk

261

帖子

0

威望

3756

eV
积分
4017

Level 6 (一方通行)

发表于 Post on 2020-1-30 21:18:39 | 显示全部楼层 Show all
liuyuje714 发表于 2020-1-30 21:17
改后缀。比如1.gjf变为1.chk

十分感谢!

本版积分规则 Credits rule

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

GMT+8, 2023-2-2 22:40 , Processed in 0.202784 second(s), 25 queries .

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