1 //////////////////////////////////////////////////////////////////////////////////////
2 // This file is distributed under the University of Illinois/NCSA Open Source License.
3 // See LICENSE file in top directory for details.
4 //
5 // Copyright (c) 2016 Jeongnim Kim and QMCPACK developers.
6 //
7 // File developed by: Jaron T. Krogel, krogeljt@ornl.gov, Oak Ridge National Laboratory
8 //                    Mark A. Berrill, berrillma@ornl.gov, Oak Ridge National Laboratory
9 //
10 // File created by: Jaron T. Krogel, krogeljt@ornl.gov, Oak Ridge National Laboratory
11 //////////////////////////////////////////////////////////////////////////////////////
12 
13 
14 #ifndef QMCPLUSPLUS_ONE_BODY_DENSITY_MATRICES_H
15 #define QMCPLUSPLUS_ONE_BODY_DENSITY_MATRICES_H
16 
17 #include "QMCHamiltonians/OperatorBase.h"
18 #include "QMCWaveFunctions/CompositeSPOSet.h"
19 #include "ParticleBase/RandomSeqGenerator.h"
20 #include "QMCWaveFunctions/WaveFunctionFactory.h"
21 
22 namespace qmcplusplus
23 {
24 class DensityMatrices1B : public OperatorBase
25 {
26 protected:
27   enum DMTimers
28   {
29     DM_eval,
30     DM_gen_samples,
31     DM_gen_sample_basis,
32     DM_gen_sample_ratios,
33     DM_gen_particle_basis,
34     DM_matrix_products,
35     DM_accumulate,
36   };
37 
38   TimerList_t timers;
39 
40 public:
41   enum
42   {
43     DIM = OHMMS_DIM
44   };
45 
46   typedef ValueType Value_t;
47   typedef GradType Grad_t;
48   typedef SPOSet::ValueVector_t ValueVector_t;
49   typedef SPOSet::GradVector_t GradVector_t;
50   typedef ParticleSet::ParticleLayout_t Lattice_t;
51   typedef Vector<Value_t> Vector_t;
52   typedef Matrix<Value_t> Matrix_t;
53   typedef std::vector<PosType> pts_t;
54   typedef std::vector<RealType> dens_t;
55 
56   enum integrators
57   {
58     uniform_grid = 0,
59     uniform,
60     density,
61     no_integrator
62   };
63 
64   enum evaluators
65   {
66     loop = 0,
67     matrix,
68     no_evaluator
69   };
70 
71   enum samplings
72   {
73     volume_based = 0,
74     metropolis,
75     no_sampling
76   };
77 
78   //data members
79   bool energy_mat;
80   CompositeSPOSet basis_functions;
81   ValueVector_t basis_values;
82   ValueVector_t basis_norms;
83   GradVector_t basis_gradients;
84   ValueVector_t basis_laplacians;
85   ValueVector_t integrated_values;
86   bool warmed_up;
87   std::vector<PosType> rsamples;
88   Vector<RealType> sample_weights;
89   std::vector<ValueType> psi_ratios;
90   RealType dens;
91   PosType drift;
92   int nindex;
93   int eindex;
94   Lattice_t& Lattice;
95   TrialWaveFunction& Psi;
96   ParticleSet& Pq;
97   const ParticleSet* Pc;
98   TraceSample<TraceReal>* w_trace;
99   TraceSample<TraceComp>* T_trace;
100   CombinedTraceSample<TraceReal>* Vq_trace;
101   CombinedTraceSample<TraceReal>* Vc_trace;
102   CombinedTraceSample<TraceReal>* Vqq_trace;
103   CombinedTraceSample<TraceReal>* Vqc_trace;
104   CombinedTraceSample<TraceReal>* Vcc_trace;
105   std::vector<Value_t> E_samp;
106   CombinedTraceSample<TraceReal>* E_trace;
107 
108 
109   bool initialized;
110   bool normalized;
111   bool volume_normed;
112   int basis_size;
113   int samples;
114   int nparticles;
115   int nspecies;
116   std::vector<int> species_size;
117   std::vector<std::string> species_name;
118   std::vector<Vector_t*> E_N;
119   std::vector<Matrix_t*> Phi_NB, Psi_NM, Phi_Psi_NB, N_BB, E_BB;
120   Matrix_t Phi_MB;
121   bool check_overlap;
122   bool check_derivatives;
123 
124 //#define DMCHECK
125 #ifdef DMCHECK
126   std::vector<Vector_t*> E_Ntmp;
127   std::vector<Matrix_t*> Phi_NBtmp, Psi_NMtmp, Phi_Psi_NBtmp, N_BBtmp, E_BBtmp;
128   Matrix_t Phi_MBtmp;
129 #endif
130 
131   integrators integrator;
132   evaluators evaluator;
133   samplings sampling;
134   int points;
135   RealType scale;
136   PosType center, rcorner;
137   RealType volume;
138   bool periodic;
139   int warmup;
140   RealType timestep;
141   bool use_drift;
142   int nmoves;
143   int naccepted;
144   RealType acceptance_ratio;
145   bool write_acceptance_ratio;
146   bool write_rstats;
147 
148   int ind_dims[DIM];
149   RealType metric;
150 
151   PosType rpcur;
152   PosType dpcur;
153   RealType rhocur;
154 
155   RandomGenerator_t* uniform_random;
156 
157 
158   //constructor/destructor
159   DensityMatrices1B(ParticleSet& P, TrialWaveFunction& psi, ParticleSet* Pcl, const WaveFunctionFactory& factory);
160   DensityMatrices1B(DensityMatrices1B& master, ParticleSet& P, TrialWaveFunction& psi);
161   ~DensityMatrices1B();
162 
163   //standard interface
164   OperatorBase* makeClone(ParticleSet& P, TrialWaveFunction& psi);
165   bool put(xmlNodePtr cur);
166   Return_t evaluate(ParticleSet& P);
167 
168   //optional standard interface
169   void get_required_traces(TraceManager& tm);
170   void setRandomGenerator(RandomGenerator_t* rng);
171 
172   //required for Collectables interface
173   void addObservables(PropertySetType& plist, BufferType& olist);
174   void registerCollectables(std::vector<observable_helper*>& h5desc, hid_t gid) const;
175 
176   //should be empty for Collectables interface
resetTargetParticleSet(ParticleSet & P)177   void resetTargetParticleSet(ParticleSet& P) {}
setObservables(PropertySetType & plist)178   void setObservables(PropertySetType& plist) {}
setParticlePropertyList(PropertySetType & plist,int offset)179   void setParticlePropertyList(PropertySetType& plist, int offset) {}
contribute_scalar_quantities()180   void contribute_scalar_quantities() {}
checkout_scalar_quantities(TraceManager & tm)181   void checkout_scalar_quantities(TraceManager& tm) {}
collect_scalar_quantities()182   void collect_scalar_quantities() {}
delete_scalar_quantities()183   void delete_scalar_quantities() {}
184 
185   //obsolete?
get(std::ostream & os)186   bool get(std::ostream& os) const { return false; }
187 
188   //local functions
189   //  initialization/finalization
190   void reset();
191   void set_state(xmlNodePtr cur);
192   void set_state(DensityMatrices1B& master);
193   void initialize();
194   void finalize();
195   void normalize();
196   //  printing
197   void report(const std::string& pad = "");
198   //  sample generation
199   void warmup_sampling();
200   void generate_samples(RealType weight, int steps = 0);
201   void generate_uniform_grid(RandomGenerator_t& rng);
202   void generate_uniform_samples(RandomGenerator_t& rng);
203   void generate_density_samples(bool save, int steps, RandomGenerator_t& rng);
204   void diffusion(RealType sqt, PosType& diff);
205   void density_only(const PosType& r, RealType& dens);
206   void density_drift(const PosType& r, RealType& dens, PosType& drift);
207   //  basis & wavefunction ratio matrix construction
208   void get_energies(std::vector<Vector_t*>& E_n);
209   void generate_sample_basis(Matrix_t& Phi_mb);
210   void generate_sample_ratios(std::vector<Matrix_t*> Psi_nm);
211   void generate_particle_basis(ParticleSet& P, std::vector<Matrix_t*>& Phi_nb);
212   //  basis set updates
213   void update_basis(const PosType& r);
214   void update_basis_d012(const PosType& r);
215   //  testing
216   void test_overlap();
217   void test_derivatives();
218   //  original loop implementation
219   void integrate(ParticleSet& P, int n);
220   Return_t evaluate_loop(ParticleSet& P);
221   //  matrix implementation
222   Return_t evaluate_check(ParticleSet& P);
223   Return_t evaluate_matrix(ParticleSet& P);
224 
225 
226   bool match(Value_t e1, Value_t e2, RealType tol = 1e-12);
227   bool same(Vector_t& v1, Vector_t& v2, RealType tol = 1e-6);
228   bool same(Matrix_t& m1, Matrix_t& m2, RealType tol = 1e-6);
229   void compare(const std::string& name, Vector_t& v1, Vector_t& v2, bool write = false, bool diff_only = true);
230   void compare(const std::string& name, Matrix_t& m1, Matrix_t& m2, bool write = false, bool diff_only = true);
231 
232 private:
233   /// reference to the sposet_builder_factory
234   const WaveFunctionFactory& wf_factory_;
235 };
236 
237 } // namespace qmcplusplus
238 
239 #endif
240