计算化学公社

标题: 利用Gaussian分解键长/键角/二面角对重组能的贡献 [打印本页]

作者
Author:
wal    时间: 2025-11-27 22:22
标题: 利用Gaussian分解键长/键角/二面角对重组能的贡献
本帖最后由 wal 于 2025-12-12 21:45 编辑

工作原理
利用Gaussian的freq=intmode功能打印冗余内坐标对正则模式的贡献;
利用FCHT的print=huangrhys计算正则模式对重组能的贡献;
整合信息来计算冗余内坐标对重组能的贡献。
这边写不了latex,公式写在GitHub的readme部分

计算部分
Gaussian计算输入:
  1. %oldchk=opt_Azulene.chk
  2. #p geom=allcheck freq(readfc,fcht,readfcht,intmode) IOp(7/75=-1)

  3. initial=source=chk final=source=chk spectroscopy=onephotonemission print=(huangrhys,matrix=JK)

  4. td_Azulene.chk
  5. opt_Azulene.chk
复制代码
Huang-Rhys因子计算的输入参考了kalinite老师的帖子,内坐标贡献计算的输入参考了此贴sob老师和坛友Uus/pMeC6H4-/キ的回答。
Tips:


常见报错
  1. ERROR: Low progression after class xxx. Total convergence = xxx%.
复制代码
这个报错可无视,只要Huang-Rhys因子正常打印就能正常分解重组能
  1. FileIo operation on non-existent file.
复制代码
这个报错九成概率是chk位置放反引起的,应首先检查initial态和final态chk文件的位置有没有放反。如果确信无误,检查chk文件是否损坏

分解重组能
运行完成后,用脚本处理log文件:
  1. ./intmode huang-rhys.log
复制代码
首先输出了所有冗余内坐标的贡献:
  1. ====================================================================================================
  2. Internal Coordinate Contribution to Reorganization Energy
  3. ====================================================================================================

  4. Detailed Contributions by Coordinate Type:
  5. ----------------------------------------------------------------------------------------------------

  6. R (Bond Lengths/Angles/Dihedrals):
  7.    Coord         Definition Contribution(%)
  8. ---------------------------------------------
  9.      R12              R(6,7)         2.9900
  10.       R6              R(3,5)         2.9900
  11.      R11             R(5,13)         2.9407
  12. ....


  13. A (Bond Lengths/Angles/Dihedrals):
  14.    Coord         Definition Contribution(%)
  15. ---------------------------------------------
  16.      A12           A(6,4,11)         3.3113
  17.       A6            A(3,2,9)         3.3041
  18.      A11           A(1,4,11)         3.2275
  19. ....


  20. D (Bond Lengths/Angles/Dihedrals):
  21.    Coord         Definition Contribution(%)
  22. ---------------------------------------------
  23.      D19         D(11,4,6,7)         0.9518
  24.      D11          D(9,2,3,5)         0.9518
  25.       D8         D(8,1,4,11)         0.8774
  26. ....


复制代码
随后求了总和,这个是饼图的数据:
  1. ====================================================================================================
  2. Summary: Total Contribution by Coordinate Type
  3. ====================================================================================================

  4.      Coordinate Type        Total %              Description
  5. ------------------------------------------------------------
  6.     R (Bond Lengths)          33.62%
  7.      A (Bond Angles)          58.93%
  8. D (Dihedral Angles)           7.24%
  9. ============================================================

  10. Total (should be ~100%): 99.79%

  11. Internal coordinate contributions saved to: intmode.dat
  12. Reorganization energy data saved to: lambda.dat
复制代码
详细内容存在intmode.dat。另外顺便把分解重组能的结果整理了一下,在lambda.dat。

用输出的dat文件里的数据就能画出来论文里常见的这种图了:
(, 下载次数 Times of downloads: 35) (, 下载次数 Times of downloads: 38)

下载
源码 (, 下载次数 Times of downloads: 35)
也可以去GitHub:bane-dysta/intmode
预编译版和例子 (, 下载次数 Times of downloads: 83)







作者
Author:
sobereva    时间: 2025-11-28 10:52
我将在下一届北京科音高级量子化学培训班讲振动电子光谱的幻灯片里提及此程序
作者
Author:
wal    时间: 2025-11-28 19:28
sobereva 发表于 2025-11-28 10:52
我将在下一届北京科音高级量子化学培训班讲振动电子光谱的幻灯片里提及此程序

哇,被高级班收录了,那我找时间修缮一下,这帖子昨天发的有点随手了
作者
Author:
谢积环    时间: 2025-12-2 15:20
wal 发表于 2025-11-28 19:28
哇,被高级班收录了,那我找时间修缮一下,这帖子昨天发的有点随手了

