Basis set for an element. More...
#include <basislibrary.h>
Public Member Functions | |
ElementBasisSet () | |
Dummy constructor. | |
ElementBasisSet (std::string sym, size_t number=0) | |
Constructor. | |
~ElementBasisSet () | |
Destructor. | |
void | add_function (FunctionShell f) |
Add a shell of functions to the basis. | |
void | sort () |
Sort the shells in decreasing angular momentum. | |
void | print () const |
Print out basis set. | |
std::string | get_symbol () const |
Get the symbol of the element. | |
size_t | get_number () const |
Get the number. | |
void | set_number (size_t num) |
Set the number. | |
bool | operator< (const ElementBasisSet &rhs) const |
Comparison operator for sorting. | |
std::vector< FunctionShell > | get_shells () const |
Get the shells. | |
std::vector< FunctionShell > | get_shells (int am) const |
Get the shells for am value. | |
void | get_primitives (arma::vec &exps, arma::mat &coeffs, int am) const |
Get exponents and contraction coefficients of angular momentum shell am. | |
void | get_primitives (arma::vec &zfree, arma::vec &zgen, arma::mat &cgen, int am) const |
Get free primitives and generally contracted part. | |
void | orthonormalize () |
Orthonormalize generally contracted functions. NB! This can screw up your computational efficiency! | |
void | P_orthogonalize (double cutoff=1e-8, double Cortho=1e-4) |
void | prune (double cutoff=0.9, bool coulomb=false) |
Prune primitives with large overlap, replacing them with the geometric average. | |
void | merge (double cutoff=0.9, bool verbose=true, bool coulomb=false) |
Merge primitives with large overlap (also decontracts basis). This routine is intended to merging core-correlation functions. | |
int | get_max_am () const |
Get maximum angular momentum used in the shells. | |
int | get_am (size_t ind) const |
Get angular momentum of i:th shell. | |
int | get_Nbf () const |
Get amount of functions. | |
void | normalize () |
Normalize coefficients. | |
void | decontract () |
Decontract set. | |
ElementBasisSet | density_fitting (int lmaxinc, double fsam) const |
ElementBasisSet | product_set (int lmaxinc, double fsam) const |
ElementBasisSet | cholesky_set (double thr, int maxam, double ovlthr) const |
Form compact Cholesky set. | |
void | augment (int naug) |
Augment the basis. | |
Private Attributes | |
std::string | symbol |
Symbol of element. | |
size_t | number |
Atom id for which the basis is for (0 for all atoms) | |
std::vector< FunctionShell > | bf |
List of shells. | |
Friends | |
class | BasisSetLibrary |
Basis set for an element.
This class defines a basis set for an element, used by the BasisSetLibrary class.
ElementBasisSet ElementBasisSet::density_fitting | ( | int | lmaxinc, |
double | fsam | ||
) | const |
Generate density fitting basis
The procedure has been documented in the article
R. Yang, A. P. Rendell and M. J. Frisch, "Automatically generated Coulomb fitting basis sets: Design and accuracy for systems containing H to Kr", J. Chem. Phys. 127 (2007), 074102.
void ElementBasisSet::P_orthogonalize | ( | double | cutoff = 1e-8 , |
double | Cortho = 1e-4 |
||
) |
P-orthogonalization of basis set. A refinement on the Davidson refinement (a.k.a. "Davidson rotation").
The procedure is based on the article
F. Jensen, "Unifying General and Segmented Contracted Basis Sets. Segmented Polarization Consistent Basis Sets", J. Chem. Theory Comput. 10, 1074 (2014).
ElementBasisSet ElementBasisSet::product_set | ( | int | lmaxinc, |
double | fsam | ||
) | const |
Generate product basis set.
Brute-force generation of orbital products, followed by merging product exponents that deviate by fsam at maximum.