计算化学公社

标题: 利用MOKIT从ORCA向其他量化程序传轨道 [打印本页]
作者
Author: zjxitcc    时间: 2024-4-23 15:43
标题: 利用MOKIT从ORCA向其他量化程序传轨道
本帖最后由 zjxitcc 于 2024-4-23 15:44 编辑

本文在2023年8月21号首发于“量子化学”公众号上，现由本人转载，后续若有更新请到gitlab上查看。本文介绍如何使用MOKIT从ORCA向其他量化程序传轨道，有以下可能的用途：

（1）在ORCA中进行了RIJK或RIJCOSX加速的大体系HF/DFT计算，想传轨道给其他程序进行后续计算，或想产生fch文件方便可视化。
（2）在ORCA中进行了CCSD等post-HF计算，并产生了（自旋）自然轨道，想产生fch文件方便可视化。
（3）在ORCA中进行了CASSCF计算，想传CASSCF轨道给其他程序进行后续计算（例如MC-PDFT）。
（4）有些复杂体系可能ORCA可以较顺利地收敛出特殊SCF解，而目标程序难以得到，可以传轨道给目标程序。

为阅读简便，本文以水分子的RHF/def2-TZVP计算为例，输入文件h2o.inp内容如下
  1. %pal nprocs 2 end
  2. %maxcore 2000
  3. ! RHF def2-TZVP VeryTightSCF RIJCOSX def2/J defgrid3
  4. * xyz 0 1
  5. O    0.74608908   -0.25872442    0.0
  6. H    1.70608908   -0.25872442    0.0
  7. H    0.42563449    0.64621141    0.0
  8. *
复制代码
VeryTightSCF是笔者的个人喜好，但做实际计算至少应采用TightSCF。defgrid3表示使用尽可能高的COSX格点。该例使用了RIJCOSX近似来加速计算，这对轨道系数影响非常小，更何况用了比较精细的格点和严格的收敛限，预期得到的轨道与一个传统RHF/def2-TZVP计算的轨道极其接近。用ORCA算完后运行
  1. orca_2mkl h2o -mkl
复制代码
转化生成h2o.mkl文件，下面我们以该文件为例展示如何传轨道，读者若使用本文功能，请务必顺带阅读文末注意事项。

1. ORCA传轨道给AMESP
运行
  1. mkl2amo h2o.mkl
复制代码
会产生h2o.amo和h2o.aip两个文件，分别含基组和轨道数据，坐标及关键词。为了让AMESP读入轨道，需运行
  1. a2m h2o.amo
复制代码
即将h2o.amo文本文件转化为二进制文件h2o.mo，其中a2m是AMESP自带的小程序。

2. ORCA传轨道给BDF
运行
  1. mkl2bdf h2o.mkl
复制代码
会产生H2O.BAS，h2o.inp和h2o.scforb三个文件，分别含基组数据，坐标及关键词，轨道数据。

3. ORCA传轨道给CFOUR
运行
  1. mkl2cfour h2o.mkl
复制代码
会产生ZMAT, GENBAS和OLDMOS三个文件，注意CFOUR所用文件名是固定的。mkl2cfour小程序目前有个缺点：在书写ORCA输入文件时需要使用CFOUR输出文件里体系的直角坐标，产生的mkl文件才能给mkl2cfour使用。这是因为CFOUR会旋转分子朝向，只有使用CFOUR输出文件里体系的直角坐标，其已经被旋转至程序标准朝向，才往往不会再被旋转了。若遇到CFOUR的SCF不收敛，意味着有输出文件，因此这点不是大问题。

4. ORCA传轨道给Dalton
运行
  1. mkl2dal h2o.mkl
复制代码
会产生h2o.dal和h2o.mol两个文件。

5. ORCA传轨道给GAMESS
运行
  1. mkl2inp h2o.mkl
复制代码
会产生h2o.inp文件，内含坐标，基组数据和一些简单的关键词。

6. ORCA传轨道给Molpro
运行
  1. mkl2com h2o.mkl
复制代码
会产生h2o.a和h2o.com文件，前者含Alpha轨道，后者含坐标，基组和关键词。如果是UHF计算，还会产生h2o.b文件（含Beta轨道）。

7. ORCA传轨道给OpenMolcas
运行
  1. mkl2inporb h2o.mkl
复制代码
会产生h2o.INPORB和h2o.input文件。

