Basis set. More...
#include <basis.h>
Public Member Functions | |
BasisSet () | |
Dummy constructor. | |
BasisSet (size_t Nat, const Settings &set) | |
Construct basis set with Nat atoms, using given settings. | |
~BasisSet () | |
Destructor. | |
BasisSet | density_fitting (double fsam=1.5, int lmaxinc=1) const |
BasisSet | exchange_fitting () const |
BasisSet | decontract (arma::mat &m) const |
Decontract basis set, m gives mapping from old functions to new ones. | |
void | add_nucleus (const nucleus_t &nuc) |
Add nucleus. | |
void | add_shell (size_t nucind, const GaussianShell &sh, bool sort=true) |
Add a shell to a nucleus and sort functions if wanted. | |
void | add_shell (size_t nucind, int am, bool uselm, const std::vector< contr_t > &C, bool sort=true) |
Add a shell to a nucleus and sort functions if wanted. | |
void | add_shells (size_t nucind, ElementBasisSet el, bool sort=true) |
Add all shells to a nucleus and sort functions if wanted. | |
void | sort () |
Sort shells in nuclear order, then by angular momentum, then by exponents. | |
void | check_numbering () |
Check numbering of basis functions. | |
void | update_nuclear_shell_list () |
Update the nuclear list of shells. | |
void | compute_nuclear_distances () |
Compute nuclear distance table. | |
void | form_unique_shellpairs () |
Form list of unique shell pairs. | |
std::vector< shellpair_t > | get_unique_shellpairs () const |
Get list of unique shell pairs. | |
std::vector< eripair_t > | get_eripairs (arma::mat &screen, double thr, double omega=0.0, double alpha=1.0, double beta=0.0, bool verbose=false) const |
Get list of ERI pairs. Screening matrix will also be calculated. | |
void | convert_contractions () |
Convert contractions from normalized primitives to unnormalized primitives. | |
void | convert_contraction (size_t ish) |
Convert contraction on given shell. | |
void | normalize (bool coeffs=true) |
Normalize contractions. If !coeffs, only cartesian factors are calculated. | |
void | coulomb_normalize () |
Normalize contractions in Coulomb norm (for density fitting) | |
void | finalize (bool convert=false, bool normalize=true) |
Do all of the above. | |
double | nuclear_distance (size_t i, size_t j) const |
Get distance of nuclei. | |
arma::mat | nuclear_distances () const |
Get nuclear distances. | |
int | get_am (size_t shind) const |
Get angular momentum of shell. | |
int | get_max_am () const |
Get maximum angular momentum in basis set. | |
size_t | get_max_Ncontr () const |
Get maximum number of contractions. | |
size_t | get_last_ind () const |
Get index of last function, throws an exception if no functions exist. | |
size_t | get_first_ind (size_t shind) const |
Get index of first function on shell. | |
size_t | get_last_ind (size_t shind) const |
Get index of last function on shell. | |
size_t | find_shell_ind (size_t find) const |
Find shell index of basis function. | |
std::vector< GaussianShell > | get_shells () const |
Get shells in basis set. | |
GaussianShell | get_shell (size_t shind) const |
Get ind:th shell. | |
size_t | get_shell_center_ind (size_t shind) const |
Get index of the center of the ind'th shell. | |
coords_t | get_shell_center (size_t shind) const |
Get coordinates of center of ind'th shell. | |
std::vector< contr_t > | get_contr (size_t ind) const |
Get exponential contraction of the ind:th shell. | |
std::vector< contr_t > | get_contr_normalized (size_t ind) const |
Get normalized exponential contraction of the ind:th shell. | |
std::vector< shellf_t > | get_cart (size_t ind) const |
Get the cartesian functions on the ind:th shell. | |
bool | is_lm_default () const |
Are spherical harmonics the default for new shells? | |
bool | lm_in_use (size_t ind) const |
Is shell ind using spherical harmonics? | |
void | set_lm (size_t ind, bool lm) |
Toggle the use of spherical harmonics on shell ind. | |
arma::mat | get_trans (size_t ind) const |
Get transformation matrix. | |
void | fill_sph_transmat () |
Fill spherical harmonics transformation table. | |
size_t | get_Nbf () const |
Get size of basis set. | |
size_t | get_Ncart () const |
Get amount of cartesian functions on shell. | |
size_t | get_Nlm () const |
Get amount of spherical harmonics on shell. | |
size_t | get_Nshells () const |
Get number of shells. | |
size_t | get_Nbf (size_t ind) const |
Get number of basis functions on shell. | |
size_t | get_Ncart (size_t ind) const |
Get number of cartesians on shell. | |
void | compute_shell_ranges (double eps=1e-10) |
std::vector< double > | get_shell_ranges () const |
Get precomputed ranges of shells. | |
std::vector< double > | get_shell_ranges (double eps) const |
Get range of shells with given value of epsilon. | |
std::vector< double > | get_nuclear_distances (size_t inuc) const |
Get distances to other nuclei. | |
size_t | get_Nnuc () const |
Get number of nuclei. | |
nucleus_t | get_nucleus (size_t inuc) const |
Get nucleus. | |
std::vector< nucleus_t > | get_nuclei () const |
Get nuclei. | |
arma::mat | get_nuclear_coords () const |
Get coordinates of all nuclei. | |
void | set_nuclear_coords (const arma::mat &coords) |
Set coordinates of all nuclei. | |
coords_t | get_nuclear_coords (size_t inuc) const |
Get coordinates of nucleus. | |
int | get_Z (size_t inuc) const |
Get charge of nucleus. | |
std::string | get_symbol (size_t inuc) const |
Get symbol of nucleus. | |
std::string | get_symbol_hr (size_t inuc) const |
Get human readable symbol of nucleus (-Bq) | |
std::vector< GaussianShell > | get_funcs (size_t inuc) const |
Get basis functions centered on a given atom. | |
std::vector< size_t > | get_shell_inds (size_t inuc) const |
Get indices of shells centered on a given atom. | |
arma::vec | eval_func (double x, double y, double z) const |
Evaluate functions at (x,y,z) | |
arma::mat | eval_grad (double x, double y, double z) const |
Evaluate gradient at (x,y,z) | |
arma::mat | eval_hess (double x, double y, double z) const |
Evaluate Hessian at (x,y,z) | |
arma::vec | eval_func (size_t ish, double x, double y, double z) const |
Evaluate functions of shell ish at (x,y,z) | |
arma::mat | eval_grad (size_t ish, double x, double y, double z) const |
Evaluate gradients of shell ish at (x,y,z) | |
arma::vec | eval_lapl (size_t ish, double x, double y, double z) const |
Evaluate laplacian of shell ish at (x,y,z) | |
arma::mat | eval_hess (size_t ish, double x, double y, double z) const |
Evaluate Hessian of shell ish at (x,y,z) | |
arma::mat | eval_laplgrad (size_t ish, double x, double y, double z) const |
Evaluate gradient of laplacian of shell ish at (x,y,z) | |
void | print (bool verbose=false) const |
Print out basis set. | |
arma::mat | cart_to_sph_trans () const |
Calculate transformation matrix from cartesians to spherical harmonics. | |
arma::mat | sph_to_cart_trans () const |
Calculate transfomration matrix from spherical harmonics to cartesians. | |
arma::mat | overlap () const |
Calculate overlap matrix. | |
arma::mat | coulomb_overlap () const |
Calculate overlap matrix in Coulomb metric. | |
arma::mat | overlap (const BasisSet &rhs) const |
Calculate overlap with another basis set. | |
arma::mat | coulomb_overlap (const BasisSet &rhs) const |
Calculate overlap in Coulomb metric with another basis set. | |
arma::mat | kinetic () const |
Calculate kinetic energy matrix. | |
arma::mat | nuclear () const |
Calculate nuclear repulsion matrix. | |
arma::mat | potential (coords_t r) const |
Calculate electric potential matrix. | |
arma::mat | eri_screening (double omega=0.0, double alpha=1.0, double beta=0.0) const |
Calculate ERI screening matrix. | |
arma::vec | nuclear_pulay (const arma::mat &P) const |
Calculate nuclear Pulay forces. | |
arma::vec | nuclear_der (const arma::mat &P) const |
Calculate nuclear Hellman-Feynman force. | |
arma::vec | kinetic_pulay (const arma::mat &P) const |
Calculate kinetic Pulay forces. | |
arma::vec | overlap_der (const arma::mat &W) const |
Calculate overlap derivative force. | |
arma::vec | nuclear_force () const |
Calculate nuclear repulsion force. | |
std::vector< arma::mat > | moment (int mom, double x=0.0, double y=0.0, double z=0.0) const |
Compute moment integral around (x,y,z) | |
arma::vec | integral () const |
Compute integrals of basis functions (used in xc-fitting) | |
int | Ztot () const |
Calculate nuclear charge. | |
double | Enuc () const |
Calculate nuclear repulsion energy. | |
void | projectMOs (const BasisSet &old, const arma::colvec &oldE, const arma::mat &oldMOs, arma::colvec &E, arma::mat &MOs, size_t nocc) const |
Project MOs to new basis set. | |
void | projectOMOs (const BasisSet &old, const arma::cx_mat &oldOMOs, arma::cx_mat &OMOs, size_t nocc) const |
Translate OMOs to new geometry, assuming same basis set. Virtuals are discarded and regenerated. | |
bool | operator== (const BasisSet &rhs) const |
Are the basis sets the same? | |
std::vector< std::vector < size_t > > | find_identical_shells () const |
Find "identical" shells in basis set. | |
Private Member Functions | |
bool | same_geometry (const BasisSet &rhs) const |
Check for same geometry. | |
bool | same_shells (const BasisSet &rhs) const |
Check for same shells. | |
Private Attributes | |
std::vector< nucleus_t > | nuclei |
Nuclei. | |
std::vector< GaussianShell > | shells |
Basis functions. | |
bool | uselm |
Use spherical harmonics by default as basis? | |
arma::mat | nucleardist |
Internuclear distances. | |
std::vector< shellpair_t > | shellpairs |
List of unique shell pairs. | |
std::vector< double > | shell_ranges |
Ranges of shells. | |
Basis set.
Class for a Gaussian basis set.
This class contains the data structures necessary for Gaussian basis sets, and functions for evaluating integrals over them.
void BasisSet::compute_shell_ranges | ( | double | eps = 1e-10 | ) |
Get range of shell (distance at which functions have dropped below epsilon)
The default parameter is from the reference R. E. Stratmann, G. E. Scuseria, and M. J. Frisch, "Achieving linear scaling in exchage-correlation density functional quadratures", Chem. Phys. Lett. 257 (1993), pp. 213-223.
BasisSet BasisSet::density_fitting | ( | double | fsam = 1.5 , |
int | lmaxinc = 1 |
||
) | 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.
BasisSet BasisSet::exchange_fitting | ( | ) | const |
Generate Coulomb and exchange fitting basis
The procedure has been documented in the article
O. Vahtras, J. Almlöf and M. W. Feyereisen, "Integral approximations for LCAO-SCF calculations", Chem. Phys. Lett. 213, p. 514 - 518 (1993).