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
Bug report
Qiming Sun osirpt.sun@gmail.com
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 …]
- for 1e integrals of shells (i,j), data are stored as
- 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 iatm[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 modelatm[i*6+3]
env offset to save the nuclear charge distribution parameteratm[i*6+4]
unusedatm[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 atombas[i*8+1]
angular momentumbas[i*8+2]
number of primitive GTO in basis ibas[i*8+3]
number of contracted GTO in basis ibas[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/2bas[i*8+5]
env offset to save exponents of primitive GTOs. e.g. 10 exponentsenv[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 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 shell | p shell | d shell | f shell |
---|---|---|---|
s | |||
s | |||
… | |||
… | … | ||
… |
- for real spheric GTO, the output data in buf are sorted as
s shell | p shell | d shell | f shell |
---|---|---|---|
s | |||
s | |||
… | |||
… | |||
… | |||
… |
- for spinor GTO, the output data in buf correspond to
… | kappa=0,p shell | kappa=1,p shell | kappa=0,d shell | … |
---|---|---|---|---|
… | ||||
… | ||||
… |
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)
- X
- 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)
- XX
Built-in function list
- Cartesian GTO integrals
CINTcgto_cart
cint1e_ovlp_cart
cint1e_nuc_cart
cint1e_kin_cart
cint1e_ia01p_cart
cint1e_irixp_cart
cint1e_ircxp_cart
cint1e_iking_cart
cint1e_iovlpg_cart
cint1e_inucg_cart
cint1e_ipovlp_cart
cint1e_ipkin_cart
cint1e_ipnuc_cart
cint1e_iprinv_cart
cint1e_rinv_cart
cint2e_cart
cint2e_ig1_cart
cint2e_ip1_cart
- Spheric GTO integrals
CINTcgto_spheric
cint1e_ovlp_sph
cint1e_nuc_sph
cint1e_kin_sph
cint1e_ia01p_sph
cint1e_irixp_sph
cint1e_ircxp_sph
cint1e_iking_sph
cint1e_iovlpg_sph
cint1e_inucg_sph
cint1e_ipovlp_sph
cint1e_ipkin_sph
cint1e_ipnuc_sph
cint1e_iprinv_sph
cint1e_rinv_sph
cint2e_sph
cint2e_ig1_sph
cint2e_ip1_sph
- Spinor GTO integrals
CINTcgto_spinor
cint1e_ovlp
cint1e_nuc
cint1e_nucg
cint1e_srsr
cint1e_sr
cint1e_srsp
cint1e_spsp
cint1e_sp
cint1e_spspsp
cint1e_spnuc
cint1e_spnucsp
cint1e_srnucsr
cint1e_sa10sa01
cint1e_ovlpg
cint1e_sa10sp
cint1e_sa10nucsp
cint1e_sa01sp
cint1e_spgsp
cint1e_spgnucsp
cint1e_spgsa01
cint1e_ipovlp
cint1e_ipkin
cint1e_ipnuc
cint1e_iprinv
cint1e_ipspnucsp
cint1e_ipsprinvsp
cint2e
cint2e_spsp1
cint2e_spsp1spsp2
cint2e_srsr1
cint2e_srsr1srsr2
cint2e_sa10sp1
cint2e_sa10sp1spsp2
cint2e_g1
cint2e_spgsp1
cint2e_g1spsp2
cint2e_spgsp1spsp2
cint2e_ip1
cint2e_ipspsp1
cint2e_ip1spsp2
cint2e_ipspsp1spsp2
cint2e_ssp1ssp2