Program Reference of libcint

libcint is an open source library for computation of Gaussian type integrals. It provides C/Fortran function to evaluate one-electron / two-electron integrals with Cartesian / real-spheric / spinor Gaussian type functions.

Features

  • Various GTO type:
    • Cartesian GTO: 6d, 10f, 15g Gaussian type functions.
    • Real-spheric GTO: 5d, 7f, 9g Gaussian type functions.
    • Spinor GTO: J-adapted spinor Gaussian functions.
  • One electron integrals.
    • Besides the regular kinetic, or nuclear attraction etc., nuclear attraction due to Gaussian nuclear model are supported.
  • Two electron integrals (value < 1e–20 are neglected) include
    • Coulomb repulsion
    • Gaunt interaction
  • Common lisp script to generate C code for new integrals.

  • Thread safe.

  • Uniform API for all kind of integrals.
    • one electron integrals
        not0 = fn1e_name(double *buf, int *atm, int natm, int *bas, int nbas, double *env);
- two electron integrals
        not0 = fn2e_name(double *buf, int *atm, int natm, int *bas, int nbas, double *env, NULL);
- the return boolean (not0) gives the summary whether the integrals
  are completely 0.
  • Minimal overhead of initialization

    • Pre-comuptation is not required. Only basic info (see previous API) of basis function need to be initialized, within a plain integer or double precision array, once. (For 2e integral, there is an optional argument called optimizer which can increase the performance by ~10%. Using NULL to swith off the optimizer does not affect the value of integral.)
  • Minimal dependence on external library.

    • BLAS is the only library needed. Normally, the performance difference due to various BLAS implementations is less than 1%.
  • Small memory usage.

    • Very few intermediate data are allocated. ~80% of the memory are used for holding the whole contracted Cartesion integrals, which typically should be less than 1 Mega bytes.

Getting libcint

The latest version is available from:

git clone http://github.com/sunqm/libcint.git

It’s very convenient to tryout Libcint with PySCF, which is a python module for quantum chemistry program

http://github.com/sunqm/pyscf.git

Installation

  • Prerequisites
    • BLAS library
    • Python version 2.5 or higher (optional, for make check)
    • clisp / SBCL (optional, for common lisp script)
  • Build libcint
    mkdir build; cd build
cmake .. <options>
make
make install

Known problems

  • On 64-bit system, “make test” stop with error:

    MKL FATAL ERROR: Cannot load libmkl_avx.so or libmkl_def.so.

This problem is caused by the conflict between Python and MKL library. It can be fixed by adding -lmkl_avx or -lmkl_mc -lmkl_def to MKL link flags to replace the default blas link flags. Be careful with the order of -lmkl_mc and -lmkl_def.

  • For basic ERIs, the code can handle highest angular momentum up to 6 (present Rys-roots functions might be numerically unstable for l=5,6). But it has to be reduced to 5 or less for derivative or high order ERI. Depending on the derivative order, reduce 1 highest angular momentum every 4 derivative order.

  • Using SSE3 instructions can increase the performance 5 ~ 50%. Please refer to qcint library (under GPL v3 license)

    https://github.com/sunqm/qcint.git
    

Interface

C routine

dim = CINTgto_cart(bas_id, bas);
dim = CINTgto_spheric(bas_id, bas);
dim = CINTgto_spinor(bas_id, bas);
f1e(buf, shls, atm, natm, bas, nbas, env);
f2e(buf, shls, atm, natm, bas, nbas, env, opt);
f2e_optimizer(&opt, atm, nbat, bas, nbas, env);
CINTdel_optimizer(&opt);
  • buf: column-major double precision array.
    • for 1e integrals of shells (i,j), data are stored as [i1j1 i2j1 ... ]
    • for 2e integrals of shells (i,j|k,l), data are stored as
      [i1j1k1l1 i2j1k1l1 ... i1j2k1l1 ... i1j1k2l1 ... ]
    • complex data are stored as two double elements, first is real, followed by imaginary, e.g. [Re Im Re Im …]
  • shls: 0-based basis/shell indices.
    • int[2] for 1e integrals
    • int[4] for 2e integrals
  • atm: int[natm*6], list of atoms. For ith atom, the 6 slots of atm[i] are
    • atm[i*6+0] nuclear charge of atom i
    • atm[i*6+1] env offset to save coordinates (env[atm[i*6+1]], env[atm[i*6+1]+1], env[atm[i*6+1]+2]) are (x,y,z)
    • atm[i*6+2] nuclear model of atom i, = 2 indicates gaussian nuclear model ρ(r)=Z(ζπ)3/2exp(-ζr2)
    • atm[i*6+3] env offset to save the nuclear charge distribution parameter ζ
    • atm[i*6+4] unused
    • atm[i*6+5] unused
  • natm: int, number of atoms, natm has no effect except nuclear attraction integrals
  • bas: int[nbas*8], list of basis. For ith basis, the 8 slots of bas[i] are
    • bas[i*8+0] 0-based index of corresponding atom
    • bas[i*8+1] angular momentum
    • bas[i*8+2] number of primitive GTO in basis i
    • bas[i*8+3] number of contracted GTO in basis i
    • bas[i*8+4] kappa for spinor GTO.
      < 0 the basis ~ j = l + 1/2.
      > 0 the basis ~ j = l - 1/2.
      = 0 the basis includes both j = l + 1/2 and j = l - 1/2
    • bas[i*8+5] env offset to save exponents of primitive GTOs. e.g. 10 exponents env[bas[i*8+5]]env[bas[i*8+5]+9]
    • bas[i*8+6] env offset to save column-major contraction coefficients. e.g. 10 primitive -> 5 contraction needs a 10×5 array
env[bas[i*8+6]  ] | env[bas[i*8+6]+10] |     | env[bas[i*8+6]+40]
env[bas[i*8+6]+1] | env[bas[i*8+6]+11] |     | env[bas[i*8+6]+41]
 .                |  .                 | ... |  .                
 .                |  .                 |     |  .                
env[bas[i*8+6]+9] | env[bas[i*8+6]+19] |     | env[bas[i*8+6]+49]
- `bas[i*8+7]` unused
  • nbas: int, number of bases, nbas has no effect
  • env: double[], save the value of coordinates, exponents, contraction coefficients

  • struct CINTOpt *opt: so called “optimizer”, it needs to be intialized

    CINTOpt *opt = NULL; intname_optimizer(&opt, atm, natm, bas, nbas, env);

every integral type has its own optimizer with the suffix _optimizer in its name, e.g. the optimizer for cint2e_sph is cint2e_sph_opimizer. “optimizer” is an optional argument for the integrals. It can roughly speed the integration by ~30% without affecting the value of integrals. However, it consumes about 5num_primitiveGTO*2 extra memory. If no optimizer is wanted, set it to NULL.

optimizer needs to be released after using.

CINTdel_optimizer(&opt);
  • if the return value equals 0, every element of the integral is 0

  • short example

#include "cint.h"
...
CINTOpt *opt = NULL;
cint2e_sph_optimizer(&opt, atm, natm, bas, nbas, env);
for (i = 0; i < nbas; i++) {
shls[0] = i;
di = CINTcgto_spheric(i, bas);
...
for (l = 0; l < nbas; l++) {
shls[3] = l;
dl = CINTcgto_spheric(l, bas);
buf = malloc(sizeof(double) * di * dj * dk * dl);
cint2e_cart(buf, shls, atm, natm, bas, nbas, env, opt);
free(buf);
}
}
CINTdel_optimizer(&opt);

Fortran routine

dim = CINTgto_cart(bas_id, bas)
dim = CINTgto_spheric(bas_id, bas)
dim = CINTgto_spinor(bas_id, bas)
call f1e(buf, shls, atm, natm, bas, nbas, env)
call f2e(buf, shls, atm, natm, bas, nbas, env, opt)
call f2e_optimizer(opt, atm, nbat, bas, nbas, env)
call CINTdel_optimizer(opt)
  • atm and bas are 2D integer array
    • atm(1:6,i) is the (charge, offset_coord, nuclear_model, unused, unused, unused) of the ith atom
    • bas(1:8,i) is the (atom_index, angular, num_primitive_GTO, num_contract_GTO, kappa, offset_exponent, offset_coeff, unused) of the ith basis
  • parameters are the same to the C function. Note that those offsets atm(2,i) bas(6,i) bas(7,i) are 0-based.
  • buf is 2D/4D double precision/double complex array
  • opt: an integer(8) to hold the address of so called “optimizer”, it needs to be intialized by

    integer(8) opt call f2e_optimizer(opt, atm, natm, bas, nbas, env)

The optimizier can be banned by setting the “optimizier” to 0_8

call f2e(buf, atm, natm, bas, nbas, env, 0_8)

To release optimizer, execute

call CINTdel_optimizer(opt);
  • short example
...
integer,external CINTcgto_spheric
integer(8) opt
call cint2e_sph_optimizer(opt, atm, natm, bas, nbas, env)
do i = 1, nbas
shls(1) = i - 1
di = CINTcgto_spheric(i-1, bas)
...
do l = 1, nbas
shls(4) = l - 1
dl = CINTcgto_spheric(l-1, bas)
allocate(buf(di,dj,dk,dl))
call cint2e_sph(buf, shls, atm, natm, bas, nbas, env, opt)
deallocate(buf)
end do
end do
call CINTdel_optimizer(opt)

Data ordering

  • for Cartesian GTO, the output data in buf are sorted as
s shellp shelld shellf shell
sxxxxxx
syxyxxy
zxzxxz
xyyxyy
yyzxyz
zzzxzz
yyy
yyz
yzz
zzz
  • for real spheric GTO, the output data in buf are sorted as
s shellp shelld shellf shell
sxxyy(3x2-y2)
syyzxyz
zz2yz2
xxzz3
yx2-y2xz2
zz(x2-y2)
x(x2-3y2)
  • for spinor GTO, the output data in buf correspond to
kappa=0,p shellkappa=1,p shellkappa=0,d shell
p1/2(-1/2)p1/2(-1/2)d3/2(-3/2)
p1/2(1/2)p1/2(1/2)d3/2(-1/2)
p3/2(-3/2)p1/2(-1/2)d3/2(1/2)
p3/2(-1/2)p1/2(1/2)d3/2(3/2)
p3/2(1/2)p1/2(-1/2)d5/2(-5/2)
p3/2(3/2)p1/2(1/2)d5/2(-3/2)
p1/2(-1/2)d5/2(-1/2)
p1/2(1/2)d3/2(-3/2)
p3/2(-3/2)d3/2(-1/2)
p3/2(-1/2)

Tensor

Integrals like Gradients have more than one components. The output array is ordered in Fortran-contiguous. The tensor component takes the biggest strides.

  • 3-component tensor
    • X buf(:,0)
    • Y buf(:,1)
    • Z buf(:,2)
  • 9-component tensor
    • XX buf(:,0)
    • XY buf(:,1)
    • XZ buf(:,2)
    • YX buf(:,3)
    • YY buf(:,4)
    • YZ buf(:,5)
    • ZX buf(:,6)
    • ZY buf(:,7)
    • ZZ buf(:,8)

