#include <energy-opt.h>
Public Member Functions | |
EnergyOptimizer (const std::string &el, bool verbose=true) | |
Constructor. | |
virtual | ~EnergyOptimizer () |
Destructor. | |
void | set_params (const arma::ivec &am, const arma::uvec &nf, const arma::uvec &npar, const arma::uvec &optsh) |
Set parameters. | |
std::string | get_el () const |
Get element. | |
void | toggle_init (bool init) |
Initialize? | |
std::vector< arma::vec > | get_exps (const arma::vec &x) const |
Get exponents. | |
virtual BasisSetLibrary | form_basis (const arma::vec &x) const |
Generate basis set. | |
virtual BasisSetLibrary | form_basis (const std::vector< arma::vec > &exps) const |
Generate basis set. | |
virtual std::vector< double > | calcE (const std::vector< BasisSetLibrary > &baslib)=0 |
Calculate energy for basis set. | |
double | optimize (arma::vec &x, size_t maxiter, double nrthr, double gthr) |
double | optimize_full (arma::vec &x, size_t maxiter, double nrthr, double gthr) |
double | scan (arma::vec &x, double xmin, double xmax, double dx) |
Look for optimal polarization exponent. | |
arma::uvec | sh_vec () const |
Get optimized shell index vector. | |
arma::uvec | idx_vec () const |
Get parameter index vector. | |
arma::uvec | idx_vec (arma::sword am) const |
Get parameter index vector. | |
void | print_info (const arma::vec &x, const std::string &msg) const |
Print parameter info. | |
Protected Member Functions | |
arma::vec | calcG (const arma::vec &x, bool check=true) |
Calculate gradient. Check enables sanity checks. | |
arma::vec | calcG (const arma::vec &x, arma::sword am, bool check=true) |
Calculate gradient. Check enables sanity checks. | |
arma::mat | calcH (const arma::vec &x, bool check=true) |
Calculate Hessian. | |
arma::mat | calcH (const arma::vec &x, arma::sword am, bool check=true) |
Calculate Hessian. | |
arma::vec | pad_vec (const arma::vec &sd) const |
Pad vector to fit into x. | |
arma::vec | pad_vec (const arma::vec &sd, arma::sword am) const |
Pad vector to fit into x. | |
Protected Attributes | |
arma::ivec | sham |
Shell angular momentum. | |
arma::uvec | nf |
Amount of functions on the shells. | |
arma::uvec | npar |
Amount of parameters. | |
arma::uvec | optsh |
Optimize shell? | |
double | fd_h |
Step size in finite difference gradient. | |
double | ls_h |
Step size in line search. | |
size_t | ntr |
Amount of consecutive trials. | |
std::string | el |
Element to optimize. | |
bool | verbose |
Verbose operation? | |
bool | init |
Initialization mode? (Sloppier convergence) | |
This file contains routines for energy-optimization of basis sets.
The algorithms are based on gradient descent methods, and use Legendre polynomials to expand the exponents in a 10-base logarithm basis; see G. Petersson, S. Zhong, J. A. Montgomery and M. J. Frish, "On the optimization of Gaussian basis sets", J. Chem. Phys. 118, 1101 (2003).
double EnergyOptimizer::optimize | ( | arma::vec & | x, |
size_t | maxiter, | ||
double | nrthr, | ||
double | gthr | ||
) |
Run optimization. Returns energy, and optimized parameters in x. Use conjugate gradients till |g|<=nrthr, after which toggle Newton-Raphson procedure. Determine convergence by |g|<gthr
double EnergyOptimizer::optimize_full | ( | arma::vec & | x, |
size_t | maxiter, | ||
double | nrthr, | ||
double | gthr | ||
) |
Run optimization. Returns energy, and optimized parameters in x. Use conjugate gradients till |g|<=nrthr, after which toggle Newton-Raphson procedure. Determine convergence by |g|<gthr