计算化学公社

标题: 如何用Gauss直接求取并输出键长/键角的一阶二阶导数 [打印本页]

作者
Author:
ZetaFunction    时间: 2023-8-18 21:20
标题: 如何用Gauss直接求取并输出键长/键角的一阶二阶导数
请问各位熟悉Gauss的大佬,Gauss是否可以直接求取能量关于结构中某个键长/键角的一阶导数和二阶导数并输出呢?能量计算用的是DFT方法,导数是有解析表达式的,而且求一阶和二阶导数是Opt和Freq的基础,程序内部肯定是可以算的,我现在就是想要输出这个结果,比如我想知道内坐标下能量关于4号原子和1号原子的距离的一阶和二阶导数的值具体是多少。

我翻阅了Gauss官网的手册,上面说Opt任务的输出文件会给出每一步的一阶导数,但是我现在需要的是一个单点结构的导数,并不想做优化(否则所有一阶导数就都变成0了)。手册上说Freq会输出二阶导数,但是我检查输出文件发现是关于笛卡尔坐标的导数,虽然也可以手动转换但还是有点麻烦最好可以直接输出。手册上还说在输入文件末尾写上“B 4 1 D”就会计算4号原子和5号原子距离的二阶导数,但是我照做后翻阅了整个输出文件都未见计算结果。

  1. %nprocshared=40
  2. %mem=160GB
  3. %chk=BD-Der.chk
  4. # M062X/gen geom=modredundant

  5. Der

  6. -1 1
  7. C                  0.71241700   -0.23209900    0.35036600
  8. O                  1.21233200   -1.20834600   -0.22999200
  9. H                  0.54006600   -0.28116200    1.44293400
  10. C                 -1.27239600   -0.05182500    0.04828100
  11. N                 -2.42123100    0.00930900   -0.09903900
  12. C                  1.08316200    1.17684700   -0.10635200
  13. H                  0.48030600    1.93470600    0.39556400
  14. H                  2.14204800    1.33286000    0.12512100
  15. H                  0.94843600    1.25766400   -1.18418300

  16. B 4 1 D
  17. A 4 1 2 D
  18. D 3 1 2 4 D
  19. A 5 4 1 D
  20. A 4 1 6 D

  21. @ma-TZVP.txt
复制代码



希望有大佬能帮我解答,如果Gauss确实做不到,有没有其他的量化程序可以实现这样的功能?

作者
Author:
mizu-bai    时间: 2023-8-19 00:46
force nosymm 可以满足你的需求
作者
Author:
sobereva    时间: 2023-8-19 05:01
不要把Gaussian写成Gauss

用freq和#P可输出,参考北京科音基础量子化学培训班(http://www.keinsci.com/workshop/KBQC_content.html)的ppt

(, 下载次数 Times of downloads: 10)


作者
Author:
ZetaFunction    时间: 2023-8-19 19:59
sobereva 发表于 2023-8-19 05:01
不要把Gaussian写成Gauss

用freq和#P可输出,参考北京科音基础量子化学培训班(http://www.keinsci.com/ ...

感谢Sobereva老师的解答,我在输入文件写上“#P freq”后逐行检查输出的log文件,还是只能找到笛卡尔坐标下的海森矩阵和一阶导数,并未找到内坐标下定义的海森矩阵。另外海森矩阵的单位是Hartree/Bohr^2,请问如何转换到Hartree/A^2?
作者
Author:
ahxb    时间: 2023-8-19 22:57
输入文件中写上freq=modredundant后,chk文件保存了内坐标下的受力和海森矩阵。计算完毕后使用formchk将chk文件转化为fchk文件,可在fchk文件中看到Internal Forces和Internal Force Constants字段(以下三角矩阵形式表示,可与输出文件中的结果比对,但有效数字更多)。冗余内坐标定义可在输出文件开头附近看到,fchk文件中也有Redundant internal coordinate indices和Redundant internal coordinate字段,前者是冗余内坐标定义,四个一组;后者是冗余内坐标数值。
也可以执行rwfdump ***.chk 511.txt 511R,在文本文件511.txt末尾会给出一个数组,第一个数是SCF收敛后的能量;第四个数开始是在内坐标中的受力,长度为冗余内坐标个数NVAR;然后是内坐标下的二阶力常数(以下三角矩阵形式表示),长度为NVAR*(NVAR+1)/2。你要是还想获得更高精度的结果,还可以从rwfdump输出的文本文件中Number为511的部分读取Base和End,然后自行写程序读取二进制文件中从Base字节到End字节的部分,与文件最后给出的结果完全相同,是Gaussian计算出的原始结果,精度只受双精度浮点数本身的限制,一般可达15~16位有效数字。
1 Hartree/Bohr^2 = 1 Hartree/(0.529 angstrom)^2 = 3.57 Hartree/angstrom^2,想要更精确的转换关系可以去查Bohr的高精度定义。
作者
Author:
sobereva    时间: 2023-8-20 04:28
ZetaFunction 发表于 2023-8-19 19:59
感谢Sobereva老师的解答,我在输入文件写上“#P freq”后逐行检查输出的log文件,还是只能找到笛卡尔坐标 ...

拿具体输入输出文件说事。ppt里该说的都说了




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