求助如何解决CFOUR的UHF收敛到不合理波函数

胡说

想用CFOUR V2的CCSD(T)解析导数来计算一个小体系（开壳层（0，2），6个原子）的电子密度和性质。
基组aug-cc-pVTZ，300基函数。（因为还要计算其阴离子，故带弥散。）
众所周知，CFOUR的SCF收敛模块并不好，再加上体系不常规，所以UHF收敛到了一个能量较高的波函数（和高斯相比，并用高斯做了波函数稳定检测）。
另外这个分子结构具有Cs对称性，但根据高斯stable=opt的结果，波函数是broken-symmetry,不过CFOUR收敛的波函数是symmetry的。

想直接用CFOUR自带的Hartree-Fock Stability Analysis（HFSTABILITY=2）来找到稳定基态波函数，关键词如下：
*CFOUR(CALC=HF,REFERENCE=UHF,BASIS=AUG-PVTZ,MEM_UNIT=GB,MEMORY_SIZE=24
HFSTABILITY=2,COORDINATES=CARTESIAN,SCF_DAMPING=500,SCF_MAXCYC=300
CHARGE=0,MULTIPLICITY=2,SYMMETRY=OFF)

但是程序一直卡在xintprc不动了，内存占用奇高，硬盘使用近300GB：linux-top如下
25.611g 5.674g   4992 S   0.0 17.8   2:28.41 xintprc
25.611g 5.674g   5000 S   0.0 17.8   2:34.53 xintprc
25.611g 5.629g   4992 S   0.0 17.6   2:11.26 xintprc
输出文件停在：
@GMOIAA-I, Processing MO integrals for spin case BB.
@GMOIAA-I, Generation of integral list completed.
                     TYPE                NUMBER
                     ----               --------
                     PPPP                     0
                     PPPH             374543984
                     PPHH              65843530
                     PHPH              33781826
                     PHHH              11876548
                     HHHH                536113

                    TOTAL             486582001
@GMOIAB-I, Processing MO integrals for spin case AB.
                     TYPE                NUMBER
                     ----               --------
                     PPPP                     0
                     PPPH1H           381365745
                     PPPH2H           371617862
                     PPHH             134073848
                     PHPH1P            33517952
                     PHPH2P            35282846
                     PHHH1P            12092873
                     PHHH2P            12404383
                     HHHH               1118810

                    TOTAL             981474319
虽然卡了十几个小时，但程序倒没有完全死绝。还在算，就是很慢很慢。

这几天使用CFOUR碰到一些问题，想请教讨论下：
1.为什么CFOUR做波函数检测如此耗费资源？速度还跟蜗牛一样（比CCSD（T）还费资源，慢）
2.请问大家对此问题有没有什么解决方法吗？
3.另外有没有办法让CFOUR读取高斯或其他程序的波函数作为初猜（我感觉没有办法）。
4.另外发现CFOUR算CC，有时候CPU利用很低，经常进程进入S，就用了两三个核，不至于内存和硬盘成为瓶颈。如此CC速度可能还不如高斯。
5.还有发现CFOUR计算产生的波函数MOLDEN，里面记录的是笛卡尔型的基函数（轨道数目等于笛卡尔型基函数数目）。但明明计算时用的是球谐型的（默认就是）。这个对后续某些Multiwfn波函数分析是有影响的吧。有办法转换成球谐型的吗？
6.还有个问题，发现输入坐标和CFOUR计算时使用坐标不一致（不是单位问题），SYM=OFF关闭对称也会如此。CFOUR也会像高斯一样自己将分子朝向改变吗，如何避免呢？7. 笛卡尔坐标输入可以正确判断出对称性，可Z坐标（gv保存的）输入却无法判断对称性了？
谢谢！

zjxitcc

chemjay 发表于 2023-9-14 09:09
用fch2cfour转换得到的ZMAT、ECPDATA、GENBAS、OLDMOS后直接运行计算，发现CFOUR并没有读入波函数，请问 ...

注意CFOUR版本要高于等于v2.1。更早的CFOUR v1.0和v2.00beta无法从OLDMOS读入轨道。

zjxitcc

chemjay 发表于 2023-9-14 09:09
用fch2cfour转换得到的ZMAT、ECPDATA、GENBAS、OLDMOS后直接运行计算，发现CFOUR并没有读入波函数，请问 ...

你是不是看到它SCF 1圈收敛，太快了，以为没计算？

（1）判断CFOUR是否正确读入产生的轨道，看输出文件里有无这两行
  starting vectors read from OLDMOS
@INITGES-I, Nonstandard initial guessing procedure.

（2）书写gjf文件时要用CFOUR输出文件里的坐标，在Coordinates used in calculation (QCOMP)这个字段后面，注意这里坐标单位是Bohr，所以我在6L举的gjf例子里写了units=au。

如果经过上述检查还搞不定，并且你的体系不涉密的话，直接把文件（做为附件）发上来我看看。另外有个问题是你要写目标计算方法，比如CALC=CCSD。CALC=SCF是算完SCF就停了。不过这个与SCF步骤是否读入轨道、几圈收敛无关。

chemjay

zjxitcc 发表于 2023-2-18 01:14
现可使用fch2cfour小程序从高斯向CFOUR传坐标、轨道和基组，CFOUR的SCF初猜极其差，如何收敛、收敛到哪里去 ...

用fch2cfour转换得到的ZMAT、ECPDATA、GENBAS、OLDMOS后直接运行计算，发现CFOUR并没有读入波函数，请问是需要对ZMAT进行修改吗？下面是我的部分输入文件。
*CFOUR(CALC=SCF,REF=UHF,SYM=OFF,CHARGE=0,MULTI=2,BASIS=SPECIAL,ECP=ON
CC_PROG=ECC,ABCDTYPE=AOBASIS)