Built-in function list

  • Cartesian GTO integrals
    • CINTcgto_cart
    • cint1e_ovlp_cart ij
    • cint1e_nuc_cart iVnucj
    • cint1e_kin_cart .5ippj
    • cint1e_ia01p_cart irr3×j
    • cint1e_irixp_cart i(r-Ri)×j
    • cint1e_ircxp_cart i(r-Ro)×j
    • cint1e_iking_cart 0.5ippiUgj
    • cint1e_iovlpg_cart iiUgj
    • cint1e_inucg_cart iiVnucUgj
    • cint1e_ipovlp_cart ij
    • cint1e_ipkin_cart 0.5ippj
    • cint1e_ipnuc_cart iVnucj
    • cint1e_iprinv_cart ir-1j
    • cint1e_rinv_cart ir-1j
    • cint2e_cart (ijkl)
    • cint2e_ig1_cart i(iUgjkl)
    • cint2e_ip1_cart (ijkl)
  • Spheric GTO integrals
    • CINTcgto_spheric
    • cint1e_ovlp_sph ij
    • cint1e_nuc_sph iVnucj
    • cint1e_kin_sph 0.5ippj
    • cint1e_ia01p_sph irr3×j
    • cint1e_irixp_sph i(rc-Ri)×j
    • cint1e_ircxp_sph i(rc-Ro)×j
    • cint1e_iking_sph 0.5ippiUgj
    • cint1e_iovlpg_sph iiUgj
    • cint1e_inucg_sph iiVnucUgj
    • cint1e_ipovlp_sph ij
    • cint1e_ipkin_sph 0.5ippj
    • cint1e_ipnuc_sph iVnucj
    • cint1e_iprinv_sph ir-1j
    • cint1e_rinv_sph ir-1j
    • cint2e_sph (ijkl)
    • cint2e_ig1_sph i(iUgjkl)
    • cint2e_ip1_sph (ijkl)
  • Spinor GTO integrals
    • CINTcgto_spinor
    • cint1e_ovlp ij
    • cint1e_nuc iVnucj
    • cint1e_nucg iVnucUgj
    • cint1e_srsr σriσrj
    • cint1e_sr σrij
    • cint1e_srsp σriσpj
    • cint1e_spsp σpiσpj
    • cint1e_sp σpij
    • cint1e_spspsp σpiσpσpj
    • cint1e_spnuc σpiVnucj
    • cint1e_spnucsp σpiVnucσpj
    • cint1e_srnucsr σriVnucσrj
    • cint1e_sa10sa01 0.5σ×rciσ×rr3j
    • cint1e_ovlpg iUgj
    • cint1e_sa10sp 0.5rc×σiσpj
    • cint1e_sa10nucsp 0.5rc×σiVnucσpj
    • cint1e_sa01sp irr3×σσpj
    • cint1e_spgsp Ugσpiσpj
    • cint1e_spgnucsp UgσpiVnucσpj
    • cint1e_spgsa01 Ugσpirr3×σj
    • cint1e_ipovlp ij
    • cint1e_ipkin 0.5ippj
    • cint1e_ipnuc iVnucj
    • cint1e_iprinv ir-1j
    • cint1e_ipspnucsp σpiVnucσpj
    • cint1e_ipsprinvsp σpir-1σpj
    • cint2e (ijkl)
    • cint2e_spsp1 (σpiσpjkl)
    • cint2e_spsp1spsp2 (σpiσpjσpkσpl)
    • cint2e_srsr1 (σriσrjkl)
    • cint2e_srsr1srsr2 (σriσrjσrkσrl)
    • cint2e_sa10sp1 0.5(rc×σiσpjkl)
    • cint2e_sa10sp1spsp2 0.5(rc×σiσpjσpkσpl)
    • cint2e_g1 (iUgjkl)
    • cint2e_spgsp1 (σpiUgσpjkl)
    • cint2e_g1spsp2 (iUgjσpkσpl)
    • cint2e_spgsp1spsp2 (σpiUgσpjσpkσpl)
    • cint2e_ip1 (ijkl)
    • cint2e_ipspsp1 (σpiσpjkl)
    • cint2e_ip1spsp2 (ijσpkσpl)
    • cint2e_ipspsp1spsp2 (σpiσpjσpkσpl)
    • cint2e_ssp1ssp2 (iσσpjkσσpl)