|
|
本帖最后由 Freeman 于 2022-7-31 17:00 编辑
【2022.07.31更新】
原来的脚本功能单一，代码不规范，应对不同的场景还要用户自己改代码，十分不便。现在gau_bag.sh升级重写了。细节如下。
gau_bag.sh
(6.73 KB, 下载次数 Times of downloads: 50)
改进
1、原来的脚本只能应用于一阶导数的场景，现在可以应用于二阶导数的场景了。用GV观看CASPT2水平的振动模式，还能读各种热力学量，非常美滋滋。
2、原来BAGEL的输出，后一帧会把前一帧覆盖掉，没法复盘整个优化过程，现在每一帧的输出文件都会保存。
3、原来的脚本每次用都要改动以适应不同的核数、基组、活性空间、电子态的需求。现在这些信息都转到从gjf文件里的external的引号里读取了，用户就可以把脚本放到一个固定的地方而不用老是改了。
4、原来的Gaussian的输入文件里除了能量、梯度什么也没有，要看多态的能量和组态还要诉诸BAGEL的输出文件，很麻烦。现在BAGEL输出文件的重要信息也会打印到Gaussian输出文件里，用户就不用辗转于两个输出文件之间了。
5、原来的代码乱七八糟的，不规范，现在规范了很多，方便用户阅读修改。
案例：基态环丁二烯互变异构过渡态
这里用环丁二烯的例子讲述新的gau_bag.sh的用法。单重态环丁二烯的极小点是长方形，其长短边互变异构的过渡态是正方形。（图一）该过渡态是多参考体系，正好用CASPT2优化。
首先，创建并进入一个新文件夹c4h4s，本节以下所有活动都在这个文件夹里进行。建模画出正方形的初猜结构，用BAGEL算个RHF单点，以提供后续CASSCF的活性空间。（文件c4h4s.bag.json）
- {"bagel":[
- {
- "title":"molecule",
- "basis":"cc-pvdz",
- "df_basis":"cc-pvdz-jkfit",
- "angstrom" : "true",
- "geometry":[
- { "atom" : "C", "xyz" : [ 0.00000000, 1.01527700, 0.00000000]},
- { "atom" : "C", "xyz" : [ 1.01527700, 0.00000000, 0.00000000]},
- { "atom" : "C", "xyz" : [ -1.01527700, 0.00000000, 0.00000000]},
- { "atom" : "C", "xyz" : [ 0.00000000, -1.01527700, 0.00000000]},
- { "atom" : "H", "xyz" : [ 0.00000000, 2.09288200, 0.00000000]},
- { "atom" : "H", "xyz" : [ 2.09288200, 0.00000000, 0.00000000]},
- { "atom" : "H", "xyz" : [ 0.00000000, -2.09288200, 0.00000000]},
- { "atom" : "H", "xyz" : [ -2.09288200, 0.00000000, 0.00000000]}
- ]
- },
- {
- "title":"hf"
- },
- {
- "title":"save_ref", // 保存波函数文件为c4h4s.bag.archive
- "file":"c4h4s.bag"
- },
- {
- "title":"print", // archive文件没法可视化，所以还要输出molden用来看轨道
- "file":"c4h4s.bag.molden",
- "orbitals":true
- }
- ]}
- {"bagel":[
- {
- "title":"load_ref", // 从c4h4s.bag.archive读取前一步的分子坐标、基组和轨道
- "file":"c4h4s.bag"
- },
- {
- "title":"casscf",
- "nact":4, // 四个活性轨道
- "nclosed":12, // 满占据的非活性轨道数
- "active":[13,14,15,20] // 活性轨道的编号，从Multiwfn看molden文件而得来的
- },
- {
- "title":"print",
- "file":"c4h4s+.bag.molden",
- "orbitals":true
- },
- {
- "title":"save_ref",
- "file":"c4h4s+.bag"
- }
- ]}
- %nprocshared=1 // 只用一个核就可以了，因为密集的计算都是BAGEL做的
- %chk=c4h4s.bag.chk // 保存最后的频率信息到chk文件里
- # opt=(ts,noeigen,calcfc,nomicro) freq nosym // nomicro和nosym关键词一定要加
- external='bash ~/gau_bag.sh 9 cc-pvdz cc-pvdz-jkfit 2 0 4 12' // 调用在家目录下的gau_bag.sh脚本，给BAGEL用9个核，基组用cc-pvdz，辅助基组用cc-pvdz-jkfit，算两个（单重）态（其实算一个态就可以，这里是为了举例），几何优化针对第1个态（0是第一个态，1是第二个态，以此类推），有4个活性轨道、12个满占据的非活性轨道
- Title Card Required
- 0 1 // 以前的脚本是不从gjf读电荷数和电子多重度的，而是要直接改写脚本内容；现在的脚本是从这里读的，所以不要乱设了
- C 0.00000000 1.01527700 0.00000000
- C 1.01527700 0.00000000 0.00000000
- C -1.01527700 0.00000000 0.00000000
- C 0.00000000 -1.01527700 0.00000000
- H 0.00000000 2.09288200 0.00000000
- H 2.09288200 0.00000000 0.00000000
- H 0.00000000 -2.09288200 0.00000000
- H -2.09288200 0.00000000 0.00000000
******************************************************************************
BAGEL自己的优化算法不太行，经常是好不容易优化到平衡位置旁边了，突然一抖，能量上升0.0几个Hartree，就前功尽弃了，对低血压有很强的治疗效果。一个办法是降低置信半径至0.1左右并固定，这样优化时就不会振荡，但是会增加无意义的耗时。Gaussian的优化算法倒是很好，因此笔者写了个将Gaussian和BAGEL联用的脚本，将BAGEL的CASPT2解析梯度和Gaussian的优化算法结合起来。
Gaussian external接口脚本的基本思路
1、Gaussian把这一步的分子坐标、电荷、自旋多重度和要算的导数阶数（优化是一阶导数，频率是二阶导数）保存在一个文件里，文件名是Gau-[PID编号].EIn。
2、*读取EIn里的信息，生成其他程序的输入文件Gau-[PID编号].inp；
3、*设置环境变量，运行其他程序，获取其输出文件Gau-[PID编号].out；
4、*把out里的能量和梯度（或频率）信息复制到Gau-[PID编号].EOu里去；（BAGEL的能量和梯度信息可以另外分别保存在ENERGY.out和FORCE_0.out文件里，更方便读取）
5、Gaussian读取EOu里的信息，判断下一步的分子坐标该变成什么样；
6、如果是优化或IRC任务，重复1-5，直到收敛；如果是频率任务，就不用重复了。
以上打*的部分是接口脚本里要写的内容。每一步，Gaussian都会调用一次这个脚本，调用方式类似于“bash script.sh Gau-[PID编号].xxx（oniom里面用的，与本文无关） Gau-[PID编号].EIn Gau-[PID编号].EOu”，也就是说$2默认就是Gau-[PID编号].EIn，而$3默认就是Gau-[PID编号].EOu。
脚本使用流程
1、用现有的分子坐标编写一个用BAGEL的CASSCF算单点能的输入文件.json，如有必要就调整轨道顺序，注意要让它输出molden文件和archive文件，后者相当于Gaussian的chk文件；
2、用Multiwfn查看molden的轨道合不合心意，没问题的话就保存成.gjf格式；
3、改写.gjf文件，设置
- %nproc=1
- # opt(nomicro) external='bash gau_bagel.sh' nosym
5、运行Gaussian。
随后，脚本会创建一个名为“Gau-[PID编号]”的文件夹，然后把第一步的reference.archive文件复制进去，命名为Gau-[PID编号].archive。之后每走一步，都会读取上一步的archive作为轨道初猜，走完一步就更新archive里的信息，并且输出ENERGY.out和FORCE_0.out文件。最后在文件夹外面出现了Gaussian的输出文件。
效果图
先运行，获得初始molden和archive，然后运行带有这个脚本的gjf。
P.S. BAGEL的输入文件格式非常变态，大家写输入文件的时候要注意一下方括号、大括号、逗号这些乱七八糟的东西有没有写全。
下期预告
测评一下包括BAGEL在内的常用来做多参考计算的程序在不同进程数下的速度。
|
评分 Rate
-
查看全部评分 View all ratings
|
发表于 Post on 2021-12-23 15:25:59
收藏 Add to favorites
楼主 Author