ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
casida_grid.h
1 /*
2  * This source code is part of
3  *
4  * E R K A L E
5  * -
6  * DFT from Hel
7  *
8  * Written by Susi Lehtola, 2010-2011
9  * Copyright (c) 2010-2011, Susi Lehtola
10  *
11  * This program is free software; you can redistribute it and/or
12  * modify it under the terms of the GNU General Public License
13  * as published by the Free Software Foundation; either version 2
14  * of the License, or (at your option) any later version.
15  */
16 
22 #ifndef ERKALE_CASIDAGRID
23 #define ERKALE_CASIDAGRID
24 
25 #include "casida.h"
26 #include "../dftgrid.h"
27 
29 class CasidaShell : public AngularGrid {
31  std::vector< std::vector< std::vector<double> > > orbs;
32 
34  std::vector<double> fx;
36  std::vector<double> fc;
37 
38  public:
40  CasidaShell(bool lobatto=false);
42  ~CasidaShell();
43 
45  void compute_orbs(const std::vector<arma::mat> & C);
46 
48  void eval_fxc(int x_func, int c_func);
49 
51  void Kxc(const std::vector< std::vector<states_pair_t> > & pairs, arma::mat & K) const;
52 
53  void free();
54 };
55 
57 class CasidaGrid {
59  std::vector<CasidaShell> wrk;
61  std::vector<angshell_t> grids;
62 
64  const BasisSet * basp;
66  bool verbose;
69 
71  void construct(const std::vector<arma::mat> & P, double tol, int x_func, int c_func);
73  void prune_shells();
74 
75  public:
76  CasidaGrid(const BasisSet * bas, bool verbose=false, bool lobatto=false);
77  ~CasidaGrid();
78 
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);
82  void print_grid() const;
83 };
84 
85 #endif
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