计算化学公社
标题: 求助ORCA DLPNO-CCSD(t)-F12计算能量异常问题 [打印本页]
作者Author: pet 时间: 2025-8-16 03:18
标题: 求助ORCA DLPNO-CCSD(t)-F12计算能量异常问题
本帖最后由 pet 于 2025-8-16 03:18 编辑
各位老师好，我在用 ORCA 计算 C₆H₁₂O₃ 阳离子体系的单点能时遇到了一些问题，请教一下各位老师。（实际上是为了计算同一分子不同构象的电离能，目前中性分子没有问题，所以这里简化为阳离子能量计算的问题）
流程如下：
先在 Gaussian 16 中用 M06-2X-D3/def2-TZVP 优化了不同构象。
再用 ORCA 计算单点能，关键词为：! DLPNO-CCSD(T)-F12 cc-pVDZ-F12 TightPNO cc-pVDZ-F12-CABS cc-pVTZ/C tightSCF SCFConvForced noautostart miniprint，其中一个构象（06i）的能量值明显异常。波函数稳定性测试显示The stability analysis shows that the wavefunction is stable。
图 中 给出了几个构象在 DFT 和 DLPNO-CCSD(T)-F12/cc-pVDZ-F12 下的能量对比，其中 06i 为异常结构。
另外，由于最近发现 M06-2X 优化结果有时不太理想，又用 B3LYP-D3BJ/6-311+G(d,p) 对 06i 重新优化，并计算单点（06i-b）。这两个结构确实有差异，DFT结果相差2.5kcal/mol，但是按说不至于导致DLPNO-CCSD(T)-F12 结果相差这么多吧。
附件中分别是m062x（06i-stable.out）和b3lyp（06i-b.out）优化完的结构在orca中计算单点的输出文件。请老师们帮忙分析一下可能的原因。
(, 下载次数 Times of downloads: 23)
(, 下载次数 Times of downloads: 4)
(, 下载次数 Times of downloads: 1)
作者Author: sobereva 时间: 2025-8-16 05:20
多一事（DLPNO）不如少一事，是我的话就用CCSD(T)/cc-pVTZ，结合MP2估计的TZ-QZ基组校正量了，这样大小的体系算起来很轻松
我不推荐M06-2X优化。但凡B3LYP-D3(BJ)优化不理想的情况（下文），推荐优先考虑wB97XD（或者ORCA里的wB97X-D3，或组合方法wB97X-3c）
谈谈量子化学研究中什么时候用B3LYP泛函优化几何结构是适当的
http://sobereva.com/557（http://bbs.keinsci.com/thread-17899-1-1.html）
作者Author: zjxitcc 时间: 2025-8-16 18:44
试试看加上
- %basis
- PCDTrimBas Overlap
- PCDTrimAuxJ Coulomb
- PCDTrimAuxJK Coulomb
- PCDTrimAuxC Coulomb
- PCDThresh -1
- end
作者Author: pet 时间: 2025-8-18 16:38
谢谢sob老师的回答，关于CCSD(T)/cc-pVTZ结合MP2估计的TZ-QZ基组校正量这种外推方式我也有点疑问，我在ORCA中测试C3H6O3，关键词
! ExtrapolateEP2(3/4,cc,MP2) tightSCF noautostart miniprint
%maxcore 5000
%pal nprocs 16 end
竟然算了三个多小时，相同的体系我用CCSD(T)/jul-cc-pV(T+d)Z 也才算了40min，总觉得是不是哪里有问题，，，这个耗时我根本不敢去算C6H12O3阳离子的体系。
(, 下载次数 Times of downloads: 0)
作者Author: pet 时间: 2025-8-18 16:39
谢谢老师的回复，结果还是一样.....
作者Author: zjxitcc 时间: 2025-8-18 16:42
是否可以上传一下06i-stable.out任务 在加了3L建议的关键词后的计算输出文件？
作者Author: pet 时间: 2025-8-18 17:01
好的，已传
(, 下载次数 Times of downloads: 2)
作者Author: zjxitcc 时间: 2025-8-18 20:00
本帖最后由 zjxitcc 于 2025-8-18 20:08 编辑
这个体系从中性单重态 失去1个电子变成 阳离子二重态时，最可能失去1个电子的地方有两处：（1）C=O双键中那个O的孤电子对；（2）O-O过氧键中某一个O的孤电子对。因此当前二重态体系大概率会有2个UHF稳定解，代表两种物理含义：
（1）C=O双键中的O失去了1个电子，剩余1个单电子仍然在这个O附近。
（2）O-O过氧键中某个O失去了1个电子，剩余1个单电子仍然在这附近。
你上传的06i-stable.out和06i-3.out文件都是收敛到了情况（2），它的UHF能量远高于（1）的UHF能量，也就是说通常情况下我们应该取UHF解（1）进行后续计算。无论你用UHF DLPNO-CCSD(T)，还是UHF UCCSD(T)，由于它们都是从UHF出发，你都可能会碰到CCSD(T)能量异常，碰到问题的概率约等于50%。
我们先构造和使用4种SCF初猜，分别进行UHF/cc-pVDZ-F12计算，4个gjf文件如下
(, 下载次数 Times of downloads: 9)
,
(, 下载次数 Times of downloads: 9)
,
(, 下载次数 Times of downloads: 8)
,
(, 下载次数 Times of downloads: 7)
为方便阅读，展示第一个gjf文件内容
- %chk=06i_uhf_gau.chk
- %mem=160GB
- %nprocshared=64
- #p UHF gen nosymm int=nobasistransform scf(xqc,maxcycle=500) stable=opt
- UHF/cc-pVDZ-F12 sp calculation
- 1 2
- C 2.20560700 3.00976100 -1.10959300
- ...
- 06i_uhf_gau.log: SCF Done: E(UHF) = -458.546387235 a.u. after 18 cycles
- 06i_frag1_uhf.log: SCF Done: E(UHF) = -458.482233098 A.U. after 31 cycles
- 06i_frag2_uhf.log: SCF Done: E(UHF) = -458.482233098 A.U. after 32 cycles
- 06i_frag3_uhf.log: SCF Done: E(UHF) = -458.546387234 A.U. after 23 cycles
- formchk 06i_uhf_gau.chk 06i_uhf_gau.fch
- fch2mkl 06i_uhf_gau.fch
- orca_2mkl 06i_uhf_gau_o -gbw
- %pal nprocs 32 end
- %maxcore 5000
- ! UHF DLPNO-CCSD(T)-F12 TightPNO cc-pVDZ-F12-CABS cc-pVTZ/C TightSCF
- ----------------------------------------D-I-I-S--------------------------------------------
- Iteration Energy (Eh) Delta-E RMSDP MaxDP DIISErr Damp Time(sec)
- -------------------------------------------------------------------------------------------
- *** Starting incremental Fock matrix formation ***
- *** Initializing SOSCF ***
- ---------------------------------------S-O-S-C-F--------------------------------------
- Iteration Energy (Eh) Delta-E RMSDP MaxDP MaxGrad Time(sec)
- --------------------------------------------------------------------------------------
- 1 -458.5463872330064987 0.00e+00 8.74e-09 7.14e-08 6.49e-09 3.1
- *** Restarting incremental Fock matrix formation ***
- 2 -458.5463872330067261 -2.27e-13 5.89e-09 3.35e-08 3.19e-08 2.9
- *** Gradient check signals convergence ***
- *****************************************************
- * SUCCESS *
- * SCF CONVERGED AFTER 2 CYCLES *
- *****************************************************
PS1. 若想查看两个UHF解的单电子分别在哪里，可运行以下Python脚本
- from mokit.lib.gaussian import uno
- uno('06i_uhf_gau.fch')
PS2. 虽然该二重态体系在UHF下有两个解，但它在UM062X下可能有两个解、也可能只有一个解，结论不可简单照搬。原则上你之前算的01i~05i也应当做这种检查。由于此体系电子结构略为复杂，我不建议你采用基组外推的方式进行计算，那会包含多个SCF步骤然后把自己绕晕进去。
PS3. 上述技巧多次用到了MOKIT，若您使用这些技巧并发表文章，除引用ORCA外也建议引用MOKIT。引用MOKIT的50篇已发表文章可以看https://jeanwsr.gitlab.io/mokit-doc-mdbook/citing.html
作者Author: pet 时间: 2025-8-18 21:47
非常感谢老师您细致的回复！ 一定认真学习并正确引用~
我能够理解整个流程，但是关于您的gjf文件还有两个疑问：
1. 06i_frag1_uhf.gjf 和 06i_frag2_uhf.gjf 应该是对应您提到的第二种失电子情况：从O-O过氧键中某一个O的孤电子对上失去电子，但是为什么gjf的片段2里只包含OH呢？
2. 我尝试把OOH包含进frag2，设置自旋多重度为1 2 1 1 0 2以及1 2 0 2 1 1 , 我理解的第二种情况对应OOH上失电子，第一种情况对应frag1失电子，但是这两种算出来都是收敛到了高能量的情况（结果上看与您的 06i_frag1_uhf.gjf 和 06i_frag2_uhf.gjf完全一致），，，这是为什么呢？为什么第一种已经排除了OOH失电子，也还是不能收敛到羰基O失电子的情况去呢？
作者Author: zjxitcc 时间: 2025-8-19 15:50
本帖最后由 zjxitcc 于 2025-8-19 17:24 编辑
存在06i_frag1_uhf.gjf和06i_frag2_uhf.gjf这两个文件，是因为当时我对“O-O过氧键中某个O失去了1个电子”没有100%的把握，我怀疑它可以进一步细分为12号O失去1个电子、13号O失去1个电子两种不同情况，所以片段2只含OH两个原子。等我算完才发现这里没有两种细分，实际上是12号O失去1个电子，说成“O-O过氧键中某个O失去了1个电子”也没问题。所以这两个gjf文件对应的都是我在8L回答中说的情况（2）。
06i_frag3_uhf.gjf对应的是我在8L回答中说的情况（1）。而06i_uhf_gau.gjf这个文件，就是想看看Gaussian默认初猜能收敛出啥来，一个常规计算，没想到它直接收敛到了较低能量的UHF解。你说把OOH划分成一个片段当然也可以，也对应情况（2）。只有像我给出的06i_frag3_uhf那样划分才对应情况（1）。
作者Author: pet 时间: 2025-8-20 20:17
好的，非常感谢您
作者Author: wjc404 时间: 2025-8-25 21:20
本帖最后由 wjc404 于 2025-8-27 20:53 编辑
这种有多个UHF稳定解的电离态，比较靠谱的诊断以及判断方式是不是IP-EOM-CCSD? 假设没遇到异常能量差的情况，可能都不知道没有找到能量最低的解。。
Update: 用ORCA算了一下IP-EOM-CCSD发现最低的电离态对应羰基的孤对电子的电离，其次才是过氧键那边的氧的电子。
(, 下载次数 Times of downloads: 11)
| 欢迎光临 计算化学公社 (http://bbs.keinsci.com/) |
Powered by Discuz! X3.3 |