8. ORCA传轨道给PSI4
运行
  1. mkl2psi h2o.mkl
复制代码
会产生h2o.A和h2o.inp文件。前者含Alpha轨道，后者含坐标，基组和关键词。如果是UHF类型计算，还会有h2o.B文件（内含Beta轨道）。

9. ORCA传轨道给PySCF
运行
  1. mkl2py h2o.mkl
复制代码
会产生h2o.fch和h2o.py文件。在运行.py文件时会从.fch文件里读取轨道。
10. ORCA传轨道给Q-Chem
运行
  1. mkl2qchem h2o.mkl
复制代码
产生h2o.in文件和一个h2o文件夹。如果机器上有定义环境变量$QCSCRATCH，则h2o文件夹会自动被移至Q-Chem临时文件目录下，屏幕上提示你可以运行
  1. qchem h2o.in h2o.out h2o
复制代码
提交Q-Chem任务，最后一个参数表示读取h2o/目录下的轨道文件。若未定义$QCSCRATCH，h2o文件夹则放在当前目录下，读者需要时自行移动。

11. ORCA传轨道给Gaussian
该功能较重要，有几种不同使用方式，此处重点介绍。

11.1 使用mkl2gjf小程序
运行
  1. mkl2gjf h2o.mkl
复制代码
会产生h2o.gjf文件，内含坐标、基组信息和轨道信息。注意Gaussian是支持从gjf文件中读取轨道的（关键词guess=cards）。

11.2 使用mkl2fch小程序
11.2.1 适用于常见HF/DFT/CASSCF轨道
  1. mkl2fch h2o.mkl
复制代码
若检测到当前目录下无h2o.fch文件，则会从零创建h2o.fch；若已存在h2o.fch文件，则仅将轨道信息写入h2o.fch。

11.2.2 适用于UNO、MP2、CCSD、CASSCF等自然轨道
此处以RI-MP2为例，输入文件计算级别部分为
  1. ! RHF def2-TZVP VeryTightSCF RIJCOSX def2/J defgrid3 RI-MP2 def2-TZVP/C
  2. %mp2
  3. density unrelaxed
  4. natorbs true
  5. end
复制代码
关键词unrelaxed表示计算RI-MP2非弛豫密度对应的自然轨道，若改成relaxed即计算弛豫密度对应的自然轨道，笔者更推荐使用前者进行波函数分析。正常计算完后有一个文件叫h2o.mp2nat，它其实是gbw格式，运行
  1. orca_2mkl h2o.mp2nat h2o_mp2no.mkl -mkl -anyorbs
  2. mkl2fch h2o_mp2no.mkl h2o_mp2no.fch -no
复制代码
即生成含MP2自然轨道和轨道占据数的h2o_mp2no.fch文件。参数-no意为将自然轨道占据数写入生成的fch文件。用GaussView打开该文件可视化轨道时，箭头旁的数值不是轨道能量，而是轨道占据数；用Multiwfn打开该文件可视化轨道时，Ene处数值为轨道占据数，Occ处数值无意义。

11.2.3 适用于UMP2、UCCSD(自旋)自然轨道
这类计算输入文件较复杂，往往涉及到对称破缺单重态和检验波函数稳定性等问题，计算实例见笔者在计算化学公社论坛上的回帖http://bbs.keinsci.com/thread-29762-1-1.html
这里我们假设已经获得了.mdci.nat文件，接着执行
  1. orca_2mkl O_singlet_o.mdci.nat O_singlet_CCSD.mkl -mkl -anyorbs
  2. mkl2fch O_singlet_CCSD.mkl O_singlet_CCSD.fch -nso
复制代码
参数-nso意为将alpha、beta自旋自然轨道占据数写入生成的fch文件。打开该文件可视化轨道时轨道占据数含义同11.2.2。自旋自然轨道有alpha、beta两套，有的读者可能认为不便分析，但自然轨道只有1套。为获得UCCSD自然轨道，可启动Python，运行
  1. from mokit.lib.lo import gen_no_from_nso
  2. gen_no_from_nso(fchname='O_singlet_CCSD.fch')
复制代码
产生UCCSD自然轨道文件O_singlet_CCSD_NO.fch，只含一套轨道。若读者不熟悉自然轨道、自然自旋轨道，可阅读如下三篇  
《从密度矩阵产生自然轨道-理论篇》
《从密度矩阵产生自然轨道_实战篇（上）》  
《在Multiwfn中基于fch产生自然轨道的方法与激发态波函数、自旋自然轨道分析实例》http://sobereva.com/403

