openmolcas做casscf如何读取轨道初猜

白也

最近尝试使用openmolcas进行casscf计算，先在hf/cc-pvdz水平下产生轨道，openmolcas读取波函数在cc-pvdz基组水平下做了casscf计算，得到了收敛的波函数casscf_pvdz.RasOrb

但接下来想让openmolcas读取这个casscf_pvdz.RasOrb波函数作为初猜，在cc-pvtz基组下作casscf计算，具体输入如下：



casscf部分的设置与在cc-pvdz基组水平下做的casscf设置完全一致，但是输出结果报错，提示：



请问我需要在额外写什么设置吗？还是说读取的参数设置有问题？麻烦有没有大神帮忙解答一下，万分感谢

wxsxjj

直接使用mokit：https://gitlab.com/jxzou/mokit

hebrewsnabla

并不是所有程序都支持大基组读小基组的轨道。有的程序支持是因为自动帮你做了基组投影。

zjxitcc

根据我有限的经验，OpenMolcas目前不支持基组投影，无法在input文件中加一两个关键词实现。不过借助开源程序MOKIT可以轻松实现，举三个例子

例1. 通过molden文件完成基组投影
以H2O分子为例，假设我们已有OpenMolcas计算CASSCF(4,4)/cc-pVDZ的molden文件，在命令行终端 里运行

1. molden2fch h2o_dz.molden -molcas

复制代码

旋即产生h2o_dz.fch文件。molden2fch是MOKIT中的一个小程序，负责将molden文件转化为fch文件。接着启动python，运行

1. from mokit.lib.gaussian import proj2target_basis
   
2. proj2target_basis(fchname='h2o_dz.fch',target_basis='cc-pVTZ',nmo=7)

复制代码

产生h2o_dz_proj.fch文件，其中包含cc-pVTZ基组下的轨道。注意这里的nmo=7，表示水分子CAS(4,4)的3个双占据轨道 + 4个活性轨道 = 7，这是多组态方法，在基组投影时应当考虑活性轨道。你还可以用GaussView/Multiwfn观看轨道，检验与cc-pVDZ下的4个活性轨道是否形状一样。接着运行

1. fch2inporb h2o_dz_proj.fch

复制代码

产生OpenMolcas的输入文件h2o_dz_proj.input和轨道文件h2o_dz_proj.INPORB，坐标、TZ基组和轨道都写好了。我们将input文件末尾4行改成CASSCF(4,4)信息，写法不唯一，例如

1. &RASSCF
   
2. CIMX= 200
   
3. Tight= 5d-8 5d-6
   
4. Spin= 1
   
5. Charge= 0
   
6. nActEl= 4
   
7. RAS2= 4
   
8. FILEORB= h2o_dz_proj.INPORB

复制代码

提交给OpenMolcas，发现收敛较快（因为我们的操作都是正确的）

1.                                          RASSCF iterations: Energy and convergence statistics
   
2.                                          ----------------------------------------------------
   
3. 
   
4.       Iter CI   SX   CI       RASSCF       Energy    max ROT     max BLB   max BLB  Level Ln srch  Step   QN   Walltime
   
5.           iter iter root      energy       change     param      element    value   shift minimum  type update hh:mm:ss
   
6. Nr of preliminary CI iterations:   1
   
7.         1   1   18    1   -75.92491703    0.00E+00  -1.24E-02    4  58 1  7.26E-02*  0.00   0.00    SX    NO    0:00:00
   
8.         2   1   16    1   -75.93122580   -6.31E-03*  3.73E-03    3  56 1  9.29E-03*  0.00   0.00    SX    NO    0:00:00
   
9.         3   1   14    1   -75.93134990   -1.24E-04*  1.53E-03    3  23 1 -2.18E-03*  0.00   0.00    SX    NO    0:00:00
   
10.         4   1   13    1   -75.93137867   -2.88E-05*  9.00E-04    2   7 1 -8.18E-04*  0.00   0.00    SX    NO    0:00:00
    
11.         5   1   13    1   -75.93138472   -6.05E-06*  2.59E-03    3  56 1 -6.52E-04*  0.00   1.40    QN   YES    0:00:00
    
12.         6   1   13    1   -75.93138679   -2.06E-06*  8.77E-04    3  56 1  1.34E-03*  0.00   0.63    QN   YES    0:00:00
    
13.         7   1   12    1   -75.93138855   -1.76E-06*  4.76E-04    3  56 1 -6.92E-04*  0.00   0.62    QN   YES    0:00:00
    
14.         8   1   10    1   -75.93138881   -2.64E-07*  1.15E-04    3  56 1  1.49E-04*  0.00   0.65    QN   YES    0:00:00
    
