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-2015 9 * Copyright (c) 2010-2015, 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 #ifndef ERKALE_ERIDIGEST 18 #define ERKALE_ERIDIGEST 19 20 #include "basis.h" 21 class dERIWorker; 22 23 /// Integral digestor 24 class IntegralDigestor { 25 public: 26 /// Constructor 27 IntegralDigestor(); 28 /// Destructor 29 virtual ~IntegralDigestor(); 30 /// Digest integral block 31 virtual void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, const std::vector<double> & ints, size_t ioff)=0; 32 }; 33 34 /// Coulomb matrix digestor 35 class JDigestor: public IntegralDigestor { 36 /// Density matrix 37 arma::mat P; 38 /// Coulomb matrix 39 arma::mat J; 40 public: 41 /// Construct digestor 42 JDigestor(const arma::mat & P); 43 /// Destruct digestor 44 ~JDigestor(); 45 46 /// Digest integrals 47 void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, const std::vector<double> & ints, size_t ioff); 48 /// Get output 49 arma::mat get_J() const; 50 }; 51 52 /// Exchange matrix digestor 53 class KDigestor: public IntegralDigestor { 54 /// Density matrix 55 arma::mat P; 56 /// Exchange matrix 57 arma::mat K; 58 public: 59 /// Construct digestor 60 KDigestor(const arma::mat & P); 61 /// Destruct digestor 62 ~KDigestor(); 63 64 /// Digest integrals 65 void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, const std::vector<double> & ints, size_t ioff); 66 /// Get output 67 arma::mat get_K() const; 68 }; 69 70 /// Complex exchange matrix digestor 71 class cxKDigestor: public IntegralDigestor { 72 /// Density matrix 73 arma::cx_mat P; 74 /// Exchange matrix 75 arma::cx_mat K; 76 public: 77 /// Construct digestor 78 cxKDigestor(const arma::cx_mat & P); 79 /// Destruct digestor 80 ~cxKDigestor(); 81 82 /// Digest integrals 83 void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, const std::vector<double> & ints, size_t ioff); 84 /// Get output 85 arma::cx_mat get_K() const; 86 }; 87 88 /// Force digestor 89 class ForceDigestor { 90 public: 91 /// Constructor 92 ForceDigestor(); 93 /// Destructor 94 virtual ~ForceDigestor(); 95 /// Digest derivative block 96 virtual void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, dERIWorker & deriw, arma::vec & f)=0; 97 }; 98 99 /// Coulomb force digestor 100 class JFDigestor: public ForceDigestor { 101 /// Density matrix 102 arma::mat P; 103 public: 104 /// Construct digestor 105 JFDigestor(const arma::mat & P); 106 /// Destruct digestor 107 ~JFDigestor(); 108 109 /// Digest integrals 110 void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, dERIWorker & deriw, arma::vec & f); 111 }; 112 113 /// Exchange force digestor 114 class KFDigestor: public ForceDigestor { 115 /// Density matrix 116 arma::mat P; 117 /// Fraction of exact exchange 118 double kfrac; 119 /// Degeneracy factor 120 double fac; 121 122 public: 123 /// Construct digestor 124 KFDigestor(const arma::mat & P, double kfrac, bool restr); 125 /// Destruct digestor 126 ~KFDigestor(); 127 128 /// Digest integrals 129 void digest(const std::vector<eripair_t> & shpairs, size_t ip, size_t jp, dERIWorker & deriw, arma::vec & f); 130 }; 131 132 #endif 133