胡说

zjxitcc 发表于 2023-2-18 01:14
现可使用fch2cfour小程序从高斯向CFOUR传坐标、轨道和基组，CFOUR的SCF初猜极其差，如何收敛、收敛到哪里去 ...

谢谢老师，当年这个问题没解决，后来就放弃了。
虽然现在不用CFOUR了，以后要用的话试试您说的这个方法。
MOKIT是个好工具，经常用它来将Gaussian轨道传给ORCA。

zjxitcc

现可使用fch2cfour小程序从高斯向CFOUR传坐标、轨道和基组，CFOUR的SCF初猜极其差，如何收敛、收敛到哪里去是个玄学，难以处理复杂体系。由于你没给出分子，我这里以一个Cs点群的H-O-Cl三原子分子为例，电荷-1，自旋多重度2，基组aug-cc-pVTZ，其UHF很有挑战性。写个hocl.gjf文件如下

1. %chk=hocl.chk
   
2. %mem=56GB
   
3. %nprocshared=28
   
4. #p UHF/aug-cc-pVTZ nosymm int=nobasistransform scf(maxcycle=128,xqc) stable=opt units=au
   
5. 
   
6. title
   
7. 
   
8. -1 2
   
9. O   -2.08303774     0.11989815     0.0
   
10. H   -2.75920696    -1.56351799     0.0
    
11. Cl   1.03231321    -0.00978036     0.0
    

复制代码

这个坐标数值是我随便跑几秒CFOUR任务，从输出文件复制出来的，因为我现在还弄不清SYM=OFF时CFOUR令人讨厌的自动对分子转向，因此这套数值是CFOUR旋转分子后的直角坐标。高斯前128圈SCF不收敛，关键词xqc让高斯自动切换为scf=qc，随后收敛，并通过波函数稳定性测试，最终E(UHF) = -534.896469981 a.u. 运行

1. fch2cfour hocl.fch

复制代码

产生ZMAT, GENBAS和OLDMOS三个文件，ZMAT文件末尾默认关键词为UHF计算所需，如果有特殊需求，可以自己添加。提交

1. xcfour >OUTPUT 2>&1 &

复制代码

CFOUR的UHF花了7圈收敛，但电子能量始终维持在-534.896470 a.u.，仅密度矩阵有微小变化，说明传轨道成功。如果不带弥散函数的话，一般传轨道过去1-2圈就收敛了。fch2cfour是开源程序MOKIT里的一个小程序，自己安装或下载即可。

胡说

sobereva 发表于 2019-3-3 00:48
看上去就是内存不够，详情可咨询开发者

1 原理上HF能量是可以精确对应的，可以达到不止5、6位。有差异 ...

明白了，多谢老师

sobereva

胡说 发表于 2019-3-2 18:00
谢谢老师。
我试了PSI4，但碰到个问题：
我按照“grad,wfn = gradient('CCSD(T)', return_wfn=True)；f ...

看上去就是内存不够，详情可咨询开发者

1 原理上HF能量是可以精确对应的，可以达到不止5、6位。有差异可能是基组定义、球谐型/笛卡尔型基函数的差异、坐标保留位数、双电子积分忽略的阈值之类细节所致

2 没有什么特殊注意的。只要载入文件后按照Multiwfn手册附录5的方法通过检测就没问题

3 没有任何限制

胡说

sobereva 发表于 2019-2-27 14:47
我没用过CFOUR的波函数稳定性测试，不妨尝试先把基组去掉弥散

PSI4也支持CCSD(T)解析导数，比CFOUR更友 ...

谢谢老师。
我试了PSI4，但碰到个问题：
我按照“grad,wfn = gradient('CCSD(T)', return_wfn=True)；fchk_writer = psi4.FCHKWriter(wfn)；fchk_writer.write('HF_CCSDpT.fchk')”这样的输入来产生波函数文件。虽然HF分子这个例子可以成功结束，但算自己体系的时候，算完CCSD(T)能量后，过了很久进入”DPD File4 Cache Listing：“，然后报错：”Error in: dpd_block_matrix: No memory left.“
换用aug-cc-pvdz也是出现这个错误，设定60GB内存，162基函数.
另外发现PSI4的CC梯度计算不支持冻核近似。

另外有几个问题想请教：
1. 不同程序HF能量会有略微不同，HF能量小数点后第五位（比如PSI4），或第六位(比如CFOUR)会有不同（以高斯为参考），请问这种情况是否可靠？
就是我用高斯的HF能量来check其他程序HF产生的能量是否靠谱，一般是小数点前四位都一样，但第五位或第六位开始会不一样，这样是否就可以认为其他程序HF参考波函数/能量是没问题的？（理论上是应该精确一样的，但实际上即使我把收敛限设成一样，最后产生的HF能量还是会有差异。）

2. 用了冻核近似的CC波函数用于Multiwfn分析，是否有什么要注意的？在手册上没找到相关说明
3. 笛卡尔型基函数的波函数用于Multiwfn分析，有什么限制吗，有什么功能是必须要求是球谐型基函数的？在手册上我好像没找到。
谢谢

sobereva

我没用过CFOUR的波函数稳定性测试，不妨尝试先把基组去掉弥散

PSI4也支持CCSD(T)解析导数，比CFOUR更友好，可以试试

CFOUR的耦合簇效率没有吹的那么神奇

Multiwfn不直接支持CFOUR产生的molden文件，需要先用beefly的molden2aim程序转换一下变成格式标准的molden文件，之后就可以给Multiwfn照常用了