对比ORCA6.0和Gaussian16的几何优化算法

ZetaFunction

之前在用ORCA算ONIOM的时候就感觉几何优化不如Gaussian稳健，差不多大小的体系ORCA经常需要比Gaussian多好几倍的优化轮次才能收敛。在Gaussian中，几何结构难收敛一般体现为能量反复震荡，而ORCA倒是少有出现震荡，但是却经常出现每轮能量下降幅度过于微小，以至于看上去四个收敛指标快要达成了但实际还要再优化几十轮甚至上百轮才能真正收敛。虽然在很多计算级别下，尤其是纯泛函，ORCA单点计算速度远超Gaussian，但是如果考虑额外所需的优化轮次，在几何优化任务中哪个效率更高就不好说了。

最近更新的ORCA6.0号称又提示了优化算法的效率，但是我实际试用下来还是经常出现能量龟速下降的情况，正好最近要优化一个含铁-镍核心的体系，可以用来对比一下两个程序优化算法的效率。体系总计122个原子，自旋极化单重态，用uTPSSh-D3优化，铁镍核心和直接与之配位的基团用def2-TZVP基组，其余用def2-SVP基组，优化时有十个原子的坐标被冻结。波函数初猜由Gaussian产生，添加“stable=opt”以确保收敛到稳定的波函数，用Multiwfn转换为.gbw格式以供ORCA读取。所有的能量和梯度计算均在ORCA中完成，ORCA使用内置的默认优化器优化，Gaussian通过将Gaussian与ORCA联用搜索过渡态、产生IRC、做振动分析 - 思想家公社的门口：量子化学·分子模拟·二次元 (sobereva.com)中的连接脚本接受ORCA计算的梯度信息用以优化，收敛标准均使用两个程序的默认收敛标准。

ORCA在第307轮时终于达到了收敛标准，在服务器上32核并行用时15小时40分钟，最终能量为-6487.37397 A.U.，体系在第200轮的能量就来到了-6487.37377 A.U.，也就是说最后一百多轮能量只下降了0.1255 kcal/mol。而Gaussian仅用73轮就收敛了，总用时3小时30分钟，最终能量为-6487.37399 A.U.，与ORCA所优化的结构的差异肉眼不可见。

可见，即便更新了6.0版本，ORCA的几何优化算法依然不如Gaussian高效稳健，所需轮次甚至可能到达Gaussian的数倍。最具效率的做法是还是通过连接脚本，各取所长，ORCA负责波函数计算，Gaussian负责几何优化，这也是目前我依然不能完全用ORCA代替Gaussian的重要原因之一。

wjc404

几何优化的算法挺复杂的，不仅涉及到优化算法本身，还有二阶导矩阵的更新算法。可以看看calcall时两者的步数差别。

大王来何

优秀

wzkchem5

可否给一下相应的输入文件？我们测一下，看看orca优化算法是否有什么地方需要改进的。
最好也可以给一下高斯的输出文件

ZetaFunction

wzkchem5 发表于 2024-8-3 16:19
可否给一下相应的输入文件？我们测一下，看看orca优化算法是否有什么地方需要改进的。
最好也可以给一下高 ...

通过百度网盘分享的文件：SIa-S1-test.rar
链接：https://pan.baidu.com/s/1PIAFZgB4mSlgDlt6LVv5lA?pwd=0803
提取码：0803
--来自百度网盘超级会员V5的分享

SIa-S1_gau.gbw是在Gaussian中计算的用作初猜的轨道

wzkchem5

ZetaFunction 发表于 2024-8-3 12:04
通过百度网盘分享的文件：SIa-S1-test.rar
链接：https://pan.baidu.com/s/1PIAFZgB4mSlgDlt6LVv5lA?pwd ...

提示页面已过期，下载不了

lanthanum

wzkchem5 发表于 2024-8-4 00:08
提示页面已过期，下载不了

SIa-S1-test.7z.001 (9 MB, 下载次数 Times of downloads: 28)

2024-8-4 00:31 上传 Uploaded

点击下载
Click to download

SIa-S1-test.7z.002 (7.78 MB, 下载次数 Times of downloads: 31)

2024-8-4 00:31 上传 Uploaded

点击下载
Click to download

wzkchem5

ZetaFunction 发表于 2024-8-3 12:04
通过百度网盘分享的文件：SIa-S1-test.rar
链接：https://pan.baidu.com/s/1PIAFZgB4mSlgDlt6LVv5lA?pwd ...