大佬,太强了。之前想做这个的时候,使用IOp(7/75=-1)内部选项没有效果,原来是要在Gaussian计算振动分辨电子光谱中设置呀。现在看看那个IOp的描述明明都提示了“正常模式位移”,属实被自己蠢到了。
作者
Author:
wal    时间: 2025-12-2 16:24
谢积环 发表于 2025-12-2 15:20
大佬,太强了。之前想做这个的时候,使用IOp(7/75=-1)内部选项没有效果,原来是要在Gaussian计算振动分辨 ...

啊其实不一定要在振动光谱里呀。这个IOp是配合freq=intmode用的,他俩单独搭配也会生效的
作者
Author:
谢积环    时间: 2025-12-3 00:03
wal 发表于 2025-12-2 16:24
啊其实不一定要在振动光谱里呀。这个IOp是配合freq=intmode用的,他俩单独搭配也会生效的

哦哦,那我再试试了
作者
Author:
Bamind    时间: 2025-12-3 10:35
老师好,我想请教一下关于Gaussian分解键长/键角/二面角对重组能的贡献,我的计算步骤如下:先将分子根据该关键词进行基态和激发态计算,#p opt freq b3lyp/6-311g(d) scrf=(solvent=toluene)、#p opt freq td=(nstates=5,root=1) cam-b3lyp/6-311g(d)
scrf=(solvent=toluene)
然后根据老师的输入进行Huang-Rhys因子计算,结果如附图1,结果显示出现错误,随后用老师写的脚本进行分解,结果显示重组能为零。如附图2,请问老师,整个过程哪里出现了问题?

作者
Author:
北大-陶豫    时间: 2025-12-3 10:53
谢积环 发表于 2025-12-2 15:20
大佬,太强了。之前想做这个的时候,使用IOp(7/75=-1)内部选项没有效果,原来是要在Gaussian计算振动分辨 ...

这个一般翻译成“简正模式”而不是“正常模式”吧
作者
Author:
wal    时间: 2025-12-3 11:37
本帖最后由 wal 于 2025-12-3 11:39 编辑
Bamind 发表于 2025-12-3 10:35
老师好,我想请教一下关于Gaussian分解键长/键角/二面角对重组能的贡献,我的计算步骤如下:先将分子根据该 ...

目前我电脑不在手边,不过根据我的经验,file io类报错大概率是td和opt的chk位置放反引起的,你可以检查一下有没有此类问题。如果确信无误,我晚上回去给你看下
作者
Author:
杲杲出日    时间: 2025-12-3 13:01
老师您例子中计算的S0与S1态结构的重组能包含了S0→S1和S1→S0两部分的弛豫能,那您的程序可以计算分解键长/键角/二面角对S1→S0弛豫能的贡献吗?
作者
Author:
Bamind    时间: 2025-12-3 15:02
wal 发表于 2025-12-3 11:37
目前我电脑不在手边,不过根据我的经验,file io类报错大概率是td和opt的chk位置放反引起的,你可以检查 ...

看了一下我的输入文件,确实放反了
作者
Author:
wal    时间: 2025-12-3 16:32
杲杲出日 发表于 2025-12-3 13:01
老师您例子中计算的S0与S1态结构的重组能包含了S0→S1和S1→S0两部分的弛豫能,那您的程序可以计算分解键长 ...

没有吧,spectroscopy=onephotonemission就是纯计算S1→S0啊,onephotonabsorption才是S0→S1
作者
Author:
杲杲出日    时间: 2025-12-3 16:45
wal 发表于 2025-12-3 16:32
没有吧,spectroscopy=onephotonemission就是纯计算S1→S0啊,onephotonabsorption才是S0→S1

我明白了,谢谢老师。那想分解anion和S0之间重组能,也是相同的操作吗
作者
Author:
Bamind    时间: 2025-12-3 16:49
wal 发表于 2025-12-3 11:37
目前我电脑不在手边,不过根据我的经验,file io类报错大概率是td和opt的chk位置放反引起的,你可以检查 ...

老师,更正chk文件的位置之后出现了新的问题,显示
ERROR: Low progression after class 2. Total convergence =  0.0%.
        The vibronic spectrum will likely be unreliable. Stopping.
随后查看了Kalinite老师的帖子 基于Gaussian FCHT计算分解重组能和计算Huang-Rhys因子的Python小脚本
http://bbs.keinsci.com/forum.php ... 8&fromuid=70892
(出处: 计算化学公社)
大家也有类似的问题,Kalinite老师回复说若只需要重组能和Huang-Rhys因子,无视该报错即可。想请问一下,如果继续用脚本处理分解键长/键角/二面角对重组能的贡献能进行下去吗?
作者
Author:
wal    时间: 2025-12-3 16:53
Bamind 发表于 2025-12-3 16:49
老师,更正chk文件的位置之后出现了新的问题,显示
ERROR: Low progression after class 2. Total conve ...

