ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
Public Member Functions | Private Member Functions | Private Attributes | List of all members
DensityFit Class Reference

Density fitting routines. More...

#include <density_fitting.h>

Public Member Functions

 DensityFit ()
 Constructor.
 
 ~DensityFit ()
 Destructor.
 
size_t fill (const BasisSet &orbbas, const BasisSet &auxbas, bool direct, double erithr, double linthr, bool hf=false)
 
size_t idx (size_t imu, size_t inu) const
 Compute index in integral table.
 
size_t memory_estimate (const BasisSet &orbbas, const BasisSet &auxbas, bool direct) const
 Compute estimate of necessary memory.
 
arma::vec compute_expansion (const arma::mat &P) const
 Compute expansion coefficients c.
 
std::vector< arma::vec > compute_expansion (const std::vector< arma::mat > &P) const
 Compute expansion coefficients c.
 
arma::mat calcJ (const arma::mat &P) const
 Get Coulomb matrix from P.
 
std::vector< arma::mat > calcJ (const std::vector< arma::mat > &P) const
 Get Coulomb matrix from P.
 
arma::vec forceJ (const arma::mat &P)
 Calculate force from P.
 
arma::mat calcK (const arma::mat &C, const std::vector< double > &occs, size_t memlimit) const
 Get exchange matrix from orbitals with occupation numbers occs.
 
arma::cx_mat calcK (const arma::cx_mat &C, const std::vector< double > &occs, size_t memlimit) const
 Get exchange matrix from orbitals with occupation numbers occs.
 
size_t get_Norb () const
 Get the number of orbital functions.
 
size_t get_Naux () const
 Get the number of auxiliary functions.
 
double get_a_munu (size_t ia, size_t imu, size_t inu) const
 Get the three-electron integral.
 
arma::mat get_ab_inv () const
 Get ab_inv.
 
void B_matrix (arma::mat &B) const
 Get B matrix (must have HF enabled)
 

Private Member Functions

void form_screening ()
 Form screening matrix.
 

Private Attributes

size_t Nbf
 Amount of orbital basis functions.
 
size_t Naux
 Amount of auxiliary basis functions.
 
bool direct
 Direct calculation? (Compute three-center integrals on-the-fly)
 
bool hf
 Hartree-Fock calculation?
 
size_t Nnuc
 Amount of nuclei.
 
int maxam
 Maximum angular momentum.
 
int maxcontr
 Maximum contractions.
 
std::vector< GaussianShellorbshells
 Orbital shells.
 
int maxorbam
 
size_t maxorbcontr
 
std::vector< GaussianShellauxshells
 Density fitting shells.
 
int maxauxam
 
size_t maxauxcontr
 
GaussianShell dummy
 Dummy shell.
 
size_t dummyind
 Index of dummy function.
 
std::vector< eripair_torbpairs
 List of unique orbital shell pairs.
 
std::vector< size_t > iidx
 Index helper.
 
arma::mat a_munu
 Integrals $ ( \alpha | \mu \nu) $.
 
arma::mat ab
 $ ( \alpha | \beta) $
 
arma::mat ab_inv
 $ ( \alpha | \beta)^-1 $
 
arma::mat ab_invh
 $ ( \alpha | \beta)^-1/2 $
 

Detailed Description

Density fitting routines.

Density fitting / RI routines.

This class contains density fitting and resolution of the identity routines used for the approximate calculation of the Coulomb and exchange operators J and K.

The RI-JK implementation is based on the procedure described in

F. Weigend, "A fully direct RI-HF algorithm: Implementation, optimised auxiliary basis sets, demonstration of accuracy and efficiency", Phys. Chem. Chem. Phys. 4, 4285 (2002).

If only RI-J is necessary, then the procedure described in

K. Eichkorn, O. Treutler, H. Öhm, M. Häser and R. Ahlrichs, "Auxiliary basis sets to approximate Coulomb potentials", Chem. Phys. Lett. 240 (1995), 283-290.

is used.

Author
Susi Lehtola
Date
2012/08/22 23:53

Member Function Documentation

size_t DensityFit::fill ( const BasisSet orbbas,
const BasisSet auxbas,
bool  direct,
double  erithr,
double  linthr,
bool  hf = false 
)

Compute integrals, use given linear dependency threshold. The HF flag here controls formation of (a|b)^{-1/2} and (a|b)^{-1}; the HF routine should be more tolerant of linear dependencies in the basis. Returns amount of significant orbital shell pairs.


The documentation for this class was generated from the following files: