|
|
本帖最后由 zjxitcc 于 2025-10-5 11:35 编辑
这个问题比较复杂。可以分为两部分来理解
（1）UCCSD(T)计算
W10这个体系有两个UHF解，在cc-pVTZ基组下电子能量分别为 (单位: a.u.)
| | solution 1 | solution 2 | | UHF | -381.911153 | -381.922008 | | UCCSD(T) | -383.495364 | -383.490805 |
第一个UHF解是不稳定的，第二个是稳定的。你用了前者进行UCCSD(T)计算，这通常是不合理的。有一个经验规律是，不稳定的UHF解通常会导致更大的相关能，所以尽管你的UHF能量偏高，但你的UCCSD(T)能量却偏低（这不是什么好事）。可能你是用Gaussian/ORCA“写一行关键词，提交1个任务获得结果”这种黑箱方式做的计算，以为这样就可以了，没有检验UHF波函数稳定性。对于电子结构较为复杂的体系，即使是单参考计算也要加倍小心，无论使用哪个传统量化软件都要仔细检查计算过程和计算结果。建议先做UHF计算，检验UHF波函数稳定性，确保波函数稳定，然后读取UHF波函数进行UCCSD(T)计算，分成两三步。
而你给出的TS11分子只有一个UHF稳定解，UHF能量是-381.910078 a.u., UCCSD(T)能量是-383.482087 a.u., 所以这一基元步骤的UCCSD(T)电子能垒是5.47 kcal/mol。这一步如果弄错了，后面对比NEVPT2能垒就容易困惑
最后算出的能垒为2.3kcal/mol(NEVPT2/cc-PVTZ，没有加零点能校正)，比CCSD(T)/cc-PVTZ低了6kcal/mol
（2）NEVPT2计算
这一基元步骤所需的最小活性空间为CAS(3,3)，活性轨道包含目标C-C键成键轨道、反键轨道，以及体系自身二重态对应的1个单占据轨道；对应的活性电子也是3个。
上图是W10的3个重要GVB轨道，可以用GaussView/Multiwfn打开文件W10_cc-pVTZ_uhf_uno_asrot2gvb22_s.fch进行可视化。在W10分子中目标C-C键还未断裂，成键轨道占据数为1.99，非常接近2.0，因此本段所说的“重要”不是根据活性轨道占据数来判断的，而是根据基元反应机理判断的。目标C-C键轨道不在32号轨道附近，需要用permute_orb()函数将其调换至31号、33号轨道位置。
理论上讲CAS(3,3)就足够了，可能你会疑惑NEVPT2(3,3)是否足够准确，这的确需要测试。我们可以使用更高阶的多参考方法例如FIC-MRCISD+Q或FIC-NEVPT4(SD)，这两个在ORCA里都有，在活性空间合理的情况下往往可以获得靠谱的结果。若活性空间增大，由于基组cc-pVTZ下空轨道数目较多，FIC-NEVPT4(SD)大概率算不动。这是我算的结果
可以看到，基于CASSCF(3,3)的FIC-NEVPT2和FIC-NEVPT4(SD)能垒与UCCSD(T)能垒很接近。就算FIC-NEVPT2结果是凑巧，还有FIC-NEVPT4(SD)可以说明问题。如果一定要追求更大活性空间，需要小心地增加的活性轨道和活性电子。使用MOKIT的automr对W10和TS11分别进行自动CASSCF/cc-pVTZ计算，自动确定的活性空间分别是(7,7)和(9,9)，显然我们需要增大/减小活性空间使二者一致。
注意，无论看xxx_gvbN_s.fch文件（在这里N=21或22）还是xxx_CASSCF_NO.fch文件，都可以看出活性轨道主要成份，前者是后者的初始轨道。但如果要调换轨道实现活性空间的一致，则需要调换第一个文件，这是因为GVB pair轨道是局域在每一根化学键上的，化学意义十分清晰，而CASSCF自然轨道是离域在整个体系中的，往往作为结果查看，而不作为可调换的过程。GVB是automr做自动CASSCF计算默认的一个中间步骤，它也很重要，不能忽略它直接去看CASSCF自然轨道。待我们看完这两个文件中的轨道
W10_cc-pVTZ_uhf_uno_asrot2gvb22_s.fch
TS11_cc-pVTZ_uhf_uno_asrot2gvb21_s.fch
可以发现乙烯基取代五元环上3根C-C pi键是比较重要的轨道，占据数偏离2.0/0.0最多，同时它们也跟自由基轨道可能共轭
所以优先将这3对轨道加入活性空间。注意每根化学键对应一对轨道（成键和反键轨道）、2个电子，因此在CAS(3,3)基础上扩增为CAS(9,9)，电子能垒变为4.38 kcal/mol（见上面表格数据），还算合理。若还想进一步扩大活性空间，应当优先考虑C-O pi键（因为它在GVB轨道里是偏离2.0/0.0次多的），从而可以构造CAS(11,11)。但这里有个额外的问题，如果进行CASSCF(11,11)轨道优化，W10体系计算会很顺利（能量丝滑下降，下降幅度小，很快收敛），但对TS11会发现CASSCF轨道变化较大，能量下降较多，这是因为在W10的GVB和CASSCF计算中C-O上始终只有一根C-O pi键和一根C-O sigma键，O的其他价电子是俩 孤电子对，在TS11的GVB计算中也是这样；但在TS11的CASSCF计算中有两根C-O pi键、一根C-O sigma键、一对O孤电子，意思就是有一对O孤电子转变成了O->C 配位pi键，成键方式发生较大变化，需要进一步扩大活性空间，将O的孤电子对、C-O sigma等加入活性空间，问题趋于复杂。有一个折中的奇技淫巧是，既然轨道比较合理，我们可以不做轨道优化，直接进行CASCI(11,11)-NEVPT2计算，把O的俩 孤电子对留在双占据空间里，不管它们怎么变化，这就是上述表格中FIC-NEVPT2(11,11) CIonly的来历，能垒为5.56 kcal/mol，也比较合理。
总结：该有机反应在单参考和多参考计算均合理的前提下，二者能垒非常接近（有1.0 kcal/mol左右的区别是正常的，但不会有6 kcal/mol的区别）。
PS1: NEVPT2方法有两个变种SC-NEVPT2和FIC-NEVPT2，后者稍微准确一点点，本回答用的也是FIC-NEVPT2方法。ORCA支持两种，PySCF目前仅支持SC-NEVPT2。
PS2: 如果使用MOKIT以及本文技巧发表文章，请记得引用MOKIT程序，点击此处可以查看引用MOKIT的53篇已发表文献。
计算涉及的输入文件和脚本见附件（没有提供GVB和CASSCF的fch文件以及输出文件，那些文件也比较大，不适合上传论坛。请阅读后自行练习计算，才能较好地理解）
C8H7O_cc-pVTZ_input.7z
(4.52 MB, 下载次数 Times of downloads: 2)
|
评分 Rate
-
查看全部评分 View all ratings
|
楼主 Author
发表于 Post on 2025-10-5 11:29:36
