计算化学公社

标题: 高斯ROHF计算单个原子收敛问题 [打印本页]
作者
Author: wenxiao    时间: 2021-9-4 11:34
标题: 高斯ROHF计算单个原子收敛问题
各位老师好，请教一个问题，我用高斯（ROHF）计算单个Ga原子，试了多种方法都收敛不了（http://sobereva.com/61)，想问一下怎样才能收敛，谢谢！


输入文件如下：


%nprocs=4
%chk=r.chk
%mem=8GB
#p roPBEPBE/def2TZVP pop=full int=ultrafine SCF=NoVarAcc

Title Card Required

0 2
  Ga     0.000 0.000 0.000



作者
Author: zjxitcc    时间: 2021-9-4 12:06
本帖最后由 zjxitcc 于 2021-9-4 12:09 编辑

办法有很多，我提供一种比较自动化的、成功率较大的计算策略：先算UHF，获得UNO，用UNO作为ROHF初猜。
第一步UHF，例如
  1. %chk=Ga.chk
  2. #p UHF/def2TZVP nosymm stable=opt

  4. title

  6. 0 2
  7. Ga   0.0   0.0   0.0

  9. --Link1--
  10. %chk=Ga.chk
  11. #p chkbasis guess(read,only,save,NaturalOrbitals) geom=allcheck nosymm
复制代码

正常结束后，将chk文件转化为fch(k)文件，然后去掉其中多余的beta轨道（因为UNO轨道只需1列，正如ROHF轨道一般）
  1. formchk Ga.chk
  2. fch_u2r Ga.fchk
  3. unfchk Ga_r.fch Ga_rohf.chk
复制代码

formchk和unfchk是高斯自带的小程序。fch_u2r是MOKIT中的一个小程序，开源免费，下载处https://gitlab.com/jxzou/mokit。然后ROHF计算，例如
  1. %chk=Ga_rohf.chk
  2. #p ROHF chkbasis guess=read geom=allcheck nosymm
复制代码

当然，也有其他策略（像慢慢调各种关键词），但我觉得不如这种人力干预少的。


作者
Author: sobereva    时间: 2021-9-4 13:09
如果没有特殊必要，甭用PBEPBE之类纯泛函
ROM06-2X、ROwB97XD直接就能轻松收敛
作者
Author: wenxiao    时间: 2021-9-5 09:09
zjxitcc 发表于 2021-9-4 12:06
办法有很多，我提供一种比较自动化的、成功率较大的计算策略：先算UHF，获得UNO，用UNO作为ROHF初猜。
第 ...

还是收敛不了
作者
Author: wenxiao    时间: 2021-9-5 09:10
sobereva 发表于 2021-9-4 13:09
如果没有特殊必要，甭用PBEPBE之类纯泛函
ROM06-2X、ROwB97XD直接就能轻松收敛

泛函不能改
作者
Author: zjxitcc    时间: 2021-9-5 14:54
本帖最后由 zjxitcc 于 2021-9-5 15:19 编辑
wenxiao 发表于 2021-9-5 09:09
还是收敛不了

哦你问的是PBE纯泛函，我举的例子是ROHF的，确实，ROHF能收敛 不代表 ROPBE也能收敛。那你直接在ROPBE下算也是可以的，例如
%chk=Ga_roPBE_s.chk
%mem=32GB
%nprocshared=16
#p roPBEPBE/def2TZVP nosymm scf(NoVarAcc,vshift)

title

0 2
Ga   0.0   0.0   0.0

这个文件我测试过可以收敛。注意算完看一下轨道形状，确认是否4s^2 4p^1.
作者
Author: wenxiao    时间: 2021-9-5 16:25
zjxitcc 发表于 2021-9-5 14:54
哦你问的是PBE纯泛函，我举的例子是ROHF的，确实，ROHF能收敛 不代表 ROPBE也能收敛。那你直接在ROPBE下 ...

嗯，是可以收敛了，多谢多谢
作者
Author: wenxiao    时间: 2021-9-6 10:44
zjxitcc 发表于 2021-9-5 14:54
哦你问的是PBE纯泛函，我举的例子是ROHF的，确实，ROHF能收敛 不代表 ROPBE也能收敛。那你直接在ROPBE下 ...

调这个有什么技巧吗，算C原子又不能收敛了
%nprocs=4
%chk=r_rohf.chk
%mem=8GB
#p roPBEPBE/def2TZVP nosymm scf(NoVarAcc,vshift)

Title Card Required

0 3
  C     0.000 0.000 0.000
