优化T1能量时发现能量值比S0态能量低

千面追风

如下面输出文件所示，
Excited State   1:      Triplet-A      -1.1324 eV    -1094.85 nm  f=-0.0000  <S**2>=2.000
     324 -> 329        0.28616
     325 -> 329       -0.55964
     325 -> 330       -0.16628
     326 -> 329       -0.23879
     325 <- 329       -0.13272
This state for optimization and/or second-order correction.
Total Energy, E(TD-HF/TD-KS) =  -3911.75446089   
Copying the excited state density for this state as the 1-particle RhoCI density.

Excited State   2:      Triplet-A      1.8447 eV  672.11 nm  f=0.0000  <S**2>=2.000
     324 -> 329        0.57285
     324 -> 330        0.16969
     325 -> 329        0.28765

Excited State   3:      Triplet-A      2.4214 eV  512.04 nm  f=0.0000  <S**2>=2.000
     325 -> 329        0.21559
     326 -> 329       -0.62582
     326 -> 330       -0.15334

当我在使用0 1 td(triplet)方法计算T1态时，能量相较于基态呈现负增长，想寻求各位老师如何解决这一问题。
附：使用0 3方法也出现负增长
附录已上传，包含基态优化结果，S1T1优化结果，带有（1）的是使用如上01方法计算T1态的输出文件。

（没加极化的原因是时间来不及了）

zjxitcc

在论坛上搜索“负激发能”、“波函数不稳定”两个词，论坛上对此问题有很多帖子（及回复）。

千面追风

zjxitcc 发表于 2024-7-25 00:26
在论坛上搜索“负激发能”、“波函数不稳定”两个词，论坛上对此问题有很多帖子（及回复）。

老师您好：在听从论坛各大帖子回复的意见后我进行了单点能计算并添加了stable=opt这一关键词，后用guess=read读取初猜进行了T1激发态的计算，但是中途计算时出现了如下结果
Excited State   1:  3.000-A      0.0991 eV 12510.42 nm  f=0.0000  <S**2>=2.000
    298A -> 361A       0.12816
    307A -> 353A       0.13572
    309A -> 330A      -0.10088
    315A -> 330A      -0.19356
    318A -> 330A      -0.14837
    324A -> 330A       1.34600
    324A -> 331A       0.13018
    324A -> 333A      -0.36753
    324A -> 346A       0.17770
    325A -> 330A      -0.95169
    325A -> 333A      -0.26763
    325A -> 353A      -0.11548
    326A -> 330A      -0.46092
    298B -> 361B      -0.12816
    307B -> 353B      -0.13572
    309B -> 330B       0.10087
    315B -> 330B       0.19357
    318B -> 330B       0.14837
    324B -> 330B      -1.34604
    324B -> 331B      -0.13015
    324B -> 333B       0.36753
    324B -> 346B      -0.17770
    325B -> 330B       0.95164
    325B -> 333B       0.26766
    325B -> 353B       0.11548
    326B -> 330B       0.46091
    298A <- 361A       0.12506
    307A <- 353A       0.13153
    315A <- 330A      -0.18135
    318A <- 330A      -0.13673
    324A <- 330A       1.24104
    324A <- 331A       0.12296
    324A <- 333A      -0.34805
    324A <- 346A       0.16725
    325A <- 330A      -0.87771
    325A <- 333A      -0.25516
    325A <- 353A      -0.11071
    326A <- 330A      -0.41966
    298B <- 361B      -0.12506
    307B <- 353B      -0.13153
    315B <- 330B       0.18135
    318B <- 330B       0.13673
    324B <- 330B      -1.24108
    324B <- 331B      -0.12293
    324B <- 333B       0.34804
    324B <- 346B      -0.16725
    325B <- 330B       0.87766
    325B <- 333B       0.25518
    325B <- 353B       0.11071
    326B <- 330B       0.41966
This state for optimization and/or second-order correction.
Total Energy, E(TD-HF/TD-KS) =  -3911.71051490   
Copying the excited state density for this state as the 1-particle RhoCI density.

Excited State   2:  3.000-A      2.3129 eV  536.05 nm  f=0.0000  <S**2>=2.000
    307A -> 330A      -0.10322
    324A -> 330A       0.37547
    324A -> 333A      -0.13873
    325A -> 330A       0.54252
    307B -> 330B       0.10322
    324B -> 330B      -0.37547
    324B -> 333B       0.13873
    325B -> 330B      -0.54256

Excited State   3:  3.000-A      2.5350 eV  489.09 nm  f=0.0000  <S**2>=2.000
    327A -> 329A       0.18829
    328A -> 329A      -0.67092
    327B -> 329B      -0.18829
    328B -> 329B       0.67092
这个激发光波长都已经到远红外线那一区域了，这可能是什么原因造成的呢？
文件已上传

zjxitcc

千面追风 发表于 2024-7-25 13:46
老师您好：在听从论坛各大帖子回复的意见后我进行了单点能计算并添加了stable=opt这一关键词，后用guess= ...

