计算化学公社

标题: ORCA的FIC-MRCC使用体验 [打印本页]

作者
Author: ionexchangeC 时间: 2022-4-23 18:05
标题: ORCA的FIC-MRCC使用体验
本帖最后由 ionexchangeC 于 2022-4-23 18:06 编辑

计算的体系是乙烯扭转的势垒。优化极小点构型和过渡态的级别是PBE0/aug-pcseg-2，优化过渡态时使用UKS做对称破缺计算。
然后选择了以下级别做单点计算计算势垒，基组都是cc-pVTZ：
多参考方法（活性空间取（2，2））：
FIC-MRCC，FIC-MRCI，FIC-NEVPT2（及其RI近似和DLPNO近似），FIC-CASPT2（及其RI近似），SORCI，CASSCF
单参考方法（过渡态用UHF做参考态）：
QCISD(TQ)，CCSD(T)，CCSD，MP2，DLPNO-CCSD(T)，DLPNO-CCSD
DFT（过渡态用UKS）：
B2PLYP，PWPB95，PBE0，M06-2X，wB97M-V，BLYP，PBE，MN15L
除了QCISD(TQ)因为ORCA不支持用Gaussian算之外，所有数据都由ORCA5计算。
结果如下：
(, 下载次数 Times of downloads: 57)
以FIC-MRCC为标准值，各方法误差按大小顺序排如下：
(, 下载次数 Times of downloads: 54) 因为只测试了一个化学反应，所以偶然误差可能大，得出的结论也比较有限。

作者
Author: ionexchangeC 时间: 2022-4-23 18:17
本帖最后由 ionexchangeC 于 2022-4-23 18:36 编辑

大致总结一下：
1.从相关能来看FIC-MRCC能比FIC-MRCI计算出更多的相关能，但是求差计算化学反应能量变化时二者差异比较小。另外必须要说的是，ORCA的FIC-MRCC效率极低，四核i5 10210M计算该体系时，两个分子FIC-MRCI用了不到10分钟，但是FIC-MRCC用了将近5天！
2.CASSCF在本例中表现良好，但这显然是运气好而已。
3.RI或者DLPNO版的NEVPT2（CASPT2）和原版几乎无差异，计算化学反应时可以放心使用，激发能的精度应该也是很好的。
4.虽然乙烯扭转到90度时的静态相关不算是很强（T1diag值小于0.02），但是DLPNO-CCSD和DLPNO-CCSD(T)都差原版相当多，这表明对成键方式怪异的分子，DLPNO-CCSD(T)不是很可靠。极小点附近DLPNO-CCSD(T)和原版CCSD(T)差别则很小。
5.QCISD(TQ)从原理上能比CCSD(T)考虑更多的静态相关，但是本体系中优势比较有限（本来静态相关也不强）。CCSD(T)依旧比较理想。CCSD不好。
6.SORCI算化学反应很差，仅好于UHF和DLPNO-CCSD。算激发能的精度我还没做过测试。
7.泛函方面，算有机体系能量精度高的M06-2X和wB97M-V不负众望。若以FIC-MRCI为高精度计算参考值，双杂化中PWPB95明显好于最老的B2PLYP。纯泛函和PBE0比较差。

作者
Author: sobereva 时间: 2022-4-23 20:24
DLPNO-CCSD(T)可以考虑改用DLPNO-CCSD(T1) tightPNO试试

作者
Author: beefly 时间: 2022-4-24 08:31
Molpro的GeCCo模块已经发布了，支持MR-CCSD和MR-CCSD(T)，可以试试效率如何。从一些用户发表的论文来看，活性轨道数已经突破了2，8个应该能算动

https://github.com/ak-ustutt/GeCCo-public

作者
Author: scf 时间: 2022-4-24 08:45
MRCC里的mrcc也没有2的限制，不过可能不是internal contracted的, 算氧分子就挺慢的-.-

作者
Author: chrinide 时间: 2022-4-24 11:33
本帖最后由 chrinide 于 2022-4-24 12:41 编辑

UHF DLPNO似乎有个非常严重的bug，如果UHF初始猜测收敛的波函数不稳定，你做STABRestartUHFifUnstable后重新收敛的正确波函数似乎不会作为post-scf的参考，后面的post-scf计算结果可能是错的。
DLPNO-CCSD(T):
E(0)                                     ... -77.872255838 （非稳定的UHF波函数？）
E(SL-UMP2 in NEV-PNO space)  ... -0.352418642
Initial E(tot)                         ... -78.224674480
<T|T>                                  ...    0.105550592

