Species.h 6.17 KB
Newer Older
Francois Gygi committed
1 2
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
3 4 5 6
// Copyright (c) 2008 The Regents of the University of California
//
// This file is part of Qbox
//
Francois Gygi committed
7 8
// Qbox is distributed under the terms of the GNU General Public License
// as published by the Free Software Foundation, either version 2 of
Francois Gygi committed
9 10 11 12 13 14
// the License, or (at your option) any later version.
// See the file COPYING in the root directory of this distribution
// or <http://www.gnu.org/licenses/>.
//
////////////////////////////////////////////////////////////////////////////////
//
Francois Gygi committed
15 16 17 18 19 20 21 22 23 24 25
// Species.h:
//
////////////////////////////////////////////////////////////////////////////////

#ifndef SPECIES_H
#define SPECIES_H

#include <iostream>
#include <string>
#include <vector>
#include <cmath>
26
#include <utility>
Francois Gygi committed
27 28 29 30

class Species
{
  private:
31

32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48
  struct ProjectorData
  {
    int l, m, n;
  };

  enum PP_type
  {
    // norm conserving pseudopotential
    NCPP,
    // norm conserving semilocal pseudopotential
    SLPP
  };

  PP_type type_; // identify type of the PP

  int nlm_;             // number of non-local projectors (including m)
  int nop_;             // number of non-local projectors (excluding m)
Francois Gygi committed
49
  int ndft_;
50

51 52 53
  std::vector<std::vector<double> > vps_spl_, vps_spl2_, phi_spl_, phi_spl2_;
  std::vector<double>               gspl_, vlocg_spl_, vlocg_spl2_;
  std::vector<std::vector<double> > vnlg_spl_, vnlg_spl2_;
54 55
  std::vector<double> wsg_;  // wsg_[l] Kleinman-Bylander weight
                             // 1/<phi|delta_V|phi>
56 57
  std::vector<double>               nlcc_spl_, nlcc_spl2_;
  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
Francois Gygi committed
63

64
  std::string symbol_;
Francois Gygi committed
65 66
  int atomic_number_;
  double mass_;        // mass in a.m.u (Carbon = 12.0)
67

68
  std::string description_; // description of the pseudopotential
Francois Gygi committed
69
  int zval_;           // valence charge
70 71
  double zcore_;       // core charge
  double ztot_;        // total charge
Francois Gygi committed
72 73 74 75 76
  int lmax_;           // largest angular momentum
  int llocal_;         // angular momentum taken as local
  int nquad_;          // number of semi-local quadrature points
  double rquad_;       // end of semi-local quadrature interval
  double deltar_;      // mesh spacing for potentials and wavefunctions
77 78
  std::vector<std::vector<double> > vps_;  // potentials for each l (input)
  std::vector<std::vector<double> > phi_;  // atomic wf for each l (input)
Francois Gygi committed
79
  double rcps_;        // cutoff radius of gaussian pseudocharge
80 81
  std::vector<double> nlcc_; // non linear core correction

82
  // map index of projector -> angular momentum
83 84 85
  std::vector<int> lmap_;
  // local potential in radial representation
  std::vector<double> vlocal_;
86
  // projector in radial representation, indices proj_[l][n][r]
87 88 89 90 91 92 93 94 95 96 97 98 99
  std::vector<std::vector<std::vector<double> > > proj_;
  // matrix D in block diagonal storage, indices d_[l,n,m]
  std::vector<std::vector<std::vector<double> > > d_;

  // initialize a norm conserving PP
  bool initialize_ncpp();
  // initialize a semilocal potential
  bool initialize_slpp();
  // non linear core correction
  void initialize_nlcc();

  // helper function that extracts l, m and n from projector index
  ProjectorData get_proj_data(int ipr);
100

Francois Gygi committed
101 102
  public:

103
  Species(std::string name);
104

105 106 107 108
  const std::string& name(void) const { return name_; }
  const std::string& symbol(void) const { return symbol_; }
  const std::string& description(void) const { return description_; }
  const std::string& uri(void) const { return uri_; }
Francois Gygi committed
109 110
  int atomic_number(void) const { return atomic_number_; }
  int zval(void) const { return zval_; }
111 112
  double zcore(void) const { return zcore_; }
  double ztot(void) const { return ztot_; }
Francois Gygi committed
113 114 115 116 117 118 119
  double mass(void) const { return mass_; }
  int lmax(void) const { return lmax_; }
  int llocal(void) const { return llocal_; }
  int nquad(void) const { return nquad_; }
  double rquad(void) const { return rquad_; }
  double deltar(void) const { return deltar_; }
  double rcps(void) const { return rcps_; }
120 121
  bool has_nlcc(void) const { return nlcc_.size() > 0; }
  bool has_dmatrix(void) const { return d_.size() > 0; }
Francois Gygi committed
122 123 124

  // number of non-local projectors sum_(l!=llocal) (2*l+1)
  int nlm(void) { return nlm_; }
125 126
  // number of non-local projectors w/o m degeneracy
  int nop(void) { return nop_; }
127
  // angular momentum of projector with index iop
128 129 130 131 132
  int l(int iop) { return lmap_[iop]; }
  // extract D matrix
  double dmatrix(int ipr, int jpr);

  bool non_local(void) { return nop_ > 0; };
Francois Gygi committed
133
  double eself(void)
134
  { return ztot_ * ztot_ / ( sqrt ( 2.0 * M_PI ) * rcps_ ); };
Francois Gygi committed
135

Francois Gygi committed
136
  void phi(int l, double r, double &val);              // phi(l,r) in r space
Francois Gygi committed
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
Francois Gygi committed
139 140
  void vlocg(double q, double &v);                     // Vloc(g)
  void dvlocg(double q, double &v, double &dv);        // Vloc(g) and dVloc/dg
141 142
  void vnlg(int iop, double q, double &v);             // Vnl(iop,g)
  void dvnlg(int iop, double q, double &v, double &dv);// Vnl(iop,g) and dVnl/dg
Francois Gygi committed
143
  double rhopsg(double q);        // pseudocharge in g space
144 145 146 147 148
  void rhocoreg(double q, double &rho);      // core correction in g space
  // core correction and its derivative in g space
  void drhocoreg(double q, double &rho, double &drho);

  double wsg(int iop) { return wsg_[iop]; };
Francois Gygi committed
149
  double rcut_loc(double epsilon); // radius beyond which potential is local
150

151 152
  const std::vector<std::vector<double> >& vps(void) const { return vps_spl_; }
  const std::vector<std::vector<double> >& phi(void) const { return phi_spl_; }
153

Francois Gygi committed
154
  bool initialize(double rcps);
155
  void info(std::ostream& os);
156
  void print(std::ostream& os, bool expanded_form);
157 158
  void print_nlcc(std::ostream& os);
  void print_radial_function(std::ostream& os, const std::vector<double>& rad_func);
159

Francois Gygi committed
160 161
  friend class SpeciesReader;
  friend class SpeciesHandler;
162

Francois Gygi committed
163
};
164
std::ostream& operator << ( std::ostream &os, Species &a );
Francois Gygi committed
165 166 167
class SpeciesInitException
{
  public:
168 169
  std::string msg;
  SpeciesInitException(std::string s) : msg(s) {}
Francois Gygi committed
170 171
};
#endif