对的 这是振动分辨光谱的报错 无需在意 分解重组能只要Huang-Rhys正常打印即可
作者
Author:
wal    时间: 2025-12-3 16:54
杲杲出日 发表于 2025-12-3 16:45
我明白了,谢谢老师。那想分解anion和S0之间重组能,也是相同的操作吗

这个我没做过 不过底层逻辑应该是类似的
作者
Author:
Bamind    时间: 2025-12-3 16:57
wal 发表于 2025-12-3 16:53
对的 这是振动分辨光谱的报错 无需在意 分解重组能只要Huang-Rhys正常打印即可

谢谢老师,继续分解得到了结果,特别好用
作者
Author:
杲杲出日    时间: 2025-12-3 17:08
wal 发表于 2025-12-3 16:54
这个我没做过 不过底层逻辑应该是类似的

好嘞,谢谢老师
作者
Author:
wxyhgk    时间: 2025-12-4 09:11
这个原来搞过,实际对开发实际分子作用不是很大,只能处理点儿很简单的分子。
作者
Author:
sxkl    时间: 2026-1-12 17:53
你好,我在window系统下采用http://bbs.keinsci.com/forum.php ... 8&fromuid=63020的进行黄昆因子相关计算,然后想采用你上面的脚本分解键长/键角/二面角对重组能的贡献,但是脚本一直编译异常,计算机这块比较弱,能否指导一下window下的脚本运行呀,万分感谢

作者
Author:
sxkl    时间: 2026-1-12 18:08
sxkl 发表于 2026-1-12 17:53
你好,我在window系统下采用http://bbs.keinsci.com/forum.php?mod=viewthread&tid=28708&fromuid=63020的 ...

已解决
作者
Author:
wal    时间: 2026-1-12 18:25
sxkl 发表于 2026-1-12 17:53
你好,我在window系统下采用http://bbs.keinsci.com/forum.php?mod=viewthread&tid=28708&fromuid=63020的 ...

需要编译一个Windows版么?我觉得这种计算一般都做成工作流了,就没做win版。
作者
Author:
wal    时间: 2026-1-13 09:27
楼上好像有人需要,mingw编译了一个window版,测试是没问题的
作者
Author:
yzh    时间: 2026-3-31 10:45
本帖最后由 yzh 于 2026-3-31 11:46 编辑

请教一下老师:
      如果需要计算类似图1过程中的λ(I)和λ(II)的重组能分解,其中λ(I)就是onephotonemission过程,初td_Azulene ,末opt_Azulene;
        λ(II)就是onephotonabsorption过程,初opt_Azulene ,末td_Azulene。      这样理解对么?http://sobereva.com/330中sober老师是直接用S0 opt freq和S1 opt freq进行的dushin分解重组能,和您这个过程是否也是一样的?
      谢谢!



作者
Author:
wal    时间: 2026-3-31 11:42
yzh 发表于 2026-3-31 10:45
请教一下老师:
      如果计算激发态弛豫的重组能(S1-S0),选择onephotonemission, 初态应该S0垂直吸收 ...

分解重组能要的是freq的信息,所以只需要两个opt freq的chk。算吸收/发射只要把相应初态/末态chk放对位置即可。
作者
Author:
yzh    时间: 2026-3-31 11:46
wal 发表于 2026-3-31 11:42
分解重组能要的是freq的信息,所以只需要两个opt freq的chk。算吸收/发射只要把相应初态/末态chk放对位置 ...

谢谢老师,我看了您的chk信息,确实是只需要S0 opt freq和S1 opt freq。你和sober老师330帖子的dushin分解重组能的过程和结果应该是一致的吧?谢谢
作者
Author:
yzh    时间: 2026-4-1 11:12
wal 发表于 2026-3-31 11:42
分解重组能要的是freq的信息,所以只需要两个opt freq的chk。算吸收/发射只要把相应初态/末态chk放对位置 ...

谢谢,请教一下老师,那开头的oldchk%算吸收也要同步更换么?吸收改为oldchk%=td_Azulene.chk?
作者
Author:
yzh    时间: 2026-4-1 11:57
杲杲出日 发表于 2025-12-3 16:45
我明白了,谢谢老师。那想分解anion和S0之间重组能,也是相同的操作吗

