ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
scf.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 
17 
18 
19 #ifndef ERKALE_SCF
20 #define ERKALE_SCF
21 
22 #include "global.h"
23 
24 #include <armadillo>
25 #include <vector>
26 class DFTGrid;
27 #include "eritable.h"
28 #include "eriscreen.h"
29 #include "erichol.h"
30 #include "density_fitting.h"
31 class Settings;
32 
33 class Checkpoint;
34 
48 typedef struct {
51  double deltaEmax;
53  double deltaPmax;
55  double deltaPrms;
57 
59 typedef struct {
61  int x_func;
63  int c_func;
64 
66  bool adaptive;
68  double gridtol;
70  int nrad;
72  int lmax;
74  bool lobatto;
75 
76  // Non-local part?
77  bool nl;
78  double vv10_b, vv10_C;
79  int nlnrad, nllmax;
80 } dft_t;
81 
83 typedef struct {
85  double Ecoul;
87  double Ekin;
89  double Enuca;
91  double Exc;
93  double Eone;
94 
96  double Eel;
98  double Enucr;
100  double Enl;
102  double Esic;
103 
105  double E;
106 } energy_t;
107 
108 
110 typedef struct {
112  arma::mat C;
114  arma::vec E;
116  arma::mat H;
118  arma::mat P;
119 
121  arma::mat J;
123  arma::mat K;
125  arma::mat XC;
126 
128  arma::cx_mat cC;
130  arma::mat P_im;
132  arma::mat K_im;
133 
136 } rscf_t;
137 
139 typedef struct {
141  arma::mat Ca, Cb;
143  arma::vec Ea, Eb;
145  arma::mat Ha, Hb;
147  arma::mat P, Pa, Pb;
148 
150  arma::mat J;
152  arma::mat Ka, Kb;
154  arma::mat XCa, XCb;
155 
157  arma::mat Ka_im, Kb_im;
158 
160  arma::cx_mat cCa, cCb;
162  arma::mat Pa_im, Pb_im;
163 
166 } uscf_t;
167 
168 
170 enum guess_t {
172  CORE_GUESS,
174  SAD_GUESS,
176  NO_GUESS,
178  GWH_GUESS
179 };
180 
182 typedef struct {
184  bool X;
186  bool C;
188  bool D;
189 } pzmet_t;
190 
191 // Parse PZ method
192 dft_t parse_pzmet(const std::string & str, const dft_t & ovmethod);
193 
195 enum pzham {
197  PZSYMM,
199  PZUNITED
200 };
201 
202 // Parse PZ hamiltonian
203 enum pzham parse_pzham(const std::string & str);
204 
205 class SCF {
206  protected:
208  arma::mat S;
209 
211  arma::mat T;
213  arma::mat Vnuc;
215  arma::mat Hcore;
216 
218  const BasisSet * basisp;
223 
225  arma::mat Sinvh;
226 
228  size_t Nbf;
229 
231  int Nel;
232 
234  int mult;
235 
237  enum guess_t guess;
238 
240  bool usediis;
242  bool diis_c1;
246  double diiseps;
248  double diisthr;
249 
251  bool useadiis;
255  bool usetrrh;
256 
258  int maxiter;
260  double shift;
262  bool verbose;
263 
265  bool direct;
267  bool decfock;
269  bool strictint;
271  double intthr;
272 
276  size_t fitmem;
278  double fitthr;
279 
281  bool cholesky;
283  double cholthr;
285  double cholshthr;
287  int cholmode;
288 
290  bool doforce;
291 
293  double Enuc;
294 
309 
313  arma::mat decconv;
314 
316  std::vector<arma::mat> freeze;
317 
318  public:
320  SCF(const BasisSet & basis, const Settings & set, Checkpoint & chkpt);
321  ~SCF();
322 
324  void RHF(rscf_t & sol, const std::vector<double> & occs, const convergence_t conv);
326  void ROHF(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, const convergence_t conv);
328  void UHF(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, const convergence_t conv);
329 
331  void RDFT(rscf_t & sol, const std::vector<double> & occs, const convergence_t conv, const dft_t dft);
333  void UDFT(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, const convergence_t conv, const dft_t dft);
334 
336  void Fock_RHF(rscf_t & sol, const std::vector<double> & occs) const;
338  void Fock_ROHF(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb) const;
340  void Fock_UHF(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb) const;
341 
343  void Fock_RDFT(rscf_t & sol, const std::vector<double> & occs, const dft_t dft, DFTGrid & grid, DFTGrid & nlgrid) const;
345  void Fock_UDFT(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, const dft_t dft, DFTGrid & grid, DFTGrid & nlgrid) const;
346 
348  void PZSIC_Fock(std::vector<arma::cx_mat> & Forb, arma::vec & Eorb, const arma::cx_mat & C, dft_t dft, DFTGrid & grid, DFTGrid & nlgrid, bool fock);
349 
351  arma::vec force_RHF(rscf_t & sol, const std::vector<double> & occs, double tol);
353  arma::vec force_ROHF(uscf_t & sol, int Nel_alpha, int Nel_beta, double tol);
355  arma::vec force_UHF(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, double tol);
357  arma::vec force_RDFT(rscf_t & sol, const std::vector<double> & occs, const dft_t dft, DFTGrid & grid, DFTGrid & nlgrid, double tol);
359  arma::vec force_UDFT(uscf_t & sol, const std::vector<double> & occa, const std::vector<double> & occb, const dft_t dft, DFTGrid & grid, DFTGrid & nlgrid, double tol);
360 
362  void set_frozen(const arma::mat & C, size_t ind);
363 
365  void set_fitting(const BasisSet & fitbas);
366 
368  void set_verbose(bool verb);
370  bool get_verbose() const;
371 
373  void do_force(bool val);
374 
376  size_t get_maxiter() const;
378  void set_maxiter(size_t maxiter);
379 
381  arma::mat get_S() const;
383  arma::mat get_Sinvh() const;
385  arma::mat get_Hcore() const;
387  Checkpoint *get_checkpoint() const;
389  bool get_strictint() const;
390 
392  void fill_rs(double omega);
393 
395  void core_guess(rscf_t & sol) const;
397  void core_guess(uscf_t & sol) const;
398 
400  void gwh_guess(rscf_t & sol) const;
402  void gwh_guess(uscf_t & sol) const;
403 
405  arma::mat exchange_localization(const arma::mat & Co, const arma::mat & Cv) const;
406 };
407 
409 void imag_lost(const rscf_t & sol, const arma::mat & S, double & d);
411 void imag_lost(const uscf_t & sol, const arma::mat & S, double & da, double & db);
412 
414 void diagonalize(const arma::mat & S, const arma::mat & Sinvh, rscf_t & sol, double shift=0.0);
416 void diagonalize(const arma::mat & S, const arma::mat & Sinvh, uscf_t & sol, double shift=0.0);
417 
418 
427 void ROHF_update(arma::mat & Fa, arma::mat & Fb, const arma::mat & P, const arma::mat & S, std::vector<double> occa, std::vector<double> occb, bool verbose=true);
428 
429 
431 void determine_occ(arma::vec & nocc, const arma::mat & C, const arma::vec & nocc_old, const arma::mat & C_old, const arma::mat & S);
432 
434 arma::mat form_density(const arma::mat & C, size_t nocc);
436 arma::mat form_density(const arma::mat & C, const arma::vec & occs);
438 void form_density(rscf_t & sol, size_t nocc);
440 void form_density(rscf_t & sol, const arma::vec & occa);
442 void form_density(uscf_t & sol, const arma::vec & occa, const arma::vec & occb);
444 arma::mat form_density(const arma::vec & E, const arma::mat & C, const std::vector<double> & nocc);
446 arma::mat purify_density(const arma::mat & P, const arma::mat & S);
448 arma::mat purify_density_NO(const arma::mat & P, const arma::mat & S);
450 arma::mat purify_density_NO(const arma::mat & P, arma::mat & C, const arma::mat & S);
451 
453 std::vector<double> atomic_occupancy(int Nel);
455 std::vector<double> get_restricted_occupancy(const Settings & set, const BasisSet & basis, bool atomic=false);
457 void get_unrestricted_occupancy(const Settings & set, const BasisSet & basis, std::vector<double> & occa, std::vector<double> & occb, bool atomic=false);
458 
460 double dip_mom(const arma::mat & P, const BasisSet & basis);
462 arma::vec dipole_moment(const arma::mat & P, const BasisSet & basis);
464 double electron_spread(const arma::mat & P, const BasisSet & basis);
465 
467 void get_Nel_alpha_beta(int Nel, int mult, int & Nel_alpha, int & Nel_beta);
468 
470 dft_t parse_dft(const Settings & set, bool init);
471 
473 void calculate(const BasisSet & basis, const Settings & set, bool doforce=false);
474 
476 typedef struct {
478  double S;
480  size_t idx;
481 } ovl_sort_t;
482 
484 bool operator<(const ovl_sort_t & lhs, const ovl_sort_t & rhs);
485 
492 arma::mat project_orbitals(const arma::mat & Cold, const BasisSet & minbas, const BasisSet & augbas);
493 
495 std::vector<int> symgroups(const arma::mat & C, const arma::mat & S, const std::vector<arma::mat> & freeze, bool verbose=false);
496 
498 void freeze_orbs(const std::vector<arma::mat> & freeze, const arma::mat & C, const arma::mat & S, arma::mat & H, bool verbose=false);
499 
501 size_t localize_core(const BasisSet & basis, int nocc, arma::mat & C, bool verbose=false);
502 
504 arma::mat interpret_force(const arma::vec & f);
505 
507 #include "checkpoint.h"
508 
509 #endif
Perdew-Zunger SIC mode.
Definition: scf.h:182
arma::mat T
Kinetic energy matrix.
Definition: scf.h:211
int diisorder
Number of DIIS matrices to use.
Definition: scf.h:244
arma::vec Ea
Orbital energies.
Definition: scf.h:143
ERItable tab
Electron repulsion table.
Definition: scf.h:296
double Eel
Total electronic energy.
Definition: scf.h:96
arma::mat H
Fock operator.
Definition: scf.h:116
ERItable tab_rs
Electron repulsion table, range separation.
Definition: scf.h:298
int c_func
Used correlation functional.
Definition: scf.h:63
void RHF(rscf_t &sol, const std::vector< double > &occs, const convergence_t conv)
Calculate restricted Hartree-Fock solution.
arma::mat P_im
Imaginary part of complex-CMO density matrix (for complex exchange contribution)
Definition: scf.h:130
double fitthr
Threshold for density fitting.
Definition: scf.h:278
arma::vec force_UHF(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, double tol)
Calculate force in unrestricted Hartree-Fock.
void Fock_UHF(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb) const
Calculate unrestricted Hartree-Fock operator.
bool X
Exchange?
Definition: scf.h:184
bool direct
Direct calculation?
Definition: scf.h:265
arma::mat J
Coulomb operator.
Definition: scf.h:150
int maxiter
Maximum number of iterations.
Definition: scf.h:258
double Enucr
Nuclear repulsion energy.
Definition: scf.h:98
arma::mat Vnuc
Nuclear attraction matrix.
Definition: scf.h:213
DFT settings.
Definition: scf.h:59
arma::vec force_UDFT(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, const dft_t dft, DFTGrid &grid, DFTGrid &nlgrid, double tol)
Calculate force in unrestricted density-functional theory.
Convergence criteria.
Definition: scf.h:49
std::vector< arma::mat > freeze
List of frozen orbitals by symmetry group. index+1 is symmetry group, group 0 contains all non-frozen...
Definition: scf.h:316
arma::mat XCa
KS-XC matrix.
Definition: scf.h:154
void set_fitting(const BasisSet &fitbas)
Set the density-fitting basis set.
Definition: scf-base.cpp:379
void RDFT(rscf_t &sol, const std::vector< double > &occs, const convergence_t conv, const dft_t dft)
Calculate restricted density-functional theory solution.
arma::mat J
Coulomb operator.
Definition: scf.h:121
void UHF(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, const convergence_t conv)
Calculate unrestricted Hartree-Fock solution.
bool C
Correlation?
Definition: scf.h:186
bool densityfit
Density fitting calculation?
Definition: scf.h:274
energy_t en
Energy information.
Definition: scf.h:165
arma::mat S
Overlap matrix.
Definition: scf.h:208
bool verbose
Verbose calculation?
Definition: scf.h:262
int Nel
Total number of electrons.
Definition: scf.h:231
arma::mat Ka_im
Imaginary exchange.
Definition: scf.h:157
Density fitting routines.
Definition: density_fitting.h:56
double cholshthr
Cholesky shell threshold (for caching)
Definition: scf.h:285
bool cholesky
Cholesky calculation?
Definition: scf.h:281
arma::mat get_Hcore() const
Get core Hamiltonian matrix.
Definition: scf-base.cpp:510
int x_func
Used exchange functional.
Definition: scf.h:61
Screening of electron repulsion integrals.
Definition: eriscreen.h:46
arma::mat P
Density matrices.
Definition: scf.h:147
int mult
Multiplicity.
Definition: scf.h:234
arma::mat Ha
Fock operators.
Definition: scf.h:145
double diiseps
Threshold of enabling use of DIIS.
Definition: scf.h:246
arma::vec force_RHF(rscf_t &sol, const std::vector< double > &occs, double tol)
Calculate force in restricted Hartree-Fock.
Energy info.
Definition: scf.h:83
Unrestricted solver info.
Definition: scf.h:139
int cholmode
Cholesky mode.
Definition: scf.h:287
arma::mat Sinvh
Basis orthogonalizing matrix.
Definition: scf.h:225
Cholesky decomposition of ERIs.
Definition: erichol.h:26
void gwh_guess(rscf_t &sol) const
Do GWH guess.
Definition: scf-base.cpp:866
DensityFit dfit
Density fitting table.
Definition: scf.h:308
double deltaPmax
Convergence criterion for maximum change of an element of the density matrix.
Definition: scf.h:53
size_t Nbf
Amount of basis functions.
Definition: scf.h:228
bool doforce
Calculate forces?
Definition: scf.h:290
double Enl
Non-local energy.
Definition: scf.h:100
void set_frozen(const arma::mat &C, size_t ind)
Set frozen orbitals in ind:th symmetry group. ind+1 is the resulting symmetry group, group 0 contains all non-frozen orbitals.
Definition: scf-base.cpp:369
void ROHF(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, const convergence_t conv)
Calculate restricted open-shell Hartree-Fock solution.
ERIscreen scr_rs
Electron repulsion screening table, range separation.
Definition: scf.h:302
Helper for sorting orbitals into maximum overlap.
Definition: scf.h:476
double E
Total energy.
Definition: scf.h:105
Restricted solver info.
Definition: scf.h:110
Settings used for a calculations.
Definition: settings.h:74
size_t idx
Index.
Definition: scf.h:480
void PZSIC_Fock(std::vector< arma::cx_mat > &Forb, arma::vec &Eorb, const arma::cx_mat &C, dft_t dft, DFTGrid &grid, DFTGrid &nlgrid, bool fock)
Helper for PZ-SIC: compute orbital-dependent Fock matrices.
Definition: scf-base.cpp:522
enum guess_t guess
Which guess to use.
Definition: scf.h:237
double Exc
Exchange(-correlation) energy.
Definition: scf.h:91
void fill_rs(double omega)
Fill range-separated integrals.
Definition: scf-base.cpp:403
size_t get_maxiter() const
Get maximum iterations.
Definition: scf-base.cpp:387
void Fock_RDFT(rscf_t &sol, const std::vector< double > &occs, const dft_t dft, DFTGrid &grid, DFTGrid &nlgrid) const
Calculate restricted density-functional theory KS-Fock operator.
double gridtol
Integration grid tolerance.
Definition: scf.h:68
arma::vec force_ROHF(uscf_t &sol, int Nel_alpha, int Nel_beta, double tol)
Calculate force in restricted open-shell Hartree-Fock.
int lmax
Maximum angular quantum number to integrate exactly (if not adaptive)
Definition: scf.h:72
ERIchol chol
Cholesky integrals.
Definition: scf.h:304
bool adaptive
Adaptive grid?
Definition: scf.h:66
arma::cx_mat cCa
Complex orbitals (for SIC)
Definition: scf.h:160
arma::mat K_im
Imaginary exchange.
Definition: scf.h:132
arma::mat get_S() const
Get overlap matrix.
Definition: scf-base.cpp:502
arma::mat Hcore
Core Hamiltonian.
Definition: scf.h:215
double intthr
Integral screening threshold.
Definition: scf.h:271
arma::mat C
Orbitals.
Definition: scf.h:112
Self-consistent field solver routines.
Definition: scf.h:205
arma::mat decconv
Conversion matrix.
Definition: scf.h:313
arma::mat Ka
Exchange operators.
Definition: scf.h:152
bool lobatto
Use Lobatto quadrature?
Definition: scf.h:74
void Fock_UDFT(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, const dft_t dft, DFTGrid &grid, DFTGrid &nlgrid) const
Calculate unrestricted density-functional theory KS-Fock operator.
bool get_verbose() const
Set verbose setting.
Definition: scf-base.cpp:395
BasisSet dfitbas
Density fitting basis.
Definition: scf.h:220
double deltaPrms
Convergence criterion for the RMS change of the density matrix.
Definition: scf.h:55
bool D
Non-local correlation?
Definition: scf.h:188
DFT quadrature grid.
Definition: dftgrid.h:475
Table of electron repulsion integrals.
Definition: eritable.h:44
Basis set.
Definition: basis.h:187
double Ecoul
Coulombic energy.
Definition: scf.h:85
int nrad
Amount of radial shells (if not adaptive)
Definition: scf.h:70
double diisthr
Threshold of enabling full use of DIIS.
Definition: scf.h:248
void set_verbose(bool verb)
Set verbose setting.
Definition: scf-base.cpp:399
bool useadiis
Use ADIIS?
Definition: scf.h:251
void set_maxiter(size_t maxiter)
Set maximum iterations.
Definition: scf-base.cpp:391
bool get_strictint() const
Using strict integrals?
Definition: scf-base.cpp:518
arma::cx_mat cC
Complex orbitals (for SIC)
Definition: scf.h:128
const BasisSet * basisp
Basis set to use (needed for DFT grid operation)
Definition: scf.h:218
void Fock_RHF(rscf_t &sol, const std::vector< double > &occs) const
Calculate restricted Hartree-Fock operator.
arma::mat K
Exchange operator.
Definition: scf.h:123
BasisSet decbas
Decontracted basis set.
Definition: scf.h:311
bool usebroyden
Use Broyden accelerator?
Definition: scf.h:253
arma::mat get_Sinvh() const
Get half-inverse overlap matrix.
Definition: scf-base.cpp:506
double S
Overlap.
Definition: scf.h:478
bool diis_c1
Use C1-DIIS instead of C2-DIIS?
Definition: scf.h:242
double deltaEmax
Convergence criterion for change of energy.
Definition: scf.h:51
Checkpoint * get_checkpoint() const
Get checkpoint file.
Definition: scf-base.cpp:514
arma::vec force_RDFT(rscf_t &sol, const std::vector< double > &occs, const dft_t dft, DFTGrid &grid, DFTGrid &nlgrid, double tol)
Calculate force in restricted density-functional theory.
Checkpoint * chkptp
Checkpoint file.
Definition: scf.h:222
SCF(const BasisSet &basis, const Settings &set, Checkpoint &chkpt)
Constructor.
Definition: scf-base.cpp:63
double shift
Level shift.
Definition: scf.h:260
arma::mat P
Density matrix.
Definition: scf.h:118
double Eone
One-electron contribution.
Definition: scf.h:93
bool strictint
Strict integrals?
Definition: scf.h:269
void do_force(bool val)
Toggle calculation of forces.
Definition: scf-base.cpp:383
arma::vec E
Orbital energies.
Definition: scf.h:114
double Enuca
Nuclear attraction energy.
Definition: scf.h:89
energy_t en
Energy information.
Definition: scf.h:135
double cholthr
Cholesky threshold.
Definition: scf.h:283
arma::mat exchange_localization(const arma::mat &Co, const arma::mat &Cv) const
Exchange localization.
Definition: scf-base.cpp:893
arma::mat Pa_im
Imaginary part of complex-CMO density matrix (for complex exchange contribution)
Definition: scf.h:162
ERIscreen scr
Electron repulsion screening table (for direct calculations)
Definition: scf.h:300
double Enuc
Nuclear repulsion energy.
Definition: scf.h:293
void core_guess(rscf_t &sol) const
Do core guess.
Definition: scf-base.cpp:851
bool usediis
Use DIIS?
Definition: scf.h:240
arma::mat XC
KS-XC matrix.
Definition: scf.h:125
arma::mat Ca
Orbitals.
Definition: scf.h:141
double Ekin
Kinetic energy.
Definition: scf.h:87
void Fock_ROHF(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb) const
Calculate restricted open-shell Hartree-Fock operator.
bool decfock
Use decontracted basis to construct Fock matrix? (Direct formation)
Definition: scf.h:267
void UDFT(uscf_t &sol, const std::vector< double > &occa, const std::vector< double > &occb, const convergence_t conv, const dft_t dft)
Calculate unrestricted density-functional theory solution.
ERIchol chol_rs
Cholesky integrals, range separation.
Definition: scf.h:306
double Esic
Self-interaction energy.
Definition: scf.h:102
bool usetrrh
Use Trust-Region Roothaan-Hall?
Definition: scf.h:255
Checkpointing class.
Definition: checkpoint.h:69
size_t fitmem
Memory allocation for density fitting.
Definition: scf.h:276