15.         9   1    9    1   -75.93138884   -2.91E-08*  1.72E-04    2   7 1 -3.29E-05   0.00   2.50    LS   YES    0:00:00
    
16.        10   1    9    1   -75.93138884    2.79E-09   3.32E-05    3  56 1 -1.07E-04*  0.00   0.51    QN   YES    0:00:00
    
17.        11   1    8    1   -75.93138885   -1.62E-08*  7.49E-05    3  56 1 -1.85E-05   0.00   1.25    QN   YES    0:00:00
    
18.        12   1    7    1   -75.93138886   -3.80E-09   7.99E-05    3  56 1 -2.03E-05   0.00   1.33    QN   YES    0:00:00
    
19.       Convergence after 12 iterations
    
20.        13   1    7    1   -75.93138886   -1.58E-09   7.99E-05    3  56 1  7.58E-06   0.00   1.33    QN   YES    0:00:00

复制代码

注意：（1）该技巧只适用于全电子基组的计算，不能使用赝势，因为molden文件中没有赝势信息。（2）不同软件的molden文件无法直接混合使用，因此这里特地强调是OpenMolcas产生的molden文件。如果是其他软件产生的molden文件，例如ORCA，那molden2fch h2o.molden -orca这样指明相应的量化程序即可。

例2. 通过RasOrb/ScfOrb文件完成基组投影
若您没有molden文件，或误删了molden文件，或计算中用到了赝势，可使用此处技巧。这里我们需要从gjf文件出发，先做一个水分子的HF/cc-pVDZ计算

1. %mem=8GB
   
2. %nprocshared=4
   
3. %chk=h2o.chk
   
4. #p RHF/cc-pVDZ nosymm int=nobasistransform scf(xqc,maxcycle=300)
   
5. 
   
6. title
   
7. 
   
8. 0 1
   
9. O      -0.23497692    0.90193619   -0.068688
   
10. H       1.26502308    0.90193619   -0.068688
    
11. H      -0.73568721    2.31589843   -0.068688
    

复制代码

这里需要加关键词nosymm int=nobasistransform。提交给Gaussian，获得h2o.chk文件。依次运行

1. formchk h2o.chk h2o.fch
   
2. fch2inporb h2o.fch
   

复制代码

获得OpenMolcas输入文件h2o.input文件和轨道文件h2o.INPORB。用这两个文件做你想做的CAS计算（轨道局域化、观看挑选轨道、手动调换轨道顺序等繁琐步骤），获得h2o.RasOrb或h2o.RasOrb.1文件。如果你懒得自己算，也可以写一个gjf文件让mokit自动做CASSCF(4,4)/cc-pVDZ计算。总之，现在我们可以运行

1. cp h2o.fch h2o_dz.fch # 备份一下
   
2. orb2fch h2o.RasOrb h2o_dz.fch

复制代码

这样CASSCF(4,4)/cc-pVDZ轨道就保存进了h2o_dz.fch文件。接下来的步骤是类似的，启动python开始基组投影

1. from mokit.lib.gaussian import proj2target_basis
   
2. proj2target_basis(fchname='h2o_dz.fch',target_basis='cc-pVTZ',nmo=7)

复制代码

产生h2o_dz_proj.fch文件，其中包含cc-pVTZ基组下的轨道。注意RasOrb/ScfOrb文件中只有轨道系数，没有体系坐标和基组数据，因此仅靠RasOrb/ScfOrb是无法完成基组投影的，这里我们借助了fch文件来中转。该技巧适用于全电子基组、赝势基组。

例3. 超越基组投影，更高技巧直达cc-pVTZ
直接进行水分子H2O的自动CASSCF(4,4)/cc-pVTZ计算，gjf文件内容如下

1. %mem=8GB
   
2. %nprocshared=4
   
3. #p CASSCF/cc-pVTZ
   
4. 
   
5. mokit{CASSCF_prog=OpenMolcas}
   
6. 
   
7. 0 1
   
8. O      -0.23497692    0.90193619   -0.068688
   
9. H       1.26502308    0.90193619   -0.068688
   
10. H      -0.73568721    2.31589843   -0.068688
    
11. 
    

复制代码

提交任务

1. automr h2o_tz.gjf >h2o_tz.out 2>&1

复制代码

1 min不到结束（体系很小）。OpenMolcas输出文件部分展示

1.                                          RASSCF iterations: Energy and convergence statistics
   
2.                                          ----------------------------------------------------
   
3. 
   
4.       Iter CI   SX   CI       RASSCF       Energy    max ROT     max BLB   max BLB  Level Ln srch  Step   QN   Walltime
   