同样情况下，如果直接用CCSD(T)，HF参考E(0)则是后面稳定的波函数
E(0)                                     ... -77.998581729
E(MP2)                                  ...    -0.289613518
Initial E(tot)                            ... -78.288195246
<T|T>                                     ...    0.088759827

RHF参考稳定波函数下 DLPNO-CCSD(T)和正则CCSD(T)都是相同E(0) ：-78.064066430。
DLPNO-CCSD(T) ：-78.438584129368  CCSD(T)：-78.438581040600 误差非常小
@wzkchem5 帮忙check下，似乎UHF参考的DLPNO有些不恰当的地方，ORCA Forum上有几个UHF下DLPNO错误比较大的悬而未决的问题.
还是说ORCA UHF DLPNO执行其实并不完美，MRCC的LNO作者之一 Peter Nagypeter给过评论，说目前Open-Shell的Local correlation执行得都不太好，包括ORCA和MOlPRO。MRCC自己Open Shell LNO也内测好久了，一直没有发布，应该还是有一些没有解决的问题。

作者
Author: ionexchangeC 时间: 2022-4-24 12:57
(, 下载次数 Times of downloads: 12)
@chrinide 确实，DLPNO对自旋极化单重态取的参考态是错的

作者
Author: wzkchem5 时间: 2022-4-24 15:11

ionexchangeC 发表于 2022-4-24 05:57
@chrinide 确实，DLPNO对自旋极化单重态取的参考态是错的

这里的问题在于DLPNO-CCSD(T)是基于QRO做的reference来做高自旋ROHF计算，而不是基于UHF轨道。所以对于存在反铁磁耦合的体系，目前的DLPNO-CCSD(T)代码不适用。
我看看要不要在orca说明书里或者在输出文件里加一个warning，提醒用户注意这一点

作者
Author: beefly 时间: 2022-4-24 18:17

scf 发表于 2022-4-24 08:45
MRCC里的mrcc也没有2的限制，不过可能不是internal contracted的, 算氧分子就挺慢的-.-

MRCC算4个活性轨道差不多就到头了，再大就很费劲了

作者
Author: mfdsrax2 时间: 2022-4-24 18:42
居然有高斯支持，ORCA不支持的方法

作者
Author: 喵星大佬 时间: 2022-4-24 21:33

mfdsrax2 发表于 2022-4-24 18:42
居然有高斯支持，ORCA不支持的方法

QCI系列方法就是瓜丝的人搞的，其他程序几乎都不支持

作者
Author: zjxitcc 时间: 2022-4-24 21:42

mfdsrax2 发表于 2022-4-24 18:42
居然有高斯支持，ORCA不支持的方法

自CC在理论化学家里兴起后，QCI系列就渐渐被抛弃了，处于比较尴尬的地位，所以时间更晚的量化程序估计基本不考虑实现QCI系列了。QCISD跟CCSD一样的标度，都是N^6，都有大小一致性，但后者理论上更准，所以都奔着CCSD去了。

作者
Author: sobereva 时间: 2022-4-25 02:21
doi: 10.1002/wcms.1131
这篇文章算是给QCI的无用盖棺定论了，坚持用QCI的都可以用这篇来怼

作者
Author: biogon 时间: 2022-4-25 09:09

beefly 发表于 2022-4-24 18:17
MRCC算4个活性轨道差不多就到头了，再大就很费劲了

mrcc为什么比mrci难算的多

作者
Author: beefly 时间: 2022-4-25 11:56

ionexchangeC 发表于 2022-4-23 18:17
大致总结一下：
1.从相关能来看FIC-MRCC能比FIC-MRCI计算出更多的相关能，但是求差计算化学反应能量变化时 ...

Gaussian过去还有CCSD(TQ)，但是不稳定，可能从09版开始就从手册里删除了。目前只有gamess和停止维护的ACES2里还有这个方法。
无论QCISD(TQ)还是CCSD(TQ)，都只是对(T)的少许改进，和T(Q)还差得远。由于是微扰理论，(TQ)的总能量在很多情况下比(T)算的还高，不如T(Q)可靠。

