22 #ifndef ERKALE_CASIDAGRID
23 #define ERKALE_CASIDAGRID
26 #include "../dftgrid.h"
31 std::vector< std::vector< std::vector<double> > >
orbs;
34 std::vector<double>
fx;
36 std::vector<double>
fc;
48 void eval_fxc(
int x_func,
int c_func);
51 void Kxc(
const std::vector< std::vector<states_pair_t> > & pairs, arma::mat & K)
const;
59 std::vector<CasidaShell>
wrk;
71 void construct(
const std::vector<arma::mat> & P,
double tol,
int x_func,
int c_func);
80 void Kxc(
const std::vector<arma::mat> & P,
double tol,
int x_func,
int c_func,
const std::vector<arma::mat> & C,
const std::vector < std::vector<states_pair_t> > & pairs, arma::mat & Kx);
std::vector< double > fc
Values of the correlation part.
Definition: casida_grid.h:36
std::vector< double > fx
Values of the exchange part of .
Definition: casida_grid.h:34
std::vector< angshell_t > grids
Angular grids.
Definition: casida_grid.h:61
void eval_fxc(int x_func, int c_func)
Evaluate fx and fc.
Definition: casida_grid.cpp:66
void construct(const std::vector< arma::mat > &P, double tol, int x_func, int c_func)
Construct grid.
Definition: casida_grid.cpp:231
bool verbose
Verbose operation?
Definition: casida_grid.h:66
void Kxc(const std::vector< std::vector< states_pair_t > > &pairs, arma::mat &K) const
Evaluate Kxc.
Definition: casida_grid.cpp:142
Perform integral over atomic grid.
Definition: casida_grid.h:29
const BasisSet * basp
Basis set.
Definition: casida_grid.h:64
void compute_orbs(const std::vector< arma::mat > &C)
Evaluate values of orbitals at grid points.
Definition: casida_grid.cpp:32
Basis set.
Definition: basis.h:187
~CasidaShell()
Destructor.
Definition: casida_grid.cpp:29
void Kxc(const std::vector< arma::mat > &P, double tol, int x_func, int c_func, const std::vector< arma::mat > &C, const std::vector< std::vector< states_pair_t > > &pairs, arma::mat &Kx)
Evaluate Kxc.
Definition: casida_grid.cpp:311
void prune_shells()
Prune shells without points.
Definition: casida_grid.cpp:305
Angular integration grid on a radial shell of an atom.
Definition: dftgrid.h:164
CasidaShell(bool lobatto=false)
Constructor. Need to set tolerance as well before using constructor!
Definition: casida_grid.cpp:26
Perform integral over molecular grid.
Definition: casida_grid.h:57
std::vector< std::vector< std::vector< double > > > orbs
Stack of values of orbitals at grid points: orbs[nspin][ngrid][norb].
Definition: casida_grid.h:31
void print_grid() const
Print out grid composition.
Definition: casida_grid.cpp:379
std::vector< CasidaShell > wrk
Work grids.
Definition: casida_grid.h:59
bool use_lobatto
Use Lobatto quadrature?
Definition: casida_grid.h:68