5.           iter iter root      energy       change     param      element    value   shift minimum  type update hh:mm:ss
   
6. Nr of preliminary CI iterations:   1
   
7.         1   1   10    1   -75.93103040    0.00E+00  -9.02E-03    3  37 1 -1.65E-02*  0.00   0.00    SX    NO    0:00:00
   
8.         2   1    9    1   -75.93135210   -3.22E-04* -2.82E-03    3  37 1 -2.86E-03*  0.00   0.00    SX    NO    0:00:00
   
9.         3   1    9    1   -75.93138182   -2.97E-05*  1.18E-03    3   8 1  9.65E-04*  0.00   0.00    SX    NO    0:00:00
   
10.         4   1    8    1   -75.93138722   -5.41E-06*  4.80E-04    2   9 1 -4.95E-04*  0.00   0.00    SX    NO    0:00:00
    
11.         5   1    7    1   -75.93138838   -1.15E-06*  3.34E-04    3   8 1  1.62E-04*  0.00   1.70    LS   YES    0:00:00
    
12.         6   1    8    1   -75.93138867   -2.95E-07*  2.69E-04    2   9 1  3.65E-04*  0.00   1.12    QN   YES    0:00:00
    
13.         7   1    7    1   -75.93138881   -1.36E-07*  1.98E-04    3   8 1 -8.56E-05   0.00   1.14    QN   YES    0:00:00
    
14.         8   1    6    1   -75.93138884   -3.56E-08*  1.52E-04    2   7 1 -1.71E-05   0.00   1.51    QN   YES    0:00:00
    
15.         9   1    5    1   -75.93138886   -1.31E-08* -1.43E-05    2   9 1 -4.15E-05   0.00   1.00    QN   YES    0:00:00
    
16.        10   1    5    1   -75.93138886   -1.09E-09  -6.14E-06    3   8 1 -1.56E-05   0.00   1.24    QN   YES    0:00:00
    
17.       Convergence after 10 iterations
    
18.        11   1    5    1   -75.93138886   -2.53E-10  -6.14E-06    2   9 1  2.39E-06   0.00   1.24    QN   YES    0:00:00

复制代码

可以看到，CASSCF轨道优化圈数比前述两种技巧的迭代圈数更少，这是因为我们的CASSCF初始轨道是在cc-pVTZ下直接产生的，而非通过cc-pVDZ基组投影。

beefly

>>COPY $CurrDir/casscf_pvdz.RasOrb $Project.RasOrb

&RASSCF
Lumorb
...

印象中openmolcas的rasscf没有基组投影功能，不同基组算的轨道文件不能互相转化

白也

zjxitcc 发表于 2024-11-14 16:11
根据我有限的经验，OpenMolcas目前不支持基组投影，无法在input文件中加一两个关键词实现。不过借助开源程 ...

非常感谢，我尝试一下。
刚才我发现手册中似乎有个关键词EXPBAS，可以调用EXPBAS模块进行基组扩张，不知是否可行

zjxitcc

白也 发表于 2024-11-14 17:12
非常感谢，我尝试一下。
刚才我发现手册中似乎有个关键词EXPBAS，可以调用EXPBAS模块进行基组扩张，不知 ...

不可行/基本没用。https://molcasforum.univie.ac.at/viewtopic.php?id=739

白也

zjxitcc 发表于 2024-11-14 17:23
不可行/基本没用。https://molcasforum.univie.ac.at/viewtopic.php?id=739

尝试了您说的第一种方法，有两点疑问还得麻烦您一下
1. 程序最后一步生成的input文件里把每个原子的basis信息都详细地给出了，看上去不免有些冗杂。改写成basis=cc-pvtz是否可行？
2. input文件中还特意写了noCD关键词，新版默认是开启RICD的，按说会加快电子积分的计算。写noCD是有什么特殊的考虑吗？如果开启会存在什么问题？（比如精度下降明显？）

zjxitcc

白也 发表于 2024-11-14 18:29
尝试了您说的第一种方法，有两点疑问还得麻烦您一下
1. 程序最后一步生成的input文件里把每个原子的basi ...

1. 不可行。不同量化程序内置的基组数据有时候不一样，格式不一样（尽管表面上名称都叫cc-pVTZ），轨道数据需要与基组数据严格对应。
2. 你可以开启，没啥问题，但那样的话可能要手动指定辅助基，我不经常用RICD。

白也

zjxitcc 发表于 2024-11-14 18:58
1. 不可行。不同量化程序内置的基组数据有时候不一样，格式不一样（尽管表面上名称都叫cc-pVTZ）
2. 你 ...

谢谢解答