高斯ROHF计算单个原子收敛问题

wenxiao

各位老师好，请教一个问题，我用高斯（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

zjxitcc

办法有很多，我提供一种比较自动化的、成功率较大的计算策略：先算UHF，获得UNO，用UNO作为ROHF初猜。
第一步UHF，例如

1. %chk=Ga.chk
   
2. #p UHF/def2TZVP nosymm stable=opt
   
3. 
   
4. title
   
5. 
   
6. 0 2
   
7. Ga   0.0   0.0   0.0
   
8. 
   
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
   

复制代码

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

sobereva

如果没有特殊必要，甭用PBEPBE之类纯泛函
ROM06-2X、ROwB97XD直接就能轻松收敛

wenxiao

zjxitcc 发表于 2021-9-4 12:06
办法有很多，我提供一种比较自动化的、成功率较大的计算策略：先算UHF，获得UNO，用UNO作为ROHF初猜。
第 ...

还是收敛不了

wenxiao

sobereva 发表于 2021-9-4 13:09
如果没有特殊必要，甭用PBEPBE之类纯泛函
ROM06-2X、ROwB97XD直接就能轻松收敛

泛函不能改

zjxitcc

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.

wenxiao

zjxitcc 发表于 2021-9-5 14:54
哦你问的是PBE纯泛函，我举的例子是ROHF的，确实，ROHF能收敛 不代表 ROPBE也能收敛。那你直接在ROPBE下 ...

嗯，是可以收敛了，多谢多谢

wenxiao

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

zjxitcc

wenxiao 发表于 2021-9-6 10:44
调这个有什么技巧吗，算C原子又不能收敛了
%nprocs=4
%chk=r_rohf.chk

把vshift改成vshift=400就行了。经验，唯手熟尔。再提醒一遍，这类计算算完一定要看电子组态是不是跟元素周期表上的一样，或者是否与自己期望的一致。

wenxiao

zjxitcc 发表于 2021-9-6 11:04
把vshift改成vshift=400就行了。经验，唯手熟尔。再提醒一遍，这类计算算完一定要看电子组态是不是跟元素 ...

好的，多谢

阿甘

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理论水平下得到的波函数就没有这个问题。不知老师可否有空帮学生查看一下。提前谢谢老师~

zjxitcc

阿甘 发表于 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)
   
5. 
   
6. title
   
7. 
   
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波函数是内部稳定的。

scf

qchem里有一种SCF Metadynamics
https://manual.q-chem.com/5.4/sec_SCF-meta-dyn.html
可能能在一定程度上搜索全局极小值

fatpig

有时RO根本就没法找到正确的态，所以也没法收敛。

例如V+的三重态基态，其实是有3个alpha电子+1个beta电子...RO在原理上就没法整出来。

scf

fatpig 发表于 2022-9-10 15:29
有时RO根本就没法找到正确的态，所以也没法收敛。

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

为什么？