注意事项
1. 目前ORCA中开启点群对称性做计算只能旋转分子朝向、识别点群和不可约表示，不能利用对称性加速计算，所以不经常使用。是否开启点群对称性做计算不影响本文功能。但传轨道后在目标程序中应关闭对称性（生成的输入文件里写了这些关键词），笔者没有让各个小程序支持点群对称性的意向。

2. 以上各个小程序都会产生目标程序输入文件，内含坐标和基组数据，不需要用户再去手动书写基组名称。强烈推荐用户使用该文件进行计算，既能免去手写基组的麻烦，也能保证传轨道时对应性更好。

3. ORCA仅支持球谐型基函数，不支持笛卡尔型基函数，因此传轨道后产生的文件采用的也是球谐型基函数进行计算。

4. 若ORCA的SCF有几百、上千a.u.的剧烈振荡，很可能是出现了基函数线性依赖导致的，此时即使侥幸收敛了能量也未必靠谱，需要在输入文件里加上
  1. %scf
  2. sthresh 1e-6
  3. end
复制代码
此阈值是Gaussian和GAMESS处理基函数线性依赖的默认值，比较靠谱。

5. mkl文件有个缺点是不含赝势信息，笔者曾在ORCA官方论坛上发帖建议未来ORCA在mkl文件中增加赝势信息，但ORCA开发者表示目前没兴趣。若读者在计算中使用全电子基组，自然无此问题；若用了赝势，按上文操作产生其他量化程序的文件不会含赝势信息，即使轨道系数正确，SCF也会剧烈振荡。这里笔者推荐一种解决办法：

Step 1. 先用Gaussian算个几秒的任务
这里以二甲基锌（ZnMe2）为例，假设我们对Zn用SDD，gjf内容如下
  1. %chk=ZnMe2.chk
  2. #p RHF genecp nosymm int=nobasistransform guess(only,save)

  4. title

  6. 0 1
  7. [ZnMe2坐标，略]

  9. H C 0
  10. cc-pVDZ
  11. ****
  12. Zn 0
  13. SDD
  14. ****

  16. Zn 0
  17. SDD

复制代码
关键词only表示构建SCF初猜后即退出程序，即使大体系耗时也只有秒级别。若读者想让Gaussian正常算完，获得收敛的RHF轨道，可以删去guess(only,save)。正常结束后获得ZnMe2.chk文件。

Step 2. 获得含赝势的fch文件和ORCA相关文件
运行
  1. fch2mkl ZnMe2.chk
复制代码
随即产生ZnMe2.fch，ZnMe2_o.inp和ZnMe2_o.mkl文件，现在我们就有了一个含赝势信息的fch文件。ZnMe2_o.mkl含有来自Gaussian的轨道，除非有特殊计算目的，否则一般推荐使用该mkl文件。运行
  1. orca_2mkl ZnMe2_o -gbw
复制代码
将mkl转化为gbw文件。我们用ZnMe2_o.inp文件在ORCA中做目标计算，最后把轨道传回已存在的fch文件
  1. orca_2mkl ZnMe2_o -mkl    # 假设ZnMe2_o.gbw已被更新，重新产生mkl
  2. cp ZnMe2.fch ZnMe2_o.fch  # 备份一下，以防覆盖
  3. mkl2fch ZnMe2_o.mkl         # 将轨道传回ZnMe2_o.fch
复制代码
后续可以使用fch2inp，fch2inporb，fch2com，bas_fch2py等小程序传给其他量化程序做计算，且均含赝势信息。mkl2fch检测到有fch文件存在时会实直接使用该文件，而非从零生成，因此避免了赝势缺失问题。注意不能跳过Step 1随便塞一个fch文件给mkl2fch小程序企图欺骗之，那样即使不报错也可能存在基组数据或轨道不对应的错误。

致谢
感谢wzkchem5 ，hebrewsnabla和niobium的修改建议。

相关阅读
《ORCA 5.0安装及运行》
《利用MOKIT从Gaussian向其他量化程序传轨道》
《利用MOKIT从PySCF向其他量化程序传轨道》
MOKIT在线手册https://jeanwsr.gitlab.io/mokit-doc-mdbook
引用MOKIT的发表文章一览https://jeanwsr.gitlab.io/mokit-doc-mdbook/citing.html




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