作者
Author: beefly 时间: 2022-4-25 11:58

biogon 发表于 2022-4-25 09:09
mrcc为什么比mrci难算的多

ci解线性方程组，cc解非线性方程组，难度更大

作者
Author: ionexchangeC 时间: 2022-4-25 12:00

beefly 发表于 2022-4-25 11:56
Gaussian过去还有CCSD(TQ)，但是不稳定，可能从09版开始就从手册里删除了。目前只有gamess和停止维护的AC ...

之后的版本里CCSD(TQ)输进去不会提示关键词语法错误，但后面会报错

作者
Author: zjxitcc 时间: 2022-4-25 12:05
本帖最后由 zjxitcc 于 2022-4-25 12:10 编辑

biogon 发表于 2022-4-25 09:09
mrcc为什么比mrci难算的多

两方面：（1）MRCC公式和代码太难了，都要用程序自动推导、自动生成代码、自动构造中间数组，还要并行，优化细节很多，随便写的话，不仅算不大，还容易爆内存；（2）MRCC代码量一般巨大，计算量也巨大。

作者
Author: Freeman 时间: 2022-4-28 14:51
本帖最后由 Freeman 于 2022-4-28 18:51 编辑

chrinide 发表于 2022-4-24 11:33
UHF DLPNO似乎有个非常严重的bug，如果UHF初始猜测收敛的波函数不稳定，你做STABRestartUHFifUnstable后重 ...

我遇到过相关问题。解决方法是，写moread读取前一步稳定波函数时，还要加上uhf关键词，否则会默认用rhf。
而且这跟dlpno无关。就算用非dlpno的post-hf，也有这个问题，也大概是同样的解决方法。
---------------------
看了你的输出文件，原来你已经加过了。不好意思忘记看了。
这才发现我这儿也有同样的问题：以不稳定波函数作初猜，scf部分的能量和后面dlpno-ccsd(t)的参考能量是不同的。

(base) [yzhangnn@login-0 light++]$ grep 'Total Energy' 15ts___.orca.out
Total Energy : -1613.21317953 Eh -43897.76233 eV
(base) [yzhangnn@login-0 light++]$ grep 'E(0)' 15ts___.orca.out
E(0) ... -1613.210319341

复制代码

这里两个能量差得不太大，我猜可能是跟自旋污染有一定联系。这个体系<S^2>=0.4，而你输出文件显示扭曲的乙烯<S^2>=1.0。可能自旋污染越大，uhf偏离rhf越明显，两能量差得越大？

作者
Author: zjxitcc 时间: 2022-9-21 10:52
本帖最后由 zjxitcc 于 2022-9-21 10:55 编辑

这几天我了解了一下，目前ORCA 5里的FIC-MRCCSD只考虑了BCH展开式中两个对易式。另外，组态系数和振幅之间没有做自洽迭代，相当于真正FIC-MRCCSD的第一圈。而MkMRCCSD只能算小活性空间。我自己开发的GVB-BCCC3b虽然活性空间中等，但目前只能考虑活性轨道，还没法算空轨道，而且精度我也不是很满意。早年很多MRCC代码没有开放，估计大多失传了，想做还得从零开始实现。哎，MRCC这个领域仍需要巨量的投入，不知多少年以后才能做多核过渡金属

作者
Author: 胡说 时间: 2024-5-18 14:23
在UHF DLPNO-CCSD(T) 计算中遇到了 6楼和19楼同样的问题，UHF收敛到了稳定波函数，但是post-HF阶段的参考态能量和SCF阶段不一样，用的是更高能的不稳定的SCF波函数，导致得到的DLPNO-CCSD(T)能量是明显不对的。
两年过去了，对于这个问题，是不是有什么解决方法？

作者
Author: zjxitcc 时间: 2024-5-18 16:08
本帖最后由 zjxitcc 于 2024-5-18 16:14 编辑

胡说发表于 2024-5-18 14:23
在UHF DLPNO-CCSD(T) 计算中遇到了 6楼和19楼同样的问题，UHF收敛到了稳定波函数，但是post-HF阶段的参考态 ...

一直以来，这不是个问题，也没有“还没解决”之说。高等级方法的使用需要了解一些原理，才能做出合理的计算。

