NESSi  v1.0.2
The NonEquilibrium Systems Simulation Library
cntr_herm_matrix_timestep_decl.hpp
Go to the documentation of this file.
1 #ifndef CNTR_HERM_MATRIX_TIMESTEP_DECL_H
2 #define CNTR_HERM_MATRIX_TIMESTEP_DECL_H
3 
5 
6 namespace cntr {
7 
8 template <typename T> class function;
10 template <typename T> class cyclic_timestep;
11 template <typename T> class herm_matrix;
12 template <typename T> class herm_matrix_timestep_view;
13 
14 // NOTE: the bose/fermi sign for the herm_matrix_timestep (needed to compute
15 // vt from tv etc.)
16 // is currently not consistently treated ... safe only for fermionic GF
17 template <typename T>
33 class herm_matrix_timestep {
34  public:
35  typedef std::complex<T> cplx;
36  typedef T scalar_type;
37  /* construction, destruction */
40  herm_matrix_timestep(int tstp,int ntau,int size1=1);
41  herm_matrix_timestep(int tstp,int ntau,int size1,int sig);
42  herm_matrix_timestep(int tstp,int ntau,int size1,int size2,int sig);
45 #if __cplusplus >= 201103L
48 #endif
49  /* resize */
50  void resize(int nt, int ntau, int size1);
51  void clear(void);
52  /* access size etc ... */
53  int size1(void) const { return size1_; }
54  int size2(void) const { return size2_; }
55  int tstp(void) const { return tstp_; }
56  int ntau(void) const { return ntau_; }
57  int sig(void) const { return sig_; }
58  void get_matrixelement(int i1, int i2, herm_matrix<T> &g);
59  void set_timestep(cyclic_timestep<T> &timestep,
60  cyclic_timestep<T> &timestep_cc, int sig);
61  // Here timestep is tstp of pseudo ("ret=gtr")
62  void set_timestep_pseudo(cyclic_timestep<T> &timestep,
63  cyclic_timestep<T> &timestep_cc, int sig);
64  // assumes hermitan:
65  void set_timestep(cyclic_timestep<T> &timestep, int sig);
66  void set_timestep_pseudo(cyclic_timestep<T> &timestep, int sig);
67  inline cplx *retptr(int i) { return data_ + i * element_size_; }
68  inline cplx *tvptr(int i) {
69  return data_ + (tstp_ + 1 + i) * element_size_;
70  }
71  inline cplx *lesptr(int i) {
72  return data_ + (tstp_ + 1 + ntau_ + 1 + i) * element_size_;
73  };
74  // preferred interfaces
75 
76  void set_timestep_zero(int tstp);
77 
78  void set_timestep(int tstp, herm_matrix<T> &g1);
80 
82  template <class Matrix> void set_ret(int i, int j, Matrix &M);
84  template <class Matrix> void set_les(int i, int j, Matrix &M);
86  template <class Matrix> void set_tv(int i, int j, Matrix &M);
88  template <class Matrix> void set_mat(int i, Matrix &M);
89 
91  template <class Matrix> void get_les(int i, int j, Matrix &M);
93  template <class Matrix> void get_ret(int i, int j, Matrix &M);
95  template <class Matrix> void get_tv(int i, int j, Matrix &M);
97  template <class Matrix> void get_mat(int i, Matrix &M);
99  template <class Matrix> void get_matminus(int i, Matrix &M);
100 
101  template <class Matrix> void density_matrix(int tstp, Matrix &M);
102 
103  // reading complex numbers:
104  // these will adress only (0,0) element for dim>1:
105 
107  inline void set_les(int i, int j, cplx &x);
109  inline void set_ret(int i, int j, cplx &x);
111  inline void set_tv(int i, int j, cplx &x);
113  inline void set_mat(int i, cplx &x);
114 
116  inline void get_les(int i, int j, cplx &x);
118  inline void get_ret(int i, int j, cplx &x);
120  inline void get_tv(int i, int j, cplx &x);
122  inline void get_mat(int i, cplx &x);
124  inline void get_matminus(int i, cplx &x);
125 
126  inline void density_matrix(int tstp, cplx &rho);
127 
128  // legacy interfaces
129 
131  template <class Matrix> void set_ret(int j, Matrix &M);
133  template <class Matrix> void set_les(int j, Matrix &M);
135  template <class Matrix> void set_tv(int j, Matrix &M);
136 
139  template <class Matrix> void get_les_t_tstp(int i, Matrix &M);
141  template <class Matrix> void get_les_tstp_t(int i, Matrix &M);
143  template <class Matrix> void get_ret_tstp_t(int j, Matrix &M);
145  template <class Matrix> void get_ret_t_tstp(int i, Matrix &M);
147  template <class Matrix> void get_tv(int j, Matrix &M);
149  template <class Matrix> void get_vt(int i, Matrix &M, int sig);
151  template <class Matrix> void get_matminus(int i, Matrix &M, int sig);
153  template <class Matrix> void get_gtr_tstp_t(int i, Matrix &M);
155  template <class Matrix> void get_gtr_t_tstp(int i, Matrix &M);
157  inline void get_vt(int i, int j, cplx &x);
158 
160  cplx density_matrix(int tstp);
162  template <class Matrix> void density_matrix(Matrix &M);
163 
165  inline cplx *matptr(int i) { return data_ + i * element_size_; }
166  // multiplication with function
168  void left_multiply(cplx *f0, cplx *ft, T weight = 1.0);
170  void right_multiply(cplx *f0, cplx *ft, T weight = 1.0);
171  void left_multiply(function<T> &ft, T weight = 1.0);
172  void right_multiply(function<T> &ft, T weight = 1.0);
173  void left_multiply(int tstp, function<T> &ft, T weight = 1.0);
174  void right_multiply(int tstp, function<T> &ft, T weight = 1.0);
175  void left_multiply_hermconj(int tstp, function<T> &ft, T weight = 1.0);
176  void right_multiply_hermconj(int tstp, function<T> &ft, T weight = 1.0);
177  void incr(herm_matrix_timestep<T> &g1, T weight);
178  void incr(herm_matrix<T> &g, T weight = 1.0);
179  void incr_timestep(int tstp, herm_matrix_timestep<T> &g1, T weight);
180  void incr_timestep(int tstp, herm_matrix<T> &g, T weight = 1.0);
181  void smul(int tstp, T weight);
182  void smul(T weight);
183  void set_matrixelement(int i1, int i2, herm_matrix_timestep_view<T> &g,
184  int j1, int j2);
185  void set_matrixelement(int tstp, int i1, int i2, herm_matrix_timestep_view<T> &g,
186  int j1, int j2);
187  void set_matrixelement(int i1, int i2, herm_matrix_timestep<T> &g, int j1,
188  int j2);
189  void set_matrixelement(int tstp, int i1, int i2, herm_matrix_timestep<T> &g, int j1,
190  int j2);
191  void set_matrixelement(int i1, int i2, herm_matrix<T> &g, int j1, int j2);
192  void set_matrixelement(int tstp, int i1, int i2, herm_matrix<T> &g, int j1, int j2);
193  void set_submatrix(int tstp, std::vector<int> &i1, std::vector<int> &i2,herm_matrix<T> &g,
194  std::vector<int> &j1,std::vector<int> &j2);
195  void set_submatrix(int tstp, std::vector<int> &i1, std::vector<int> &i2,herm_matrix_timestep<T> &g,
196  std::vector<int> &j1,std::vector<int> &j2);
197 // MPI UTILS
198 #if CNTR_USE_MPI == 1
199  // preferred interfaces
200  void Reduce_timestep(int tstp, int root);
201  void Bcast_timestep(int tstp, int root);
202  void Send_timestep(int tstp, int dest, int tag);
203  void Recv_timestep(int tstp, int root, int tag);
204 
205  // legacy interfaces
207  void Reduce_timestep(int root);
209  void Bcast_timestep(int tstp, int ntau, int size1, int root);
211  void Send_timestep(int tstp, int ntau, int size1, int dest, int tag);
213  void Recv_timestep(int tstp, int ntau, int size1, int root, int tag);
214 #endif
215 // HDF5 I/O
216 #if CNTR_USE_HDF5 == 1
217  void write_to_hdf5(hid_t group_id);
218  void write_to_hdf5(hid_t group_id, const char *groupname);
219  void write_to_hdf5(const char *filename, const char *groupname);
220  // does a resize:
221  void read_from_hdf5(hid_t group_id, const char *groupname);
222  void read_from_hdf5(hid_t group_id);
223  void read_from_hdf5(const char *filename, const char *groupname);
224 #endif
225 
227  cplx *data_;
229 
230  int tstp_;
232 
233  int ntau_;
235 
236  int size1_;
238 
239  int size2_;
241 
242  int element_size_;
244 
245  int total_size_;
247 
248  int sig_;
249 };
250 
251 } // namespace cntr
252 
253 #endif // CNTR_HERM_MATRIX_TIMESTEP_DECL_H
herm_matrix_timestep & operator=(const herm_matrix_timestep &g)
void smul(int tstp, T weight)
Multiply herm_matrix_timestep with scalar weight.
Class herm_matrix_timestep_view serves for interfacing with class herm_matrix_timestep without copyi...
void density_matrix(int tstp, Matrix &M)
Returns the density matrix at given time step.
void Recv_timestep(int tstp, int root, int tag)
Recevies the herm_matrix at a given time step from a specific task.
Class herm_matrix_timestep deals with contour objects at a particular timestep .
void right_multiply_hermconj(int tstp, function< T > &ft, T weight=1.0)
Right-multiplies the herm_matrix_timestep with the hermitian conjugate of a contour function...
void set_timestep_zero(int tstp)
Sets all components at time step tstp to zero.
void resize(int nt, int ntau, int size1)
Resizes herm_matrix_timestep object with respect to the number of points on the Matsubara branch and...
void set_timestep(cyclic_timestep< T > &timestep, cyclic_timestep< T > &timestep_cc, int sig)
Class function for objects with time on real axis.
void Send_timestep(int tstp, int dest, int tag)
Sends the herm_matrix_timestep at a given time step to a specific task.
void write_to_hdf5(hid_t group_id)
Write herm_matrix_timestep to hdf5 group given by group_id
void set_matrixelement(int i1, int i2, herm_matrix_timestep_view< T > &g, int j1, int j2)
Set the matrix element of for each component from
Class herm_matrix for two-time contour objects with hermitian symmetry.
void set_submatrix(int tstp, std::vector< int > &i1, std::vector< int > &i2, herm_matrix< T > &g, std::vector< int > &j1, std::vector< int > &j2)
Sets a (sub-) matrix of this contour object at a given time step to a (sub-) matrix of a given herm_...
void get_matrixelement(int i1, int i2, herm_matrix< T > &g)
void left_multiply_hermconj(int tstp, function< T > &ft, T weight=1.0)
Left-multiplies the herm_matrix_timestep with the hermitian conjugate of a contour function...
void set_timestep_pseudo(cyclic_timestep< T > &timestep, cyclic_timestep< T > &timestep_cc, int sig)
void clear(void)
Sets all values to zero.
void Bcast_timestep(int tstp, int root)
Broadcasts the herm_matrix_timestep at a given time step to all ranks.
void read_from_hdf5(hid_t group_id, const char *groupname)
Read herm_matrix_timestep from hdf5 group given by group_id and group name groupname ...
void incr_timestep(int tstp, herm_matrix_timestep< T > &g1, T weight)
Increase the value of the herm_matrix_timestep by weight
void Reduce_timestep(int tstp, int root)
MPI reduce for the herm_matrix_timestep
void incr(herm_matrix_timestep< T > &g1, T weight)
Increase the value of the herm_matrix_timestep by weight