作者
Author: zjxitcc    时间: 2021-9-6 11:04
本帖最后由 zjxitcc 于 2021-9-6 11:05 编辑
wenxiao 发表于 2021-9-6 10:44
调这个有什么技巧吗，算C原子又不能收敛了
%nprocs=4
%chk=r_rohf.chk

把vshift改成vshift=400就行了。经验，唯手熟尔。再提醒一遍，这类计算算完一定要看电子组态是不是跟元素周期表上的一样，或者是否与自己期望的一致。
作者
Author: wenxiao    时间: 2021-9-6 12:38
zjxitcc 发表于 2021-9-6 11:04
把vshift改成vshift=400就行了。经验，唯手熟尔。再提醒一遍，这类计算算完一定要看电子组态是不是跟元素 ...

好的，多谢
作者
Author: 阿甘    时间: 2022-9-9 10:00
本帖最后由 阿甘 于 2022-9-9 11:19 编辑
zjxitcc 发表于 2021-9-6 11:04
把vshift改成vshift=400就行了。经验，唯手熟尔。再提醒一遍，这类计算算完一定要看电子组态是不是跟元素 ...
看回A_25[PBE]诊断的原文时候发现，自己没仔细看清楚A_25[PBE]诊断用的是UKS而不是ROKS……


老师，有个问题想请教您。
我最近在做A_25[PBE]诊断的时候需要计算碳原子在ROPBE下的单点能。用Gaussian16在ROPBE/cc-pVDZ理论水平下得到了收敛的波函数，但是查看其分子轨道发现，那个空的2p轨道的轨道能比两个单占据的2p轨道的轨道能还低，这是不是说明此时波函数不稳定呢？还是不需要理会它呀？随后我又测试了在ROHF/cc-pVDZ理论水平下得到的波函数就没有这个问题。不知老师可否有空帮学生查看一下。提前谢谢老师~

作者
Author: zjxitcc    时间: 2022-9-9 15:56
本帖最后由 zjxitcc 于 2022-9-9 19:12 编辑
阿甘 发表于 2022-9-9 10:00
看回A_25诊断的原文时候发现，自己没仔细看清楚A_25诊断用的是UKS而不是ROKS……
你观察力挺敏锐的。。。ROPBE/cc-pVDZ算三重态C原子结果确实很奇怪，但波函数是稳定的。目前没有更有效的程序用来收敛ROHF/ROKS波函数看是不是有更低的稳定解存在，我只能说没有特殊情况的话建议使用CASSCF或MRCI算开壳层原子。。。如果你想检验ROPBE波函数稳定性也是可以的，在另一个帖子里我介绍过用PSI4检验ROHF波函数稳定性，但是PSI4这个功能对ROKS不支持。不过PySCF支持，比如我们先用高斯算一个gjf
  1. %chk=C_cc-pVDZ_ROPBE.chk
  2. %mem=16GB
  3. %nprocshared=8
  4. #p ROPBEPBE/cc-pVDZ nosymm int=nobasistransform scf(NoVarAcc,vshift=400)

  6. title

  8. 0 3
  9. C   0.0   0.0   0.0
复制代码

然后传轨道给PySCF，即执行
  1. bas_fch2py C_cc-pVDZ_ROPBE.chk
复制代码

上述命令中换成fch(k)文件也行。打开C_cc-pVDZ_ROPBE.py文件，在末尾加上如下几行
  1. dm = mf.make_rdm1()
  2. mf = dft.ROKS(mol)
  3. mf.verbose = 4
  4. mf.xc = 'pbe'
  5. mf.grids.atom_grid = (99,590) # ultrafine
  6. mf.kernel(dm0=dm)
  7. mf.stability()
复制代码

提交给pyscf，2圈收敛，结果为
  1. <class 'pyscf.dft.roks.ROKS'> wavefunction is stable in the internal stability analysis
复制代码

可见该ROPBE波函数是内部稳定的。

作者
Author: scf    时间: 2022-9-10 07:20
qchem里有一种SCF Metadynamics
https://manual.q-chem.com/5.4/sec_SCF-meta-dyn.html
可能能在一定程度上搜索全局极小值
作者
Author: fatpig    时间: 2022-9-10 15:29
有时RO根本就没法找到正确的态，所以也没法收敛。

例如V+的三重态基态，其实是有3个alpha电子+1个beta电子...RO在原理上就没法整出来。
作者
Author: scf    时间: 2022-9-14 19:50
fatpig 发表于 2022-9-10 15:29
有时RO根本就没法找到正确的态，所以也没法收敛。

例如V+的三重态基态，其实是有3个alpha电子+1个beta电 ...

为什么？



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