ERKALE
ERKALE - DFT from Hel
 All Classes Functions Variables Friends Pages
emd.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 
29 #ifndef ERKALE_RADEMD
30 #define ERKALE_RADEMD
31 
32 #include "../basis.h"
33 #include "../global.h"
34 #include "../emd/spherical_expansion.h"
35 
36 #include <complex>
37 #include <vector>
38 #include <armadillo>
39 
42  protected:
44  int l;
45 
46  public:
48  RadialFourier(int l);
50  virtual ~RadialFourier();
51 
53  int getl() const;
54  // Print expansion
55  virtual void print() const = 0;
56 
61  virtual std::complex<double> get(double p) const = 0;
62 };
63 
65 typedef struct {
67  int l;
69  int lp;
70 
71  // couple to
72 
74  int L;
76  int M;
77 
79  std::complex<double> c;
81 
83 bool operator<(const coupl_coeff_t & lhs, const coupl_coeff_t & rhs);
85 bool operator==(const coupl_coeff_t & lhs, const coupl_coeff_t & rhs);
86 
88 typedef struct {
90  int L;
92  int M;
94  std::complex<double> c;
96 
98 bool operator<(const total_coupl_t & lhs, const total_coupl_t & rhs);
100 bool operator==(const total_coupl_t & lhs, const total_coupl_t & rhs);
101 
103 void add_coupling_term(std::vector<total_coupl_t> & v, total_coupl_t & t);
104 
106 typedef struct {
108  int l;
110  std::complex<double> f;
111 } radf_val_t;
112 
114 typedef struct {
116  size_t i;
118  size_t j;
119 } noneqradf_t;
120 
127  std::vector< std::vector<size_t> > idfuncs;
129  std::vector< std::vector<coupl_coeff_t> > cc;
130 
132  std::vector<size_t> loc;
133 
135  size_t Nat;
137  std::vector<double> dist;
139  std::vector< std::vector< std::complex<double> > > YLM;
140 
142  arma::cx_mat P;
143 
145  int Lmax;
146 
148  void distance_table(const std::vector<coords_t> & coord);
149 
151  void compute_coefficients(const std::vector< std::vector<ylmcoeff_t> > & clm, int lp, int mp);
152 
154  void add_coupling(size_t ig, size_t jg, coupl_coeff_t c);
155 
157  void get_coupling(size_t ig, size_t jg, int l, int lp, std::vector<total_coupl_t> & c) const;
158 
160  std::vector<radf_val_t> get_radial(size_t ig, double p) const;
161 
163  void get_total_coupling(size_t ig, size_t jg, double p, std::vector<total_coupl_t> & c, std::vector<total_coupl_t> & tmp) const;
164 
165  protected:
173  std::vector< std::vector<RadialFourier *> > rad;
174 
175  public:
177  EMDEvaluator();
178 
193  EMDEvaluator(const std::vector< std::vector<size_t> > & idfuncsv, const std::vector< std::vector<ylmcoeff_t> > & clm, const std::vector<size_t> & locv, const std::vector<coords_t> & coord, const arma::cx_mat & Pv, int lp=0, int mp=0);
194 
196  ~EMDEvaluator();
197 
199  void print() const;
201  void check_norm() const;
202 
204  std::complex<double> get(double p) const;
205 };
206 
208 arma::mat bessel_array(const std::vector<double> & args, int lmax);
209 
210 
212 typedef struct {
214  double p;
216  double d;
217 } emd_t;
218 
219 
237 class EMD {
239  std::vector<emd_t> dens;
241  void add4(size_t ind);
242 
244  int l;
246  int m;
247 
248  protected:
250  double Nel;
251 
255  std::complex<double> poscoef;
256 
260  std::complex<double> negcoef;
261 
262  public:
264  EMD(const EMDEvaluator * poseval, const EMDEvaluator * negeval, double Nel, int l, int m);
266  ~EMD();
267 
269  double eval(double p) const;
270 
272  void initial_fill(bool verbose=true);
274  void complete_fill();
275 
277  void find_electrons(bool verbose=true, double tol=1e-4);
279  void optimize_moments(const std::vector<int> & moms, bool verbose=true, double tol=1e-8);
281  void optimize_moments(bool verbose=true, double tol=1e-8);
282 
292  void fixed_fill(bool verbose=true, double h0=1e-3, double l0=3.0, double hfac=2.0, double lfac=2.0);
293 
295  std::vector<emd_t> get() const;
296 
298  void save(const char * fname) const;
299 
301  arma::mat moments() const;
303  void moments(const char * fname) const;
304 
306  arma::mat compton_profile() const;
308  void compton_profile(const char * raw) const;
310  void compton_profile_interp(const char * interp) const;
311 };
312 
313 #endif
Functions for evaluating properties of the electron momentum density.
Definition: emd.h:237
std::vector< double > dist
The distances between the functions&#39; origins (Nat x Nat)
Definition: emd.h:137
size_t i
first index
Definition: emd.h:116
void optimize_moments(const std::vector< int > &moms, bool verbose=true, double tol=1e-8)
Optimize given moments of EMD within tolerance.
Definition: emd.cpp:824
double d
Electron momentum density at p.
Definition: emd.h:216
int lp
l&#39;
Definition: emd.h:69
double Nel
Number of electrons.
Definition: emd.h:250
int l
l
Definition: emd.h:67
size_t Nat
The number of centers.
Definition: emd.h:135
std::vector< size_t > loc
The locations of the functions on the atoms (Nbas)
Definition: emd.h:132
void find_electrons(bool verbose=true, double tol=1e-4)
Continue filling until number of electrons is reproduced within tolerance.
Definition: emd.cpp:744
const EMDEvaluator * poseval
Positive evaluator.
Definition: emd.h:253
std::vector< emd_t > dens
List of radial densities.
Definition: emd.h:239
std::complex< double > c
Coefficient.
Definition: emd.h:94
void add_coupling(size_t ig, size_t jg, coupl_coeff_t c)
Add coupling coefficient.
Definition: emd.cpp:329
void compute_coefficients(const std::vector< std::vector< ylmcoeff_t > > &clm, int lp, int mp)
Computes the coupling coefficients.
Definition: emd.cpp:197
size_t j
second index
Definition: emd.h:118
List of identical functions.
Definition: emd.h:114
std::complex< double > f
Value.
Definition: emd.h:110
std::vector< radf_val_t > get_radial(size_t ig, double p) const
Computes the ig:th radial function.
Definition: emd.cpp:315
std::vector< std::vector< RadialFourier * > > rad
Definition: emd.h:173
void check_norm() const
Check norms of radial functions.
Definition: emd.cpp:474
int L
L.
Definition: emd.h:74
Class for (basis set independent) normalized radial wfs.
Definition: emd.h:41
double eval(double p) const
Evaluate density at p.
Definition: emd.cpp:694
int l
l value
Definition: emd.h:244
int getl() const
Get l value.
Definition: emd.cpp:62
void save(const char *fname) const
Save values of momentum density.
Definition: emd.cpp:975
int M
M.
Definition: emd.h:92
Structure for holding radial EMD.
Definition: emd.h:212
EMD(const EMDEvaluator *poseval, const EMDEvaluator *negeval, double Nel, int l, int m)
Constructor.
Definition: emd.cpp:670
void get_total_coupling(size_t ig, size_t jg, double p, std::vector< total_coupl_t > &c, std::vector< total_coupl_t > &tmp) const
Get the total coupling (incl. radial function)
Definition: emd.cpp:416
arma::mat compton_profile() const
Calculate Compton profile.
Definition: emd.cpp:1037
void add4(size_t ind)
Add 4 points at ind.
Definition: emd.cpp:722
int M
M.
Definition: emd.h:76
const EMDEvaluator * negeval
Negative evaluator.
Definition: emd.h:258
Value of radial function.
Definition: emd.h:106
arma::cx_mat P
The density matrix.
Definition: emd.h:142
void initial_fill(bool verbose=true)
Initial filling of grid.
Definition: emd.cpp:708
void fixed_fill(bool verbose=true, double h0=1e-3, double l0=3.0, double hfac=2.0, double lfac=2.0)
Definition: emd.cpp:916
int l
l
Definition: emd.h:108
std::vector< std::vector< std::complex< double > > > YLM
Spherical harmonics values, complex conjugated [Nat x Nat] [(L,M)].
Definition: emd.h:139
std::vector< std::vector< size_t > > idfuncs
Definition: emd.h:127
void distance_table(const std::vector< coords_t > &coord)
Computes the distance table.
Definition: emd.cpp:147
void complete_fill()
Fill regions where density changes by huge amounts.
void print() const
Print the evaluator.
Definition: emd.cpp:463
std::vector< std::vector< coupl_coeff_t > > cc
The coupling coefficients of the nonequivalent functions.
Definition: emd.h:129
double p
Radial momentum.
Definition: emd.h:214
Coupling list.
Definition: emd.h:88
~EMDEvaluator()
Destructor.
Definition: emd.cpp:312
std::complex< double > c
The coupling coefficient.
Definition: emd.h:79
virtual ~RadialFourier()
Destructor.
Definition: emd.cpp:59
int m
m value
Definition: emd.h:246
int l
l value
Definition: emd.h:44
void compton_profile_interp(const char *interp) const
Save Compton profile in interpolated form.
Definition: emd.cpp:1083
void get_coupling(size_t ig, size_t jg, int l, int lp, std::vector< total_coupl_t > &c) const
Get the coupling constants for L=|l-lp|, ..., l+lp.
Definition: emd.cpp:361
RadialFourier(int l)
Constructor.
Definition: emd.cpp:55
Radial EMD evaluator.
Definition: emd.h:122
arma::mat moments() const
Calculate moments of momentum density.
Definition: emd.cpp:982
Coupling coefficient.
Definition: emd.h:65
~EMD()
Destructor.
Definition: emd.cpp:691
int Lmax
Maximum value of L.
Definition: emd.h:145
EMDEvaluator()
Dummy constructor.
Definition: emd.cpp:128
std::complex< double > poscoef
Coefficient.
Definition: emd.h:255
std::complex< double > negcoef
Coefficient.
Definition: emd.h:260
int L
L.
Definition: emd.h:90