您好,用这个工具进行电子重组能或者空穴重组能分解能成功么?我尝试了空穴重组能(阳离子和S0),chk文件检测都是successfully,初态和模态都进行了调整。始终报错:
FileIO operation on non-existent file.
FileIO: IOper= 2 IFilNo(1)=-40770 Len=           0 IPos=           0 Q=   23421002834232
Error termination in NtrErr:
NtrErr Called from FileIO.

调整初态和模态chk位置,以及更换oldchk%的信息。
谢谢!                        
作者
Author:
wal    时间: 2026-4-1 18:04
yzh 发表于 2026-4-1 11:12
谢谢,请教一下老师,那开头的oldchk%算吸收也要同步更换么?吸收改为oldchk%=td_Azulene.chk?

需要。
作者
Author:
wal    时间: 2026-4-1 18:04
yzh 发表于 2026-3-31 11:46
谢谢老师,我看了您的chk信息,确实是只需要S0 opt freq和S1 opt freq。你和sober老师330帖子的dushin分 ...

不同程序的输出可能存在一定差异。
作者
Author:
yzh    时间: 2026-4-1 22:13
wal 发表于 2026-4-1 18:04
需要。

谢谢
作者
Author:
yzh    时间: 2026-4-1 22:14
wal 发表于 2026-4-1 18:04
不同程序的输出可能存在一定差异。

好的,非常感谢,我这边尝试电子和空穴重组能得分解,好像还是不行,我再多尝试一下。
作者
Author:
yzh    时间: 2026-4-2 13:23
wal 发表于 2025-12-3 16:53
对的 这是振动分辨光谱的报错 无需在意 分解重组能只要Huang-Rhys正常打印即可

感谢老师,S0-S1的吸收和发光重组能分解自己的结构测试成功了,但是total 只有 98.41%,这个怎么提升到99.5%以上。
谢谢!
     Coordinate Type        Total %              Description
------------------------------------------------------------
    R (Bond Lengths)          15.45%
     A (Bond Angles)          41.49%
D (Dihedral Angles)          41.47%
============================================================

Total (should be ~100%): 98.41%

作者
Author:
wal    时间: 2026-4-2 13:56
yzh 发表于 2026-4-2 13:23
感谢老师,S0-S1的吸收和发光重组能分解自己的结构测试成功了,但是total 只有 98.41%,这个怎么提升到99 ...

加iop(9/40=3)了吗 如果加了 提到4 或者5
作者
Author:
yzh    时间: 2026-4-2 15:08
wal 发表于 2026-4-2 13:56
加iop(9/40=3)了吗 如果加了 提到4 或者5

好的,谢谢,我没有加
作者
Author:
wal    时间: 2026-4-2 15:32
yzh 发表于 2026-4-2 15:08
好的,谢谢,我没有加

不对我搞混了,这里不是9/40,是7/75,就是主楼写的7/75=-1,这里加9/40应该没啥用
作者
Author:
yzh    时间: 2026-4-3 09:07
wal 发表于 2026-4-2 15:32
不对我搞混了,这里不是9/40,是7/75,就是主楼写的7/75=-1,这里加9/40应该没啥用

谢谢 ,应该是在激发态S1优化中加的
作者
Author:
yzh    时间: 2026-4-8 17:25
老师,再次请教一下
提示ERROR: One imaginary frequency in initial state. You can use the option "ClearImFreq" to remove it  automatically.
     我再intmode.gjf的关键词输入RedDim=ClearImFreq,提示语法错误。
这一般如何解决?
谢谢
作者
Author:
wal    时间: 2026-4-8 18:15
yzh 发表于 2026-4-8 17:25
老师,再次请教一下
提示ERROR: One imaginary frequency in initial state. You can use the option "Cle ...

http://sobereva.com/223
作者
Author:
yzh    时间: 2026-4-9 10:44
wal 发表于 2026-4-8 18:15
http://sobereva.com/223

尝试了,还是要重新优化去掉S1的虚频,虽然只有一个。谢谢
作者
Author:
adreamy    时间: 2026-4-11 14:42
请问老师,您这里用spectroscopy=onephotonemission计算S1→S0的Huang-Rhys因子,这个HR因子是基于基态频率计算得到的吧?为什么不是用spectroscopy=onephotonabsorption来计算S1→S0的HR因子,然后分解得到重组能。我不是很懂这个问题,麻烦老师有空可以告知,谢谢。
作者
Author:
wal    时间: 2026-4-11 23:04
adreamy 发表于 2026-4-11 14:42
请问老师,您这里用spectroscopy=onephotonemission计算S1→S0的Huang-Rhys因子,这个HR因子是基于基态频率 ...