（1）如果有结构优化，新结构下可能基态波函数不稳定。
（2）没有涉及结构优化，但错误地读取了 算好的稳定波函数的轨道。
（3）ghost态，见这两篇博文
《Multiwfn已支持Ghost态诊断指数，以及CM5电荷》http://bbs.keinsci.com/thread-6707-1-1.html
《Multiwfn支持的电子激发分析方法一览》http://bbs.keinsci.com/thread-11406-1-1.html
逐一排查即可。

千面追风

zjxitcc 发表于 2024-7-25 16:01
（1）如果有结构优化，新结构下可能基态波函数不稳定。
（2）没有涉及结构优化，但错误地读取了 算好的 ...

谢谢老师答复！

千面追风

zjxitcc 发表于 2024-7-25 16:01
（1）如果有结构优化，新结构下可能基态波函数不稳定。
（2）没有涉及结构优化，但错误地读取了 算好的 ...

如果是波函数不稳定的话，我能做的解决方案是否是将tddft计算的结果作为初猜反复进行stable=opt，guess=read的步骤直至没有出现负激发能并且波长也在可控范围内？

zjxitcc

千面追风 发表于 2024-7-25 17:58
如果是波函数不稳定的话，我能做的解决方案是否是将tddft计算的结果作为初猜反复进行stable=opt，guess=r ...

你这个问题的答案不是简单的 是/否。我怀疑你这个结果有点问题，如果你读取了稳定基态波函数，并且基态波函数<S**2>大于0，那你每个激发态都会有自旋污染，并且难以分辨什么是T1，且此时TD()括号中的triplet关键词不起作用。而你在3L展示的结果却都是三重态自旋纯态，这有两种可能：（1）你在3L展示的结果 对应的计算是错的，由于没展示计算细节，无法判断哪里错了；（2）关键词可能没错，但激发态结构优化几步之后，基态就变成了闭壳层单重态，<S**2>等于0，激发态也就没有了自旋污染，这需要看输出文件判断。这些信息你都没交待，因此无法回答（注意，如果你要回复我，应该上传输出文件，并且自己先看输出文件计算过程，有一个初步结论，而不是甩给我一两个log文件让我帮你看，因为我不是你导师）。

千面追风

zjxitcc 发表于 2024-7-25 20:18
你这个问题的答案不是简单的 是/否。我怀疑你这个结果有点问题，如果你读取了稳定基态波函数，并且基态波 ...

老师您好：对应的输出/输入文件已全部上传至1楼。
至于结果问题：我发现是因为没有成功读取到oldchk文件导致的问题，十分抱歉占用了老师这一段的分析时间
在仔细修改输入文件并查看优化的输出文件后，我认为本次的输入结果（还是出现了负激发态）的可能性大概率为老师您所说的第二种情况（激发态仅优化一步后<S**2>=0）

老师您觉得出现这种现象的原因是因为一开始stable=opt的波函数没找对，还是因为所使用的泛函基组所含HF成分较低？
我查找了论坛较多的帖子，总结出了三个可能的解决方案:(前两个是对ghost态的解决方法)
1.用M062X和wB97XD等HF成分较高的基组进一步进行计算
2.用TDA再算一下（排除波函数不稳定的因素）
3.重复使用stable=opt, guess=read步骤直至不出现负激发态

老师认为以上这三个步骤是否可行？以及我这个结果是否可以归类到ghost态里面去？万分感谢老师的仔细答复！

zjxitcc

千面追风 发表于 2024-7-26 09:38
老师您好：对应的输出/输入文件已全部上传至1楼。
至于结果问题：我发现是因为没有成功读取到oldchk文件 ...

你在1L的附件压缩包里给出的stable=opt输入文件.gjf, stable=opt输入文件.log
里面关键词不一致，这不是个好习惯。上传的文件应具有一致性，避免绕晕读者。假设我们以log文件为准，你的关键词写的没有大问题。基组6-31G是不应该，时间不够不能作为借口，还是得加极化函数；由于基态是开壳层单重态，opt()中的triplets不起作用（写了也白写，但不算错误）。结构优化几步之后体系波函数变成了闭壳层单重态，波函数不稳定，所以又有负激发能了。

你说的三个措施均不太可行，且有描述错误。首先要把基组换成合理的，例如6-31G(d,p)。
1. M062X和wB97XD是泛函，不是基组。我不对换泛函做评论，这是你的自由选择。
2. 碰到TDDFT基态波函数不稳定 就换TDA，不能“排除波函数不稳定的因素”，是有些人不想解决问题，选择逃避或假装看不见的措施。
3. 在合理的泛函和基组下，这是一条可行的措施。

如果体系基态是单重态，你想优化T1态极小点结构，最便捷的方式是进行三重态UDFT计算，不涉及TDDFT。可能你被审稿人要求使用TDDFT的td(triplet)再算一下、比较一下结果，但目前的计算过程看来不是很乐观。若基组提升至6-31G(d,p)后还有这种负激发能行为，得换用MRSF-TDDFT方法对T1进行结构优化。