从gaussian和gamess-ag真人国际厅网站

获得电子积分是自己写量化程序的关键。在网上能找到一些零散的从常用的gaussian或gamess-us程序中提取单、双电子积分的帖子,但是往往语焉不详,说得不确切、不全面,有的做法不雅,有的还会令人白绕弯路。这里把从gaussian和gamess-us中提取电子积分的方式清楚、详细、准确地说一下。

无论输出的是单电子积分还是双电子积分,都是基函数已归一化后的值。

双电子积分形式是(ij|kl),i,j基函数是r1坐标,k,l基函数是r2坐标。由于积分有对称性,即ijkl=ijlk=jikl=jilk=klij=lkij=klji=lkji,所以程序输出双电子积分的时候都只输出对称唯一的部分。

1 从gaussian中提取积分

以下内容对gaussian09 d.01亲测有效。

由于gaussian会自动调整朝向,会影响积分值,建议加上nosymm,这样也避免因对称性等价的积分不被输出。

计算时加上iop(3/33=1)即可,输出的包括

*** overlap ***:重叠积分

*** kinetic energy ***:动能积分

***** potential energy *****:核吸引势积分(没考虑电子带负电荷,所以数值皆为正)

****** core hamiltonian ******:核哈密顿矩阵。其数值就是动能积分矩阵元减去核吸引势积分矩阵元

multipole matrices:输出三次,ix= 1、ix= 2、ix= 3分别对应x,y,z方向的偶极矩积分矩阵。

gaussian可以通过l316模块双电子积分。用这些关键词即可:iop(3/33=3) extralinks(l316) scf=conventional noraff。同时上一节的单电子积分也会输出。

输出信息诸如:

 intcnt=   4238596 itotal=   4238596 nwiib=    262144 isym2e=0

 i= 44 j= 29 k= 43 l= 10 int=  0.147872132379d-07

 i= 44 j= 43 k= 29 l= 10 int= -0.396977105787d-08

 i= 45 j= 29 k= 43 l= 10 int= -0.403826333553d-07

这里intcnt是实际输出的积分总数。注意并行计算时每次输出的积分顺序可能不同。

gaussian计算并输出双电子积分的值阈值可以用iop(3/27=n)调节,小于10^-n的积分不会被计算和输出。默认n=10,精度已经很高了,若要计算并输出数值更小的可以比如设n=12。

用32bit版本时这种方式输出积分时基函数最多127个。虽然64bit版没这个限制,但基函数多于这个时输出文件会巨大。

2 从gamess-us中提取积分

以下方法对2014-dec版亲测有效。相同基组下(dunning基组不算,因gaussian会自动会去掉其中冗余的gtf)输出的结果和gaussian中的一致。

输出单电子积分和双电子积分都有改代码和不改代码两种做法,结果一样,前者省事,但后者可以自由控制输出格式,而且还能同时用nprint=-5来避免计算时输出一堆翔一样的大量无用的信息。

输出积分时绝对不要用并行运算,否则输出内容会癫狂。

$contrl里写上nprint=3即可在计算时输出基函数间的各种积分矩阵,会看到这样的段落

          ********************

          1 electron integrals

          ********************

下面输出的包括:

overlap matrix:重叠矩阵

bare nucleus hamiltonian integrals (h=t v):核哈密顿矩阵

kinetic energy integrals:动能积分矩阵

程序没直接给出核吸引势积分,这只需自行把核哈密顿矩阵减去动能积分矩阵即可。

在gamess/source/int1.src中—– optional debug printout —–有两处,在第一处的下面加上

      write(*,”(4d18.10)”) s

      write(*,”(/,’ kinetic energy integrals, num:’,i8)”) ll2

      write(*,”(4d18.10)”) t

      write(*,”(/,’ core hamiltonian matrix, num:’,i8)”) ll2

      write(*,”(4d18.10)”) h

      write(*,”(/,’ potential energy integrals, num:’,i8)”) ll2

      write(*,”(4d18.10)”) h-t

      write(*,*)

可见,改过之后可以直接把核吸引势积分输出出来,已经考虑电子是带负电荷。

重新编译此文件并链接成新的可执行文件使之生效,即在gamess-us目录下执行./comp int1;./lked gamess 00。瞬间就完事,用不着把整个gamess-us都重头编译一遍。

之后一般计算时屏幕上就会输出各种单电子积分,上述几种积分矩阵都是对称矩阵,这里依次输出的是下三角部分的元素。如重叠积分:

 overlap integrals, num:      28

  0.1000000000d 01  0.2367039205d 00  0.1000000000d 01  0.0000000000d 00

  0.0000000000d 00  0.1000000000d 01  0.0000000000d 00  0.0000000000d 00

  0.0000000000d 00  0.1000000000d 01  0.2005800301d-16  0.2537204019d-17

  0.0000000000d 00  0.0000000000d 00  0.1000000000d 01  0.5362008406d-01

  0.4729709355d 00  0.0000000000d 00 -0.3205127531d 00  0.2265160217d 00

  0.1000000000d 01  0.5362008406d-01  0.4729709355d 00  0.0000000000d 00

  0.3205127531d 00  0.2265160217d 00  0.2327636225d 00  0.1000000000d 01

这里num后面是输出的元素数目。

一定要用$scf dirscf=.f. $end使用conventional的scf方式,否则不会输出这些积分。

$contrl里写上nprint=4即可在计算时输出双电子积分,会看到如

 ii,jst,kst,lst =  1  1  1  1 nrec =         1 intloc =    1

 ii,jst,kst,lst =  2  1  1  1 nrec =         1 intloc =    2

   1   1   1   1  1.0      4.785065752    2   1   1   1  1.0      0.741380321

   2   2   1   1  1.0      1.118946840    3   3   1   1  1.0      1.115813819

   4   4   1   1  1.0      1.115813819    5   5   1   1  1.0      1.115813819

   2   1   2   1  1.0      0.136873367    3   1   3   1  1.0      0.024477411

 total number of nonzero two-electron integrals =                 228

最后一行就是实际输出的双电子积分数,当前是228个。默认是小于10^-9的积分不被输出,可以在$contrl里用icut=n来把输出阈值调为10^-n,显然n越大输出得越多。

对(ij|kl)积分,有8种等价形式,程序自动把指标调换成i>=j、k>=l、i>=k的形式。

在gamess/source/int2a.src的i4 = locl l下面加上

write(*,”(4i6,d25.16)”) i1,i2,i3,i4,val

重新编译此文件并链接成新的可执行文件使之生效,即在gamess-us目录下执行./comp int2a;./lked gamess 00。

之后一般计算时屏幕上就会输出双电子积分,诸如

     1     1     1     1   0.4785065752035100d 01

     2     1     1     1   0.7413803207944081d 00

     2     2     1     1   0.1118946840438162d 01

     3     3     1     1   0.1115813818563188d 01

     4     4     1     1   0.1115813818563188d 01

     5     5     1     1   0.1115813818563188d 01

如上一节所述,非常小的积分不会输出,阈值可以用icut来设。如果想输出所有积分,直接把int2a.src中if(abs(val).lt.cutoff) go to 200这一行注释掉即可。

如果想输出nprint=4那样的将序号调换为i>=j、k>=l、i>=k后的积分,把前述输出语句加到if (out) call intout(i1,i2,i3,i4,qq4,ijkl_index,val)的上面一行。

原文链接:http://sobereva.com/343

网络摘文,本文作者:15h,如若转载,请注明出处:https://www.15cov.cn/2023/08/27/从gaussian和gamess-us中提取电子积分/

发表评论

邮箱地址不会被公开。 必填项已用*标注

网站地图