ORCA中的DLPNO方法对于UHF计算，采用的参考波函数不是UHF，而是QRO单行列式。所谓QRO就是UHF自然轨道（简称UNO），但扔掉UNO轨道占据数、不用这些分数占据数。所以UHF-based DLPNO-CCSD的参考态能量注定不可能等于UHF能量，而且一般高于UHF能量，这不是用关键词调节调节就能“解决”的。只有极少数的体系，UHF与ROHF能量基本相等，才有可能出现UHF-based DLPNO-CCSD的参考态能量十分接近UHF能量。

如果你的体系是开壳层单重态，不应该使用目前的UHF-based DLPNO-CCSD方法，由于它扔掉了参考波函数的多组态特征（或者说轨道的分数占据），理论上就不对了。注意，这不是把UHF轨道换成DFT轨道就可以解决的。如果你的体系是开壳层高自旋，那么QRO就是UNO轨道构成的ROHF行列式，占据数取作{2,1,0}。这时候方法是否work，取决于这个ROHF行列式描述这个体系是否定性正确。若定性正确，则DLPNO-CC合理；若体系存在两个简并的最低二重态，则不能用此方法。

作者
Author: 胡说 时间: 2024-5-18 16:38

zjxitcc 发表于 2024-5-18 16:08
一直以来，这不是个问题，也没有“还没解决”之说。高等级方法的使用需要了解一些原理，才能做出合理的计算 ...

谢谢老师的回复。

我算的体系是开壳层二重态，那么如何判断ROHF对于该体系是否定性正确呢？
计算发现，UHF稳定波函数下能量比UHF-based DLPNO-CCSD中参考态能量低了0.13 eV左右，那么这么大的差距是否可以认为ROHF描述这个体系可能定性不正确？

另外，如果ROKS可以定性正确描述该体系，那是不是可以基于UKS来计算DLPNO-CCSD(T)？

作者
Author: wjc404 时间: 2024-5-20 12:28
本帖最后由 wjc404 于 2024-5-20 13:11 编辑

胡说发表于 2024-5-18 16:38
谢谢老师的回复。

我算的体系是开壳层二重态，那么如何判断ROHF对于该体系是否定性正确呢？

先在UHF级别下找到稳定波函数并生成UNO。如果<S**2>明显高于理论值（二重态的体系是0.75）然后UNO的占据数偏离ROHF的偏差较大（二重态的ROHF的轨道占据应该是一个1.0而其他都是2.0或0.0）那就说明你的体系可能有明显multireference特征。如果这两者都正常，可以放心做DLPNO-CCSD(T)计算。

如果UNO占据数偏离较大，可以把占据数在0.1-1.9之间的轨道纳入活性空间做CASSCF计算，查看收敛后的活性空间轨道没问题后再做DLPNO-NEVPT2计算。

当然有不少体系是有弱的multireference特性但还不至于到非用multireference方法的地步，比如一些单核过渡金属配合物，UHF会高估它们的自旋极化，但后续用Brueckner-CCD一算就发现轨道优化后的<S**2>恢复正常且波函数双激发不明显。此时用ROHF的轨道作为reference，如果CCSD(T)出来的t1diag和最大双激发幅度都正常，就可以用。当然如果你的体系可能属于此类但Brueckner-CCD又算不动，可以把包含体系结构的输入文件以zip附件形式传上来让别人帮忙看看。

作者
Author: 胡说 时间: 2024-5-23 10:53

wjc404 发表于 2024-5-20 12:28
先在UHF级别下找到稳定波函数并生成UNO。如果明显高于理论值（二重态的体系是0.75）然后UNO的占据数偏离R ...

感谢回复。
我检查了下UHF稳定波函数，<S**2>为0.76，与理想值0.75差别不大。且基于UHF生成的UNO的占据数也和ROHF差别很小，都大于1.99或者小于0.01。所以从这一点上看，UHF稳定波函数是合理，且多参考态特征不强。

但ORCA作开壳层的UHF DLPNO-CCSD(T)是先把UHF轨道转成QRO轨道，再以QRO轨道组成的ROHF reference为基础做的，那么问题是如何来判断QRO轨道组成的ROHF reference是对于某体系来说是合理的呢？
比如这里QRO轨道能量顺序是正常的, 但SOMO处显示的alpha和beta轨道能量差别很大（alpha= -14.572 beta= -3.457）, 不知道这个差异大说明什么？

欢迎光临计算化学公社 (http://bbs.keinsci.com/)