1 //////////////////////////////////////////////////////////////////////////////// 2 // 3 // Copyright (c) 2008 The Regents of the University of California 4 // 5 // This file is part of Qbox 6 // 7 // Qbox is distributed under the terms of the GNU General Public License 8 // as published by the Free Software Foundation, either version 2 of 9 // the License, or (at your option) any later version. 10 // See the file COPYING in the root directory of this distribution 11 // or <http://www.gnu.org/licenses/>. 12 // 13 //////////////////////////////////////////////////////////////////////////////// 14 // 15 // Species.h: 16 // 17 //////////////////////////////////////////////////////////////////////////////// 18 19 #ifndef SPECIES_H 20 #define SPECIES_H 21 22 #include <iostream> 23 #include <string> 24 #include <vector> 25 #include <cmath> 26 #include <utility> 27 28 class Species 29 { 30 private: 31 32 struct ProjectorData 33 { 34 int l, m, n; 35 }; 36 37 enum PP_type 38 { 39 // norm conserving pseudopotential 40 NCPP, 41 // norm conserving semilocal pseudopotential 42 SLPP 43 }; 44 45 PP_type type_; // identify type of the PP 46 47 int nlm_; // number of non-local projectors (including m) 48 int nop_; // number of non-local projectors (excluding m) 49 int ndft_; 50 51 std::vector<std::vector<double> > vps_spl_, vps_spl2_, phi_spl_, phi_spl2_; 52 std::vector<double> gspl_, vlocg_spl_, vlocg_spl2_; 53 std::vector<std::vector<double> > vnlg_spl_, vnlg_spl2_; 54 std::vector<double> wsg_; // wsg_[l] Kleinman-Bylander weight 55 // 1/<phi|delta_V|phi> 56 std::vector<double> nlcc_spl_, nlcc_spl2_; 57 std::vector<double> nlccg_spl_, nlccg_spl2_; 58 59 std::vector<double> rps_spl_; // radial linear mesh (same for all l) 60 61 std::string name_; // name used in current application 62 std::string uri_; // uri of the resource defining the pseudopotential 63 64 std::string symbol_; 65 int atomic_number_; 66 double mass_; // mass in a.m.u (Carbon = 12.0) 67 68 std::string description_; // description of the pseudopotential 69 int zval_; // valence charge 70 double zcore_; // core charge 71 double ztot_; // total charge 72 int lmax_; // largest angular momentum 73 int llocal_; // angular momentum taken as local 74 int nquad_; // number of semi-local quadrature points 75 double rquad_; // end of semi-local quadrature interval 76 double deltar_; // mesh spacing for potentials and wavefunctions 77 std::vector<std::vector<double> > vps_; // potentials for each l (input) 78 std::vector<std::vector<double> > phi_; // atomic wf for each l (input) 79 double rcps_; // cutoff radius of gaussian pseudocharge 80 std::vector<double> nlcc_; // non linear core correction 81 82 // map index of projector -> angular momentum 83 std::vector<int> lmap_; 84 // local potential in radial representation 85 std::vector<double> vlocal_; 86 // projector in radial representation, indices proj_[l][n][r] 87 std::vector<std::vector<std::vector<double> > > proj_; 88 // matrix D in block diagonal storage, indices d_[l,n,m] 89 std::vector<std::vector<std::vector<double> > > d_; 90 91 // initialize a norm conserving PP 92 bool initialize_ncpp(); 93 // initialize a semilocal potential 94 bool initialize_slpp(); 95 // non linear core correction 96 void initialize_nlcc(); 97 98 // helper function that extracts l, m and n from projector index 99 ProjectorData get_proj_data(int ipr); 100 101 public: 102 103 Species(std::string name); 104 name(void)105 const std::string& name(void) const { return name_; } symbol(void)106 const std::string& symbol(void) const { return symbol_; } description(void)107 const std::string& description(void) const { return description_; } uri(void)108 const std::string& uri(void) const { return uri_; } atomic_number(void)109 int atomic_number(void) const { return atomic_number_; } zval(void)110 int zval(void) const { return zval_; } zcore(void)111 double zcore(void) const { return zcore_; } ztot(void)112 double ztot(void) const { return ztot_; } mass(void)113 double mass(void) const { return mass_; } lmax(void)114 int lmax(void) const { return lmax_; } llocal(void)115 int llocal(void) const { return llocal_; } nquad(void)116 int nquad(void) const { return nquad_; } rquad(void)117 double rquad(void) const { return rquad_; } deltar(void)118 double deltar(void) const { return deltar_; } rcps(void)119 double rcps(void) const { return rcps_; } has_nlcc(void)120 bool has_nlcc(void) const { return nlcc_.size() > 0; } has_dmatrix(void)121 bool has_dmatrix(void) const { return d_.size() > 0; } 122 123 // number of non-local projectors sum_(l!=llocal) (2*l+1) nlm(void)124 int nlm(void) { return nlm_; } 125 // number of non-local projectors w/o m degeneracy nop(void)126 int nop(void) { return nop_; } 127 // angular momentum of projector with index iop l(int iop)128 int l(int iop) { return lmap_[iop]; } 129 // extract D matrix 130 double dmatrix(int ipr, int jpr); 131 non_local(void)132 bool non_local(void) { return nop_ > 0; }; eself(void)133 double eself(void) 134 { return ztot_ * ztot_ / ( sqrt ( 2.0 * M_PI ) * rcps_ ); }; 135 136 void phi(int l, double r, double &val); // phi(l,r) in r space 137 void vpsr(int l, double r, double &v); // Vps(l,r) in r space 138 void dvpsr(int l, double r, double &v, double &dv); // Vps and dVps/dr 139 void vlocg(double q, double &v); // Vloc(g) 140 void dvlocg(double q, double &v, double &dv); // Vloc(g) and dVloc/dg 141 void vnlg(int iop, double q, double &v); // Vnl(iop,g) 142 void dvnlg(int iop, double q, double &v, double &dv);// Vnl(iop,g) and dVnl/dg 143 double rhopsg(double q); // pseudocharge in g space 144 void rhocoreg(double q, double &rho); // core correction in g space 145 // core correction and its derivative in g space 146 void drhocoreg(double q, double &rho, double &drho); 147 wsg(int iop)148 double wsg(int iop) { return wsg_[iop]; }; 149 double rcut_loc(double epsilon); // radius beyond which potential is local 150 vps(void)151 const std::vector<std::vector<double> >& vps(void) const { return vps_spl_; } phi(void)152 const std::vector<std::vector<double> >& phi(void) const { return phi_spl_; } 153 154 bool initialize(double rcps); 155 void info(std::ostream& os); 156 void print(std::ostream& os, bool expanded_form); 157 void print_nlcc(std::ostream& os); 158 void print_radial_function(std::ostream& os, const std::vector<double>& rad_func); 159 160 friend class SpeciesReader; 161 friend class SpeciesHandler; 162 163 }; 164 std::ostream& operator << ( std::ostream &os, Species &a ); 165 class SpeciesInitException 166 { 167 public: 168 std::string msg; SpeciesInitException(std::string s)169 SpeciesInitException(std::string s) : msg(s) {} 170 }; 171 #endif 172