优化的结构做频率分析了吗，有没有虚频？
由ORCA输出文件里的“Lowest eigenvalues of augmented Hessian”可以看出，ORCA的incrementally built Hessian一直有一个小虚频，而高斯的没有。ORCA后期收敛缓慢可能是因为ORCA想尝试消掉这个小虚频。所以现在问题就是这个虚频客观上是不是存在，如果不存在，那就说明是ORCA的默认初始Hessian或者Hessian update的算法不如高斯好；如果存在，那就说明高斯把问题想简单了，没有意识到这个体系有一个虚频需要消。

喵星大佬

倪都在起始结构做个振动分析，用精确的其实hessian去优化，看看着两个程序更新hessian的步骤和走步线搜索算法的区别

ZetaFunction

wzkchem5 发表于 2024-8-4 03:29
优化的结构做频率分析了吗，有没有虚频？
由ORCA输出文件里的“Lowest eigenvalues of augmented Hessia ...

学校的服务器在排队，频率还没算。不过我在使用Gaussian调用ORCA的时候还发现有时在优化的最后几步会出现梯度还在逐渐降低但是能量居然会略有上升的情况，我怀疑可能是因为ORCA用的是自适应格点所以Gaussian每次调用时的积分格点分配可能会不同导致影响计算结果。ORCA中有无可以固定格点数量的关键词？

wzkchem5

ZetaFunction 发表于 2024-8-4 09:06
学校的服务器在排队，频率还没算。不过我在使用Gaussian调用ORCA的时候还发现有时在优化的最后几步会出现 ...

梯度降低本来就不代表能量也必须降低啊

ZetaFunction

wzkchem5 发表于 2024-8-4 16:16
梯度降低本来就不代表能量也必须降低啊

我还发现ORCA的QM/XTB每次在Preparing阶段不仅会用XTB计算低层的电荷，还会用XTB计算高层的电荷，而且似乎都是先用XTB优化结构再计算电荷，第一次计算电荷结束后几乎就不会随着优化更新，这导致如果几何结构变化较大，拿最后一帧重新计算单点根本得不到一致的结果。低层的电荷可以指定ORCAFFfile来固定，但是高层的电荷我没找到可以关闭自动计算的方法。

wzkchem5

ZetaFunction 发表于 2024-8-4 09:41
我还发现ORCA的QM/XTB每次在Preparing阶段不仅会用XTB计算低层的电荷，还会用XTB计算高层的电荷，而且似 ...

这个应该是有意为之，如果不这么做的话，就需要计算XTB电荷对结构的导数，才能做结构优化，但是XTB电荷对结构的解析导数实现起来比较麻烦，也会增加不少计算量

ZetaFunction

wzkchem5 发表于 2024-8-4 17:43
这个应该是有意为之，如果不这么做的话，就需要计算XTB电荷对结构的导数，才能做结构优化，但是XTB电荷对 ...

我之前就一直好奇ONIOM如果考虑静电嵌入的话梯度要怎么处理，因为高层的能量不再是和低层坐标无关的函数了，像Gaussian一样对低层用固定点电荷也许可以获得解析梯度，我怀疑ORCA的QM/XTB也是基于第一帧计算的点电荷，在优化中不会更新。但是这么做的一大问题就是最后一帧如果重新算单点无法重现优化的结果，虽然差异很微小，但是会导致频率分析或IRC计算等需要计算级别严格一致的任务出问题。根据我的尝试，就算指定了相同的ORCAFFfile也不行，因为ORCA还是会在计算开始前去拟合高层的点电荷用来处理共价键截断。要么ORCA可以设法同时指定高层和低层的电荷，要么在需要求过渡态频率或者IRC的时候只能用Mechanical Embedding。

ZetaFunction

wzkchem5 发表于 2024-8-4 03:29
优化的结构做频率分析了吗，有没有虚频？
由ORCA输出文件里的“Lowest eigenvalues of augmented Hessia ...

频率算完了，Gaussian所优化的结构确实有一个-9的小虚频，当前计算级别下算一次频率要三个小时。ORCA所优化的结构还没算频率，二者差异极其微小。即使ORCA的结构没有虚频，但是优化所额外需要的轮次还是超出Gaussian太多了，这么看来也许可以先用Gaussian优化后再用ORCA优化，应该能节约不少时间。