波函数稳定性分析的计算时间

zha23

我在用我的程序使用6-31g*计算9个碳的正烷烃使发现，如果使用SAD作为初始猜测，会收敛到错误的波函数，这意味这我需要加入稳定性分析来检测稳定性。我阅读了相关的文献遇到了一个问题：波函数稳定性分析，需要计算(ia||jb)这样的电子积分，其时间复杂度为N^5，通过直接方法使用Davidson对角化，可以将时间复杂度降为N^4，但是需要迭代求解，每次都要计算电子积分，如此计算是否过于耗时？

常轩豪

L1yuu 发表于 2024-10-23 09:25
感谢您的回复，最后还是看贴主 @zha23 在二楼发的这篇 J. Chem. Phys. 1993, 99 (2), 1262–1270 了解了 ...

嗯嗯, 这篇文章里的 11-15 确实是这个问题direct起来实现的核心

L1yuu

常轩豪 发表于 2024-10-23 09:05
波函数的稳定性理论要追溯到 Cizek J., Paldus J. 的工作 [J. Chem. Phys., 1967, 47(10), p.3976–3985. ...

感谢您的回复，最后还是看贴主 @zha23 在二楼发的这篇 J. Chem. Phys. 1993, 99 (2), 1262–1270 了解了具体做法，主要是文章的式 (11)-(15)，把 Hessian 和试探向量的乘法结果推导出了类似双电子 Fock 矩阵的形式，所以能直接用 HF 的迭代做。最开始没看懂还是对 Davidson 对角化的理解没到位。

常轩豪

L1yuu 发表于 2024-5-10 11:31
楼主和各位老师好，想请教一下楼主提到的 “通过直接方法使用Davidson对角化，可以将时间复杂度降为N^4”是 ...

波函数的稳定性理论要追溯到 Cizek J., Paldus J. 的工作 [J. Chem. Phys., 1967, 47(10), p.3976–3985.], 这篇比pople的好读一点.

核心问题是orbital Hessian的正定性, 也就是看最小特征值的正负.
orbital hessian 是包含了(ia||jb), (ij||ab), 需要类似 Fock build 一样 Direct 计算出来, 然后用标准 Davidson 迭代求解, 每一步迭代矩阵构建复杂度是 N^4, 跟scf一样. Pople后来提出有 real, complex instability 以及triplet instability, 这些从数学和算法角度都是一样的.

你也可以建立整个矩阵存起来直接对角化, 计算和存储复杂度应该至少是 N^6 把...

wzkchem5

L1yuu 发表于 2024-5-10 04:31
楼主和各位老师好，想请教一下楼主提到的 “通过直接方法使用Davidson对角化，可以将时间复杂度降为N^4”是 ...

但凡说看了别人发的文献不理解，必须截图说明哪些地方不理解，这样别人才能针对性地给你讲

L1yuu

楼主和各位老师好，想请教一下楼主提到的 “通过直接方法使用Davidson对角化，可以将时间复杂度降为N^4”是如何实现的呢？看了楼主贴出的三篇文献后我还是不太理解，望各位老师指点！

wzkchem5

zjxitcc 发表于 2022-6-5 05:38
大体系的波函数稳定性分析耗时是超过SCF迭代的，不止是“差不多”。要快的话，可以在SCF和stable这些步骤 ...

这个也和用的是什么软件有关吧，以及和SCF迭代次数有关。一般需要检查稳定性的场合，SCF经常比较难收敛，所以SCF用时本来就比正常的SCF要多，而稳定性分析虽然也需要迭代（Davidson迭代），但是一般不至于收敛很困难，所以迭代次数可能会少些

zjxitcc

zha23 发表于 2022-6-5 11:41
嗯，应该是我画结构时手抖了。我现在主要的困惑是，做一个稳定性分析，其时间复杂度都快和完成一个自洽场 ...

大体系的波函数稳定性分析耗时是超过SCF迭代的，不止是“差不多”。要快的话，可以在SCF和stable这些步骤里全用上RI加速

zha23

zjxitcc 发表于 2022-6-5 10:59
不合理的结构就不用详细讨论了，不同程序算出不一样的结果很正常，因为异常结构一般都有多个SCF解。

嗯，应该是我画结构时手抖了。我现在主要的困惑是，做一个稳定性分析，其时间复杂度都快和完成一个自洽场差不多了，这其中是否有什么技巧，所以想向各位请教一下。

zjxitcc

不合理的结构就不用详细讨论了，不同程序算出不一样的结果很正常，因为异常结构一般都有多个SCF解。

Daniel_Arndt

你应该先看看是不是有两个氢原子相互靠得太近了，这种情况用高斯算的话很容易出现不对头的结果。估计用其他软件处理这种明显不合理的结构也不大可能算出有意义的结果。

wzkchem5

zha23 发表于 2022-6-4 17:17
我打开mol文件看了一下，确实不合理，乱了，但是所有软件都是通过这个mol文件生成的，为什么能量会有 ...

有可能。如果结构不合理的话，本来就不容易收敛到最稳定解，这时一个程序收敛到了最稳定解、另一个没有，不能说后者一定比前者差。

zha23

wzkchem5 发表于 2022-6-5 00:10
可以把坐标也发一下
看看是不是坐标单位搞错了，或者构象不合理

我打开mol文件看了一下，确实不合理，乱了，但是所有软件都是通过这个mol文件生成的，为什么能量会有差异，是因为初始猜测不同导致的吗？

wzkchem5

zha23 发表于 2022-6-4 16:55
这个是测试结果

可以把坐标也发一下
看看是不是坐标单位搞错了，或者构象不合理

zha23

zjxitcc 发表于 2022-6-4 23:34
要不你把这“九个碳的烷烃”坐标和各个程序的能量发上来看看。。。

Created by GaussView 5.0.9
31 29  0  0  0  0  0  0  0  0  0    0
   -1.6367    1.0686    1.0393 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9934    0.5642    0.1656 C   0  0  0  0  0  0  0  0  0  0  0  0
   -3.0634    0.5643    0.1656 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4801   -0.8877    0.1656 C   0  0  0  0  0  0  0  0  0  0  0  0
   -0.4101   -0.8877    0.1656 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8368   -1.3921   -0.7081 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9934   -1.6136    1.4230 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6366   -1.1094    2.2967 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.0634   -1.6135    1.4231 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4804   -3.0657    1.4229 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8382   -3.5702    0.5498 H   0  0  0  0  0  0  0  0  0  0  0  0
   -0.4104   -3.0661    1.4217 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9926   -3.7910    2.6810 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6385   -4.8007    2.6797 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6327   -3.2880    3.5542 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.0626   -3.7884    2.6835 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4353   -4.8553    3.1394 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.0786   -5.8641    3.1394 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.0786   -4.3509    4.0130 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.0786   -4.3509    2.2657 H   0  0  0  0  0  0  0  0  0  0  0  0
   -2.5053   -4.8553    3.1394 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4801    1.2902   -1.0918 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8352    0.7847   -1.9655 H   0  0  0  0  0  0  0  0  0  0  0  0
   -0.4101    1.2919   -1.0908 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.9957    2.7413   -1.0932 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.6410    3.2467   -0.2193 H   0  0  0  0  0  0  0  0  0  0  0  0
   -3.0657    2.7396   -1.0947 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.4817    3.4675   -2.3502 C   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8364    2.9621   -3.2242 H   0  0  0  0  0  0  0  0  0  0  0  0
   -1.8400    4.4757   -2.3512 H   0  0  0  0  0  0  0  0  0  0  0  0
   -0.4117    3.4692   -2.3487 H   0  0  0  0  0  0  0  0  0  0  0  0
  1  2  1  0  0  0  0
  1  3  1  0  0  0  0
  1  4  1  0  0  0  0
  1 22  1  0  0  0  0
  4  5  1  0  0  0  0
  4  6  1  0  0  0  0
  4  7  1  0  0  0  0
  7  8  1  0  0  0  0
  7  9  1  0  0  0  0
  7 10  1  0  0  0  0
10 11  1  0  0  0  0
10 12  1  0  0  0  0
10 13  1  0  0  0  0
13 14  1  0  0  0  0
13 15  1  0  0  0  0
13 16  1  0  0  0  0
17 18  1  0  0  0  0
17 19  1  0  0  0  0
17 20  1  0  0  0  0
17 21  1  0  0  0  0
22 23  1  0  0  0  0
22 24  1  0  0  0  0
22 25  1  0  0  0  0
25 26  1  0  0  0  0
25 27  1  0  0  0  0
25 28  1  0  0  0  0
28 29  1  0  0  0  0
28 30  1  0  0  0  0
28 31  1  0  0  0  0

这个是我计算用的mol文件