Pyxaid2
 All Classes
ElectronicStructure.h
1 /***********************************************************
2  * Copyright (C) 2013 Alexey V. Akimov
3  * This file is distributed under the terms of the
4  * GNU General Public License as published by the
5  * Free Software Foundation; either version 3 of the
6  * License, or (at your option) any later version.
7  * http://www.gnu.org/copyleft/gpl.txt
8 ***********************************************************/
9 
10 #ifndef ElectronicStructure_H
11 #define ElectronicStructure_H
12 
13 //#include "matrix.h"
14 #include "state.h"
15 #include "units_pyxaid.h"
16 #include "InputStructure.h"
17 #include <boost/python.hpp>
18 using namespace boost::python;
19 using namespace std;
20 
21 #include "liblibra_core.h"
22 using namespace liblibra::librandom;
23 using namespace liblibra::liblinalg;
24 
25 
26 class ElectronicStructure{
27 
28  // For DISH
29  void update_decoherence_times(CMATRIX& rates);
30  void project_out(int i);
31  void decohere(int i);
32 
33  // For integrators
34  void rot1(double phi,int i,int j);
35  void rot2(double phi,int i,int j);
36  void rot(complex<double> Hij,double dt,int i,int j);
37  void phase(complex<double> Hii,double dt,int i);
38 
39 public:
40 
41  //=========== Members ===============
42  int num_states; // number of adiabatic states
43  int curr_state; // current adiabatic state
44 
45  // Wavefunction
46  CMATRIX* Ccurr;
47  CMATRIX* Cprev;
48  CMATRIX* Cnext;
49  CMATRIX* A; // density matrix - populations and coherences
50 
51  // Hamiltonian: meaning of "current" and "next" depends on the algorithm to solve TD-SE
52  CMATRIX* Hcurr; // current Hamiltonian Hij = H[i*num_states+j]
53  CMATRIX* Hprev; // Hamiltonian on previous nuclear time step ---
54  CMATRIX* Hnext; // Hamiltonian on next nuclear time step ---
55  CMATRIX* dHdt; // slope of Hamiltonian - for interpolation scheme ---
56 
57  CMATRIX* Hprimex;
58  CMATRIX* Hprimey;
59  CMATRIX* Hprimez;
60 
61  vector<double> g; // num_states x num_states matrix, reshaped in 1D array
62 
63  // DISH variables:
64  vector<double> tau_m; // times since last decoherence even for all PES (actually rates, that is inverse times)
65  vector<double> t_m; // time counters for each PES
66 
67 
68 
69  //========== Methods ================
70  // Constructors
71  ElectronicStructure(int n){
72  num_states = n;
73 
74  complex<double> tmp(0.0,0.0);
75  vector<double> tmp1(n,0.0);
76  Ccurr = new CMATRIX(n,1); *Ccurr = tmp;
77  Cprev = new CMATRIX(n,1); *Cprev = tmp;
78  Cnext = new CMATRIX(n,1); *Cnext = tmp;
79 
80  g = std::vector<double>(n*n,0.0); // g[i*n+j] ~=g[i][j] - probability of i->j transition
81 
82  A = new CMATRIX(n,n); *A = tmp;
83 
84  Hcurr = new CMATRIX(n,n); *Hcurr = tmp;
85  Hprev = new CMATRIX(n,n); *Hprev = tmp;
86  Hnext = new CMATRIX(n,n); *Hnext = tmp;
87  dHdt = new CMATRIX(n,n); *dHdt = tmp;
88 
89  Hprimex = new CMATRIX(n,n); *Hprimex = tmp;
90  Hprimey = new CMATRIX(n,n); *Hprimey = tmp;
91  Hprimez = new CMATRIX(n,n); *Hprimez = tmp;
92 
93 
94  tau_m = std::vector<double>(n,0.0);
95  t_m = std::vector<double>(n,0.0);
96  }
97 
98  ElectronicStructure(const ElectronicStructure& es){ // Copy constructor
99  int n = es.num_states;
100  num_states = n;
101  curr_state = es.curr_state;
102 
103  complex<double> tmp(0.0,0.0);
104  vector<double> tmp1(n,0.0);
105  Ccurr = new CMATRIX(n,1);
106  Cprev = new CMATRIX(n,1);
107  Cnext = new CMATRIX(n,1);
108 
109  g = std::vector<double>(n*n,0.0); // g[i*n+j] ~=g[i][j] - probability of i->j transition
110 
111  A = new CMATRIX(n,n);
112 
113  Hcurr = new CMATRIX(n,n);
114  Hprev = new CMATRIX(n,n);
115  Hnext = new CMATRIX(n,n);
116  dHdt = new CMATRIX(n,n);
117 
118  Hprimex = new CMATRIX(n,n);
119  Hprimey = new CMATRIX(n,n);
120  Hprimez = new CMATRIX(n,n);
121 
122 
123  tau_m = es.tau_m;
124  t_m = es.t_m;
125 
126  *Ccurr = *es.Ccurr; *Cprev = *es.Cprev; *Cnext = *es.Cnext;
127  g = es.g; *A = *es.A;
128  *Hcurr = *es.Hcurr; *Hprev = *es.Hprev; *Hnext = *es.Hnext; *dHdt = *es.dHdt;
129  *Hprimex = *es.Hprimex; *Hprimey = *es.Hprimey; *Hprimez = *es.Hprimez;
130  }
131  // Destructor
132 
133  ~ElectronicStructure(){
134  if(g.size()>0) {g.clear();}
135  if(Ccurr!=NULL) { delete Ccurr;}// Ccurr = NULL;}
136  if(Cprev!=NULL) {delete Cprev;}// Cprev = NULL;}
137  if(Cnext!=NULL) {delete Cnext;}// Cnext = NULL;}
138  if(A!=NULL) { delete A;} // A = NULL;}
139  if(Hcurr!=NULL){ delete Hcurr;} // Hcurr = NULL;}
140  if(Hprev!=NULL){ delete Hprev;}// Hprev = NULL;}
141  if(Hnext!=NULL){ delete Hnext;}// Hnext = NULL;}
142  if(dHdt!=NULL){ delete dHdt;} // dHdt = NULL;}
143  if(Hprimex!=NULL){ delete Hprimex; }
144  if(Hprimey!=NULL){ delete Hprimey; }
145  if(Hprimez!=NULL){ delete Hprimez; }
146  if(tau_m.size()>0){ tau_m.clear(); }
147  if(t_m.size()>0){ t_m.clear(); }
148  }
149 
150 
151  ElectronicStructure operator=(ElectronicStructure es){
152  num_states = es.num_states;
153  curr_state = es.curr_state;
154  *Ccurr = *es.Ccurr; *Cprev = *es.Cprev; *Cnext = *es.Cnext;
155  g = es.g; *A = *es.A;
156  *Hcurr = *es.Hcurr; *Hprev = *es.Hprev; *Hnext = *es.Hnext;
157  *Hprimex = *es.Hprimex; *Hprimey = *es.Hprimey; *Hprimez = *es.Hprimez;
158  *dHdt = *es.dHdt;
159  tau_m = es.tau_m; t_m = es.t_m;
160  return *this;
161  }
162 
163  ElectronicStructure operator<<(ElectronicStructure es){
164 // This is basically the same as operator=, but keeps parameters the same
165 // copied only wavefunction and state
166 
167  num_states = es.num_states;
168  curr_state = es.curr_state;
169  *Cprev = *es.Cprev;
170  *Ccurr = *es.Ccurr;
171  *Cnext = *es.Cnext;
172 
173 // Also need to update the DISH timescales
174  tau_m = es.tau_m;
175  t_m = es.t_m;
176 
177  return *this;
178  }
179 
180 
181 
182  // Other methods
183  void set_state(int indx){
184  for(int i=0;i<num_states;i++){
185  if(i==indx){ Ccurr->M[i] = complex<double>(1.0,0.0); } else{ Ccurr->M[i] = 0.0; }
186  }
187  curr_state = indx;
188  }
189  double energy(); // calculate the total energy
190  double norm(); // calculate the norm of the wavefunction
191 
192  void update_populations(); // update matrix A from Ccurr
193  void update_hop_prob(double dt,int is_boltz_flag,double Temp, CMATRIX& Ef);
194 
195 
196  void update_hop_prob_fssh(double dt,int is_boltz_flag,double Temp, CMATRIX& Ef,double Ex, CMATRIX&);
197  void update_hop_prob_mssh(double dt,int is_boltz_flag,double Temp, CMATRIX& Ef,double Ex, CMATRIX&);
198  void update_hop_prob_gfsh(double dt,int is_boltz_flag,double Temp, CMATRIX& Ef,double Ex, CMATRIX&);
199  void init_hop_prob1();
200 
201  void check_decoherence(double dt,int boltz_flag,double Temp, CMATRIX& rates, Random& rnd); // practically DISH correction
202 
203  void propagate_coefficients(double dt, CMATRIX& Ef); // Trotter factorization
204  void propagate_coefficients(double dt, CMATRIX& Ef, CMATRIX&); // Trotter factorization with purostat
205  void propagate_coefficients1(double dt,int opt, CMATRIX& Ef); // Finite difference
206  void propagate_coefficients2(double dt, CMATRIX& Ef); // "Exact"
207 
208 };
209 
210 
211 #endif // ElectronicStructure_H
212 
ElectronicStructure
Definition: ElectronicStructure.h:26
Generated on Sat Aug 26 2017 08:29:34 for Pyxaid2 by   doxygen 1.8.6