以圆锥交叉点为起点，优化极小点失败，为什么？

Freeman

大家好。
我在用ORCA 5的新功能，SF-TDDFT，优化激发态。首先我找到了S0/S1的交叉点（虽然自旋污染比较严重，但是勉强能辨认），输入文件如下：

1. ! ci-opt pbe0 d3bj def2-tzvp rijcosx def2-tzvp/c miniprint # 话说def2-tzvp/c能加速tddft吗？输出文件里面好像并没有体现tddft应用了/C辅助基组。
   
2. %maxcore 4000
   
3. %pal
   
4. nproc 40
   
5. end
   
6. %tddft
   
7. sf true # 使用SF-TDDFT。
   
8. nroots 5
   
9. iroot 2 # SF-TDDFT的基态是1号，第一激发态是2号。
   
10. jroot 1
    
11. end
    
12. %geom
    
13. maxstep 0.05 # 经验表明，优化交叉点，步长不能太长，否则容易振荡。
    
14. trust -0.05
    
15. end
    
16. *xyzfile 0 3 e3.orca.xyz # SF-TDDFT 的参考态必须设为高自旋态。我的体系是单重态，所以参考态要设为三重态（应该也可以是五重或更高，没试过）。
    

复制代码

由下图所示，由于组态相互作用，交叉点差不多就是绝热激发态的极小点（右图2的位置）。

于是我以CI为起点，优化激发态的极小点，输入文件如下：

1. ! opt freq pbe0 d3bj def2-tzvp rijcosx def2-tzvp/c miniprint
   
2. %moinp "15ci++.orca.gbw"
   
3. %maxcore 4000
   
4. %pal
   
5. nproc 40
   
6. end
   
7. %tddft
   
8. sf true
   
9. nroots 5
   
10. iroot 2
    
11. end
    
12. %geom
    
13. maxstep 0.05
    
14. trust -0.05
    
15. calc_hess true
    
16. end
    
17. *xyzfile 0 3 15ci++.orca.xyz
    

复制代码

优化时发生了振荡，具体情况是：1态和2态位置总是交换，导致一会儿1态变成优化的目标，一会儿2态又变成优化的目标了。我就很奇怪阿，经过组态相互作用之后，1态和2态不应该变换位置阿。难道说SF-TDDFT用的是透热态的势能面吗？
还有一个疑点：上一个任务交叉点，两态能量十分接近，能量差只有0.000002Eh。如果是绝热态，能量差不该这么小吧！
由此我得出结论，SF-TDDFT没有考虑两个势能面的耦合，导致激发态极小点优化不出来（优化的结构有直接从交叉点沿着上左图的透热势能面滑下来的趋势）。请问大家，我的理解对不对呢？如果对的话，那我该如何优化出激发态极小点呢（如有必要）？

pwzhou

优化激发态最小点从基态构型出发就可以了
从CI点出发优化激发态最小点就是自己给自己找不痛快，没有这么做的

wzkchem5

对，这里SF-TDDFT算出来的是透热态。但是对于这个问题没办法在绝热态势能面上优化，因为导致S0、S1态耦合的是非绝热耦合向量，而一旦在计算里考虑非绝热耦合向量，就破坏了Born-Oppenheimer近似，所以势能面的概念就不成立了。所以此时你可以讨论整个体系的电子-核总波函数，但是不能讨论某个具体核坐标下的S0、S1能量，因此也就不存在S1极小值点的概念了。
你给的图只适用于耦合矩阵元不显含核算符的情形，比如电子转移矩阵元、能量转移矩阵元、SOC矩阵元等。对于耦合矩阵元是非绝热耦合向量的情形，虽然存在绝热态，但是不能画成势能面的形式，因为此时绝热态波函数的核部分和电子部分不是分离的。

Freeman

wzkchem5 发表于 2022-3-2 18:55
对，这里SF-TDDFT算出来的是透热态。但是对于这个问题没办法在绝热态势能面上优化，因为导致S0、S1态耦合的 ...

您说得真详细，非常感谢！
我还有一些疑问。下图是环丁二烯互变反应的势能面，来自https://doi.org/10.1063/1.2222366，是用多参考方法算出来的结果。

看上去，两个1Ag发生了组态相互作用，上面1Ag的左边和下面1Ag的右边应该属于同一个透热态势能面。您的回答中说，势能面的非绝热耦合有时不适合BO近似的框架，因此CI（及其邻近的几何构型）是不能被程序判定为激发态极小点的；可是这个图却明明白白地展示出了耦合后的势能面图，图上的CI（正方形构型）对应着激发态的极小点。
为什么会这样呢？难道这个互变反应涉及电子转移、能量转移或旋归耦合，因而BO近似依然适用吗？或者是源于多参考方法和TDDFT本身的不同吗？

wzkchem5

Freeman 发表于 2022-3-2 14:23
您说得真详细，非常感谢！
我还有一些疑问。下图是环丁二烯互变反应的势能面，来自https://doi.org/10.1 ...

因为绝热态和透热态不是绝对的，而是相对于你考虑的是哪个非对角Hamiltonian矩阵元而言的。你这个帖子给的图是对角化了CAS的Hamiltonian，所以考虑了两个Kekule价键式的共振作用对应的耦合矩阵元，但是没有考虑非绝热耦合向量对应的耦合矩阵元，所以这个图是对于价键共振绝热，对于非绝热耦合向量透热。而你主题帖的体系也是对于价键共振绝热，对于非绝热耦合向量透热。只不过当两个势能面出现避免交叉的时候，你更容易看到绝热的那一面；当两个势能面出现锥形交叉的时候，你更容易看到透热的那一面。
换句话说，环丁二烯这两个价键式是纯粹的透热态，它们之间有两个耦合项（价键共振和非绝热耦合），如果两者都考虑是绝热态，但是你这里只考虑了前者而没有考虑后者，所以不能脱离语境笼统地说它是绝热态还是透热态，只能说对于价键共振绝热，对于非绝热耦合向量透热。