取决于你的目的。如果你需要计算absorption,你可自行更改输入。
作者
Author:
2182356089    时间: 2026-4-15 09:41
老师您好,我想请问您一下这种通过正则振动分解计算重组能与四点法计算出来的重组能结果会有多大区别呢,我看文献大多都是用四点法计算的结果。那这两个计算方法哪个更标准呢,或者这两个方法计算出来的结果在是否都应该是在理论上一致的。
作者
Author:
wal    时间: 2026-4-15 12:45
本帖最后由 wal 于 2026-4-15 13:28 编辑
2182356089 发表于 2026-4-15 09:41
老师您好,我想请问您一下这种通过正则振动分解计算重组能与四点法计算出来的重组能结果会有多大区别呢,我 ...

如果目的只是求重组能的话,四点法更准。
分解重组能的目的通常不是求总重组能,比如我做这个是为了研究振动模式和电子跃迁的耦合程度,这个四点法给不出结论。
作者
Author:
2182356089    时间: 2026-4-15 14:39
wal 发表于 2026-4-15 12:45
如果目的只是求重组能的话,四点法更准。
分解重组能的目的通常不是求总重组能,比如我做这个是为了研究 ...

老师,我还想在请教您一下。那我能不能这样理解,就是虽然也可以通过正则振动的方式计算出重组能的值(既您给的程序或者Dushin),但是该方法主要是为了计算Huang-Rhys因子。而正常计算重组能的话,还是应该通过四点法。
作者
Author:
yzh    时间: 2026-4-15 17:20
2182356089 发表于 2026-4-15 14:39
老师,我还想在请教您一下。那我能不能这样理解,就是虽然也可以通过正则振动的方式计算出重组能的值(既 ...

和四点法的计算结果差异不大的,我试过了。这个主要算S0-S1和S1-S0的重组能,S1-T1不行,UDFT优化的T1提示自旋多重度不一致和S1。
作者
Author:
wal    时间: 2026-4-15 18:35
2182356089 发表于 2026-4-15 14:39
老师,我还想在请教您一下。那我能不能这样理解,就是虽然也可以通过正则振动的方式计算出重组能的值(既 ...

对的,一般计算重组能用四点法就行,走这个在某些情况下有概率产生较大误差
作者
Author:
wal    时间: 2026-4-15 18:36
yzh 发表于 2026-4-15 17:20
和四点法的计算结果差异不大的,我试过了。这个主要算S0-S1和S1-S0的重组能,S1-T1不行,UDFT优化的T1提 ...

小心涉及基团旋转的模式,此类模式有概率导致计算出的重组能偏离四点法的
作者
Author:
2182356089    时间: 2026-4-15 20:17
本帖最后由 2182356089 于 2026-4-15 20:25 编辑
wal 发表于 2026-4-15 18:35
对的,一般计算重组能用四点法就行,走这个在某些情况下有概率产生较大误差

明白了,谢谢老师。然后老师我还有一个不理解的地方,您说“这个计算方法某些情况下有概率产生较大误差”,那如果说该方法在计算重组能上可能产生误差,那这不就是说明最后计算出的振动模式对重组能的贡献上可能就会不准?那最后重组能的分解结果与Huang-Rhys 因子不也是有问题的吗?


作者
Author:
wal    时间: 2026-4-15 22:09
2182356089 发表于 2026-4-15 20:17
明白了,谢谢老师。然后老师我还有一个不理解的地方,您说“这个计算方法某些情况下有概率产生较大误差” ...

这个属于比较底层的问题,我目前对原理的研究没那么深入,就讲讲个人见解吧,各位老师如果看到这里觉得有误请指正
FCHT计算引入的误差分两类,有一部分误差是系统的,每个正则模式都有,所以最终计算出的数值虽然可能有误差,但仍能正确反映贡献比例情况;
还有一种,是由于谐振模型没法描述基团旋转等模式引入的,这种误差属于当前方法的缺陷,如果以后可以进行非谐振的FCHT计算,可能避免此类误差。

作者
Author:
2182356089    时间: 2026-4-15 22:53
wal 发表于 2026-4-15 22:09
这个属于比较底层的问题,我目前对原理的研究没那么深入,就讲讲个人见解吧,各位老师如果看到这里觉得有 ...

好的好的,我明白了,十分感谢老师您的回复。
作者
Author:
yzh    时间: 2026-4-16 10:40
wal 发表于 2026-4-15 18:36
小心涉及基团旋转的模式,此类模式有概率导致计算出的重组能偏离四点法的

非常感谢老师指导




欢迎光临 计算化学公社 (http://bbs.keinsci.com/) Powered by Discuz! X3.3