ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
density_fitting.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-2012
9  * Copyright (c) 2010-2012, 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 
49 #ifndef ERKALE_DENSITYFIT
50 #define ERKALE_DENSITYFIT
51 
52 #include "global.h"
53 #include "basis.h"
54 
56 class DensityFit {
58  size_t Nbf;
60  size_t Naux;
62  bool direct;
64  bool hf;
65 
67  size_t Nnuc;
69  int maxam;
71  int maxcontr;
72 
74  std::vector<GaussianShell> orbshells;
75  int maxorbam;
76  size_t maxorbcontr;
78  std::vector<GaussianShell> auxshells;
79  int maxauxam;
80  size_t maxauxcontr;
83 
85  size_t dummyind;
86 
88  std::vector<eripair_t> orbpairs;
89 
91  std::vector<size_t> iidx;
92 
94  arma::mat a_munu;
95 
97  arma::mat ab;
99  arma::mat ab_inv;
101  arma::mat ab_invh;
102 
104  void form_screening();
105 
106  public:
108  DensityFit();
110  ~DensityFit();
111 
118  size_t fill(const BasisSet & orbbas, const BasisSet & auxbas, bool direct, double erithr, double linthr, bool hf=false);
120  size_t idx(size_t imu, size_t inu) const;
121 
123  size_t memory_estimate(const BasisSet & orbbas, const BasisSet & auxbas, bool direct) const;
124 
126  arma::vec compute_expansion(const arma::mat & P) const;
128  std::vector<arma::vec> compute_expansion(const std::vector<arma::mat> & P) const;
129 
131  arma::mat calcJ(const arma::mat & P) const;
133  std::vector<arma::mat> calcJ(const std::vector<arma::mat> & P) const;
134 
136  arma::vec forceJ(const arma::mat & P);
137 
139  arma::mat calcK(const arma::mat & C, const std::vector<double> & occs, size_t memlimit) const;
141  arma::cx_mat calcK(const arma::cx_mat & C, const std::vector<double> & occs, size_t memlimit) const;
142 
144  size_t get_Norb() const;
146  size_t get_Naux() const;
148  double get_a_munu(size_t ia, size_t imu, size_t inu) const;
150  arma::mat get_ab_inv() const;
151 
153  void B_matrix(arma::mat & B) const;
154 };
155 
156 
157 #endif
bool direct
Direct calculation? (Compute three-center integrals on-the-fly)
Definition: density_fitting.h:62
std::vector< size_t > iidx
Index helper.
Definition: density_fitting.h:91
A shell of Gaussian basis functions of a given angular momentum, sharing the same exponential contrac...
Definition: basis.h:482
void form_screening()
Form screening matrix.
GaussianShell dummy
Dummy shell.
Definition: density_fitting.h:82
size_t dummyind
Index of dummy function.
Definition: density_fitting.h:85
size_t Naux
Amount of auxiliary basis functions.
Definition: density_fitting.h:60
arma::mat ab_invh
Definition: density_fitting.h:101
Density fitting routines.
Definition: density_fitting.h:56
arma::mat ab_inv
Definition: density_fitting.h:99
~DensityFit()
Destructor.
Definition: density_fitting.cpp:29
std::vector< eripair_t > orbpairs
List of unique orbital shell pairs.
Definition: density_fitting.h:88
size_t idx(size_t imu, size_t inu) const
Compute index in integral table.
Definition: density_fitting.cpp:178
void B_matrix(arma::mat &B) const
Get B matrix (must have HF enabled)
Definition: density_fitting.cpp:1177
arma::mat a_munu
Integrals .
Definition: density_fitting.h:94
DensityFit()
Constructor.
Definition: density_fitting.cpp:26
size_t fill(const BasisSet &orbbas, const BasisSet &auxbas, bool direct, double erithr, double linthr, bool hf=false)
Definition: density_fitting.cpp:33
double get_a_munu(size_t ia, size_t imu, size_t inu) const
Get the three-electron integral.
Definition: density_fitting.cpp:1161
size_t Nnuc
Amount of nuclei.
Definition: density_fitting.h:67
arma::mat get_ab_inv() const
Get ab_inv.
Definition: density_fitting.cpp:1173
std::vector< GaussianShell > orbshells
Orbital shells.
Definition: density_fitting.h:74
Basis set.
Definition: basis.h:187
arma::vec compute_expansion(const arma::mat &P) const
Compute expansion coefficients c.
Definition: density_fitting.cpp:208
size_t Nbf
Amount of orbital basis functions.
Definition: density_fitting.h:58
std::vector< GaussianShell > auxshells
Density fitting shells.
Definition: density_fitting.h:78
bool hf
Hartree-Fock calculation?
Definition: density_fitting.h:64
size_t memory_estimate(const BasisSet &orbbas, const BasisSet &auxbas, bool direct) const
Compute estimate of necessary memory.
Definition: density_fitting.cpp:184
size_t get_Naux() const
Get the number of auxiliary functions.
Definition: density_fitting.cpp:1157
arma::mat calcJ(const arma::mat &P) const
Get Coulomb matrix from P.
Definition: density_fitting.cpp:414
arma::vec forceJ(const arma::mat &P)
Calculate force from P.
Definition: density_fitting.cpp:579
arma::mat ab
Definition: density_fitting.h:97
size_t get_Norb() const
Get the number of orbital functions.
Definition: density_fitting.cpp:1153
arma::mat calcK(const arma::mat &C, const std::vector< double > &occs, size_t memlimit) const
Get exchange matrix from orbitals with occupation numbers occs.
Definition: density_fitting.cpp:766
int maxam
Maximum angular momentum.
Definition: density_fitting.h:69
int maxcontr
Maximum contractions.
Definition: density_fitting.h:71