计算化学公社

标题: 求助Gaussian中使用B3LYP与UHF得到的自旋密度差异 [打印本页]

作者
Author: qiqi7 时间: 2022-11-15 17:03
标题: 求助Gaussian中使用B3LYP与UHF得到的自旋密度差异
各位老师好！我在进行自由基反应相关计算时，首先使用Gaussian，级别为UB3LYP-D3(bj)/6-311+G(d,p)优化中间的一些结构（单重态,C14H10），并使用stable=opt检验波函数已经稳定，且<S**2>为0是闭壳层体系，符合本人预期。接下来我使用Gaussian读取上一步的优化结构，级别为UHF/cc-pVQZ去得到稳定的波函数，预备作为orca高精度单点能的初猜。两个计算结束后，我将两个fch文件读入Multiwfn检查自旋密度图，此时发现两个级别下的自旋密度图存在的不小差异，请问造成这种事情的原因可能有哪些，该如何解决？

(, 下载次数 Times of downloads: 2)
(, 下载次数 Times of downloads: 3)

作者
Author: wzkchem5 时间: 2022-11-15 17:19
HF因为没有相关项，会高估波函数发生对称性破缺的概率（倾向于让波函数发生自旋对称性破缺降低能量，来弥补相关能的缺失）。所以当UHF自旋破缺而UDFT不破缺时，一般后者更能代表真实的电子结构

作者
Author: qiqi7 时间: 2022-11-15 17:25

wzkchem5 发表于 2022-11-15 17:19
HF因为没有相关项，会高估波函数发生对称性破缺的概率（倾向于让波函数发生自旋对称性破缺降低能量，来弥补 ...

感谢老师解答！请问这种情况是否就可以理解为，该结构是闭壳层态的，当我要使用gaussian计算HF波函数作为orca的初猜时，应该使用RHF波函数？

作者
Author: qiqi7 时间: 2022-11-15 17:26

wzkchem5 发表于 2022-11-15 17:19
HF因为没有相关项，会高估波函数发生对称性破缺的概率（倾向于让波函数发生自旋对称性破缺降低能量，来弥补 ...

另外，请问老师有相关文献推荐吗？想多了解一些

作者
Author: wzkchem5 时间: 2022-11-15 20:09

qiqi7 发表于 2022-11-15 10:25
感谢老师解答！请问这种情况是否就可以理解为，该结构是闭壳层态的，当我要使用gaussian计算HF波函数作为 ...

取决于你用orca做什么计算。
如果是DFT（包括双杂化泛函），那orca用什么泛函高斯就用什么泛函，如果高斯没有这个泛函就用相近泛函代替；
如果是CCSD(T)，那么用UHF一般稍好一些，因为UCCSD(T)满足size consistency而RCCSD(T)不满足。但假如你不研究键解离等问题，那RHF也未必错误。
如果是DLPNO-CCSD(T)，个人建议用RHF，因为DLPNO-CCSD(T)暂时不支持以broken symmetry态为参考态，虽然程序也给你算，但是结果和RCCSD(T)相似。所以不如一开始就用RHF参考态。

另外说实话，这种分子建议直接用orca算。绝对不要把“用高斯收敛以后作为orca的初猜”作为默认方法，orca的SCF收敛并不差，尤其从orca 5.0开始更是如此。即使orca 5.0收敛的概率真的仍然不如高斯高（目前没见到有人测过），你先用高斯再用orca的话，等于是一定要跑至少两个计算（还要中间转化一次文件格式），直接用orca倒有很大概率是一次计算可以收敛，其实平均上讲用高斯准备初猜很可能反倒浪费时间。

作者
Author: qiqi7 时间: 2022-11-15 20:55

wzkchem5 发表于 2022-11-15 20:09
取决于你用orca做什么计算。
如果是DFT（包括双杂化泛函），那orca用什么泛函高斯就用什么泛函，如果高 ...

感谢老师的解答。因为结构较多，我是想用DLPNO-CCSD(T)计算，然后计算单点能过程中，大多数分子能量合理，对于UB3LYP下<S**2>值大于0的体系出现了明显不合理的结果，所以才尝试使用gaussian进行对称破缺态的UHF计算，然后再使用orca的DLPNO-CCSD(T) 结合UHF关键词计算DLPNO-UCCSD(T)。然后顺手选了几个闭壳层分子想测试一下，才发现了上述现象。个人接下来可能放弃DLPNO计算，直接使用CCSD(T)结合小一点的基组计算该体系。感谢老师，学到了很多平时会忽略的地方。

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

qiqi7 发表于 2022-11-15 13:55
感谢老师的解答。因为结构较多，我是想用DLPNO-CCSD(T)计算，然后计算单点能过程中，大多数分子能量合理 ...

其实可以尝试一下用orca的NEVPT2算一下，那个苯环的轨道可以不纳入活性空间，只把剩下8个碳原子的垂直于平面的p轨道纳入活性空间。这个计算可以用RI加速，虽然也可以用DLPNO，但这个大小的体系DLPNO-NEVPT2估计不比RI-NEVPT2快。有条件的话可以再加F12校正，来解决基组误差问题。
这个计算结果应该比UCCSD(T)靠谱，毕竟是多参考态方法，既size consistent又没有自旋污染。

作者
Author: qiqi7 时间: 2022-11-16 09:32

wzkchem5 发表于 2022-11-15 21:24
其实可以尝试一下用orca的NEVPT2算一下，那个苯环的轨道可以不纳入活性空间，只把剩下8个碳原子的垂直于 ...

好的，感谢老师，我去尝试一下

作者
Author: wjc404 时间: 2022-11-18 10:59
本帖最后由 wjc404 于 2023-2-16 18:43 编辑

同意二楼的看法，这个分子的对称破缺被HF高估了。
尝试了以下计算过程：
（1）用m06-2X/def2-TZVP优化几何结构，得到平面无扭曲的结构，无虚频，波函数稳定（无自旋极化）。接下来的步骤在该几何结构下进行，基组不变。
（2）用HF得到稳定的波函数，为自旋极化波函数，S2 = 1.03。用stable=opt关键词，或先设置两个片段（茚环和茂环各为一片段）生成片段组合波函数然后优化，结果一样。
（3）用Brückner-CCD优化上一步UHF所得的轨道，发现过程中参考波函数<S**2>逐渐降至0.01以下，最后的最大双激发幅度为0.053。
（4）上一步优化的轨道投影到spin-restricted模式，然后Brückner-CCD计算出密度矩阵，从自然轨道看，1以上的最低三个占据数为1.91/1.91/1.90，1以下的最高三个占据数为0.09/0.07/0.06。
以上结果应该能说明，这个分子在最优结构下的静态相关不强，HF无法合适描述它的可能原因是缺乏动态相关的考虑。HF的这一问题在大共轭体系中不罕见。
计算结果见https://pan.baidu.com/s/1S88wjWKVaF7fI-H7SJXVLw，提取码44ik。

作者
Author: qiqi7 时间: 2022-11-20 10:18

wjc404 发表于 2022-11-18 10:59
同意二楼的看法，这个分子的对称破缺被HF高估了。
尝试了以下计算过程：
（1）用m06-2X/def2-TZVP优化几 ...

感谢大佬，给我以后遇到量化问题时自己分析问题提供了一些思路

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