Freeman

wzkchem5 发表于 2022-3-2 23:58
因为绝热态和透热态不是绝对的，而是相对于你考虑的是哪个非对角Hamiltonian矩阵元而言的。你这个帖子给 ...

好深奥啊！看了三遍才有点明白了。我再用自己的语言来表达一下，您看看我理解得对不对。

势能面交叉或避免交叉点处有很多种耦合项，可以分成两大类：一类是非绝热耦合向量（我不明白它为什么是个“向量”，但是不影响全局理解），另一类是其他耦合项（包括旋轨耦合和价键耦合等等）。非绝热耦合是必定存在的，很多时候算反应、画势能面时不考虑它，因为一旦考虑它，BO近似就不能用了，甚至连势能面都不存在了；其他耦合项的存在与否是与体系相关的，例如单重态/三重态交叉点有旋轨耦合，环丁二烯互变反应有价键耦合。
如果不考虑非绝热耦合向量，只考虑其他耦合项，就会发现势能面发生了组态相互作用，出现了避免交叉点；如果连其他耦合项都不考虑，或者体系中根本就不存在其他耦合项，势能面就会交叉，而不发生组态相互作用。我一楼的体系就属于不存在其他耦合项的情况，因为两态能量差几乎为零，有交叉点而非避免交叉。
您看我说得对么？

我还有两个问题：
1、您说的“价键耦合”是不是指电子组态之间的耦合？
即透热态1波函数=组态1，透热态2波函数=组态2；
绝热态1波函数=c1*组态1+c2*组态2，绝热态2波函数=c2*组态1-c1*组态2。
该文献用的是CAS，而非GVB，所以您说“价键耦合”让我懵比了好一会儿。
2、我一楼的体系，可不可以把CI直接当成极小点对待？
我算出来的反应路径是，一个分子被激发到激发态，经过一个过渡态，然后来到CI，发生系间窜越。过渡态虚频的振动方向正好指向CI的几何结构，连弯儿都不带拐的。不过，过渡态一般连着两个极小点，所以我试图在CI附近找极小点，却遇到了一楼的问题。在这种情况下，我可否说，过渡态就是连着CI的，所以就不再优化极小点了。毕竟我从您的解答中明白了，这个体系下，CI附近是不能成功优化出极小点的。

wzkchem5

Freeman 发表于 2022-3-2 17:59
好深奥啊！看了三遍才有点明白了。我再用自己的语言来表达一下，您看看我理解得对不对。

势能面交叉或 ...

差不多，不过你一楼的体系其实也可以解释为虽然存在耦合项，但是在锥形交叉点附近耦合项恰好等于0。如果耦合项处处等于0，你会得到intersection seam而不是conical intersection。至于非绝热耦合向量相关的概念可以参见https://pubs.acs.org/doi/10.1021/acs.accounts.1c00312
后面这两个问题：
1. 对，固然你直接用CAS的话计算就不直接涉及价键的问题了，但是如果你要把你算出来的态看作绝热态，想找相应的透热态的话，你会发现透热态就是类似GVB的价键理论波函数，也就是你图里D2h结构下的那两个共振式。
2. 可以说不需要优化极小点，但是我感觉应该不能说CI是极小点。假如你只需要知道激发态PES上的那个过渡态对应的势垒，一旦越过那个过渡态之后发生的过程你不关心，或者已知一旦越过那个过渡态之后的产物是唯一的，那就找到CI就行了。但如果经过CI以后会生成不止一种产物，那就需要跑非绝热动力学了。

Freeman

wzkchem5 发表于 2022-3-3 04:13
差不多，不过你一楼的体系其实也可以解释为虽然存在耦合项，但是在锥形交叉点附近耦合项恰好等于0。如果 ...

你解答了我长期以来算激发态反应的困惑！非常感谢！

Sui1999

您好，我想问一下自旋设置成3，用SF-TDDFT的时候显示第一个激发能是正常的吗？

wzkchem5

Sui1999 发表于 2024-5-21 05:31
您好，我想问一下自旋设置成3，用SF-TDDFT的时候显示第一个激发能是正常的吗？

先定义什么叫正常

Sui1999

wzkchem5 发表于 2024-5-21 15:06
先定义什么叫正常

不好意思老师，少打字了，是想问第一个激发态的激发能是负数正不正常，在优化的过程中开始是负数，后面变正数了，正数应该才是合理的吧。

wzkchem5

Sui1999 发表于 2024-5-21 08:09
不好意思老师，少打字了，是想问第一个激发态的激发能是负数正不正常，在优化的过程中开始是负数，后面变 ...

如果参考态是T1的话，SF-TDDFT的第一激发能是S0能量减去T1能量。大部分有机分子S0低于T1，所以第一激发能一般是负数，正数反倒是少见的情形（而且此时S0不应该叫S0，而应该叫S1；T1不应该叫T1，而应该叫T0）。
“激发能一般是正数”是自旋守恒TDDFT的结论，不要忽略前提条件而直接套到SF-TDDFT上

Sui1999

wzkchem5 发表于 2024-5-21 16:27
如果参考态是T1的话，SF-TDDFT的第一激发能是S0能量减去T1能量。大部分有机分子S0低于T1，所以第一激发能 ...

好的，谢谢老师。