NESSi  v1.0.2
The NonEquilibrium Systems Simulation Library
cntr_herm_matrix_decl.hpp
Go to the documentation of this file.
1 #ifndef CNTR_HERM_MATRIX_DECL_H
2 #define CNTR_HERM_MATRIX_DECL_H
3 
5 
6 namespace cntr {
7 
8 template <typename T> class function;
9 template <typename T> class herm_matrix_timestep;
10 
11 template <typename T>
44 class herm_matrix {
45  public:
46  typedef std::complex<T> cplx;
47  typedef T scalar_type;
48 
49  /* construction, destruction */
50  herm_matrix();
51  ~herm_matrix();
52  herm_matrix(int nt, int ntau, int size1 = 1, int sig = -1);
53  herm_matrix(int nt, int ntau, int size1, int size2, int sig);
54  herm_matrix(const herm_matrix &g);
56 #if __cplusplus >= 201103L
57  herm_matrix(herm_matrix &&g) noexcept;
58  herm_matrix &operator=(herm_matrix &&g) noexcept;
59 #endif
60  /* resize */
61  void resize_discard(int nt, int ntau, int size1);
62  void resize_nt(int nt);
63  void resize(int nt, int ntau, int size1);
64  void clear(void);
65  /* access size etc ... */
67  int element_size(void) const { return element_size_; }
68  int size1(void) const { return size1_; }
69  int size2(void) const { return size2_; }
70  int ntau(void) const { return ntau_; }
71  int nt(void) const { return nt_; }
72  int sig(void) const { return sig_; }
73  void set_sig(int sig) { sig_ = sig; }
74  /* conversion from other types */
75  // herm_matrix<T> & herm_matrix(const matrix<T> &g1);
76  // herm_matrix<T> & herm_matrix(const scalar<T> &g1);
77  // raw pointer to elements ... to be used with care
79  inline cplx *lesptr(int i, int j);
81  inline cplx *retptr(int i, int j);
83  inline cplx *tvptr(int i, int j);
85  inline cplx *matptr(int i);
87  inline const cplx *lesptr(int i, int j) const;
89  inline const cplx *retptr(int i, int j) const;
91  inline const cplx *tvptr(int i, int j) const;
93  inline const cplx *matptr(int i) const;
94  // reading basic and derived elements to any Matrix type
96  template <class Matrix>
97  void get_les(int i, int j, Matrix &M) const;
99  template <class Matrix>
100  void get_gtr(int i, int j, Matrix &M) const;
102  template <class Matrix>
103  void get_ret(int i, int j, Matrix &M) const;
105  template <class Matrix>
106  void get_vt(int i, int j, Matrix &M) const;
108  template <class Matrix>
109  void get_tv(int i, int j, Matrix &M) const;
111  template <class Matrix>
112  void get_mat(int i, Matrix &M) const;
114  template <class Matrix>
115  void get_matminus(int i, Matrix &M) const;
116  // reading complex numbers:
117  // these will adress only (0,0) element for dim>1:
119  inline void get_les(int i, int j, cplx &x) const;
121  inline void get_gtr(int i, int j, cplx &x) const;
123  inline void get_ret(int i, int j, cplx &x) const;
125  inline void get_vt(int i, int j, cplx &x) const;
127  inline void get_tv(int i, int j, cplx &x) const;
129  inline void get_mat(int i, cplx &x) const;
131  inline void get_matminus(int i, cplx &x) const;
133  cplx density_matrix(int i) const;
134  // should work with eigen
135  template <class Matrix>
136  void density_matrix(int tstp, Matrix &M) const;
137 
138  // writing basic elements:
140  template <class Matrix>
141  void set_les(int i, int j, Matrix &M);
143  template <class Matrix>
144  void set_ret(int i, int j, Matrix &M);
146  template <class Matrix>
147  void set_tv(int i, int j, Matrix &M);
149  template <class Matrix>
150  void set_mat(int i, Matrix &M);
152  template<class Matrix> void set_mat_real(int i,Matrix &M);
154  void set_mat_herm(int i);
156  void set_mat_herm(void);
158  inline void set_les(int i, int j, cplx x);
160  inline void set_ret(int i, int j, cplx x);
162  inline void set_tv(int i, int j, cplx x);
164  inline void set_mat(int i, cplx x);
165  // INPUT/OUTPUT
166  void print_to_file(const char *file, int precision = 16) const;
167  void read_from_file(const char *file);
168  // read up to timestep nt1: NO RESIZE
169  void read_from_file(int nt1, const char *file);
170 // HDF5 I/O
171 #if CNTR_USE_HDF5 == 1
172  void write_to_hdf5(hid_t group_id);
173  void write_to_hdf5(hid_t group_id, const char *groupname);
174  void write_to_hdf5(const char *filename, const char *groupname);
175  void read_from_hdf5(hid_t group_id);
176  void read_from_hdf5(hid_t group_id, const char *groupname);
177  void read_from_hdf5(const char *filename, const char *groupname);
178  // read up to timestep nt1: NO RESIZE
179  void read_from_hdf5(int nt1, hid_t group_id);
180  void read_from_hdf5(int nt1, hid_t group_id, const char *groupname);
181  void read_from_hdf5(int nt1, const char *filename, const char *groupname);
182  // writing some other compressed format: just a collection of timeslices
183  // -1,0,dt,2*dt,...
184  void write_to_hdf5_slices(hid_t group_id, int dt);
185  void write_to_hdf5_slices(hid_t group_id, const char *groupname, int dt);
186  void write_to_hdf5_slices(const char *filename, const char *groupname,
187  int dt);
188  void write_to_hdf5_tavtrel(hid_t group_id, int dt);
189  void write_to_hdf5_tavtrel(hid_t group_id, const char *groupname, int dt);
190  void write_to_hdf5_tavtrel(const char *filename, const char *groupname,
191  int dt);
192 #endif
193  // simple operations on the timesteps (n<0: Matsubara)
194  void set_timestep_zero(int tstp);
195  void set_timestep(int tstp, herm_matrix &g1);
196  void set_timestep(int tstp, herm_matrix_timestep<T> &timestep);
197  void get_timestep(int tstp, herm_matrix_timestep<T> &timestep) const;
198  void get_timestep(int tstp, herm_matrix<T> &timestep) const;
199  void incr_timestep(int tstp, herm_matrix_timestep<T> &timestep,
200  cplx alpha);
201  void incr_timestep(int tstp, herm_matrix_timestep<T> &timestep);
202  void incr_timestep(int tstp, herm_matrix<T> &g, cplx alpha);
203  void incr_timestep(int tstp, herm_matrix<T> &g);
204  // check what this is actually doing
205  void incr(herm_matrix<T> &g, cplx alpha);
206  void incr(herm_matrix<T> &g);
207 
208  void set_matrixelement(int tstp, int i1, int i2,
209  herm_matrix_timestep<T> &timestep);
210  void set_matrixelement(int tstp, int i1, int i2, herm_matrix<T> &g);
211  // this(i1,i2) <- g(j1,j2)
212  void set_matrixelement(int tstp, int i1, int i2,
213  herm_matrix_timestep<T> &g, int j1, int j2);
214  void set_matrixelement(int tstp, int i1, int i2, herm_matrix<T> &g,
215  int j1, int j2);
216  void set_matrixelement(int i1,int i2,herm_matrix<T> &g);
217  void set_matrixelement(int i1,int i2,herm_matrix<T> &g,int j1,int j2);
218 
219  void set_submatrix(std::vector<int> &i1,std::vector<int> &i2,
220  herm_matrix<T> &g,std::vector<int> &j1,std::vector<int> &j2);
221  void set_submatrix(int tstp, std::vector<int> &i1,std::vector<int> &i2,
222  herm_matrix<T> &g,std::vector<int> &j1,std::vector<int> &j2);
223 
224  // multiplication with function
225 
226  // check if multiplication with pointer can be removed (there might be
227  // internal dependencies)
229  void
230  left_multiply(int tstp, cplx *f0, cplx *ft,
231  T weight = 1.0); // version with function object is safer
233  void
234  right_multiply(int tstp, cplx *f0, cplx *ft,
235  T weight = 1.0); // version with function object is safer
236  void left_multiply(int tstp, function<T> &ft, T weight = 1.0);
237  // - check why underscore (private routine?)
238  // - check what std::integral_constant is doing
239  // - add interface to multiply with a matrix from left/right
241  template <int SIZE>
242  void left_multiply_(int tstp, function<T> &ft, T weight,
243  std::integral_constant<int, SIZE>);
244  void right_multiply(int tstp, function<T> &ft, T weight = 1.0);
246  template <int SIZE>
247  void right_multiply_(int tstp, function<T> &ft, T weight,
248  std::integral_constant<int, SIZE>);
250  int tstp, function<T> &ft,
251  T weight = 1.0); // G(t,t') -> weight*f(t)^dag * G(t,t')
253  int tstp, function<T> &ft,
254  T weight = 1.0); // G(t,t') -> weight*G(t,t') * f(t')^dag
255  void smul(int tstp, T weight);
256  void smul(int tstp, cplx weight);
257 // MPI UTILS
258 #if CNTR_USE_MPI == 1
259  void Reduce_timestep(int tstp, int root);
260  void Bcast_timestep(int tstp, int root);
261  void Send_timestep(int tstp, int dest, int tag);
262  void Recv_timestep(int tstp, int root, int tag);
263 #endif
264  private:
265  int les_offset(int t, int t1) const;
266  int ret_offset(int t, int t1) const;
267  int tv_offset(int t, int tau) const;
268  int mat_offset(int tau) const;
269 
270  private:
272 
273  cplx *les_;
275 
276  cplx *ret_;
278 
279  cplx *tv_;
281 
282  cplx *mat_;
284 
285  int nt_;
287 
288  int ntau_;
290 
291  int size1_;
293 
294  int size2_;
296 
297  int element_size_;
299 
300  int sig_; // Bose = +1, Fermi =-1
301 };
302 
303 } // namespace cntr
304 
305 #endif // CNTR_HERM_MATRIX_DECL_H
void clear(void)
Sets all values to zero.
void get_timestep(int tstp, herm_matrix_timestep< T > &timestep) const
std::complex< T > cplx
Class herm_matrix_timestep deals with contour objects at a particular timestep .
void resize_nt(int nt)
Resizes herm_matrix object with respect to the number of time points nt. Internal routine; see top-l...
void resize_discard(int nt, int ntau, int size1)
Discards the herm_matrix object and resize with respect to the number of time points nt...
void left_multiply_hermconj(int tstp, function< T > &ft, T weight=1.0)
Left-multiplies the herm_matrix with the hermitian conjugate of a contour function.
void smul(int tstp, T weight)
Multiplies all component of the herm_matrix with a real scalar.
void right_multiply_hermconj(int tstp, function< T > &ft, T weight=1.0)
Right-multiplies the herm_matrix with the hermitian conjugate of a contour function.
void set_timestep_zero(int tstp)
Sets all components at time step tstp to zero.
void incr_timestep(int tstp, herm_matrix_timestep< T > &timestep, cplx alpha)
Adds a herm_matrix_timestep with given weight to the herm_matrix.
void set_matrixelement(int tstp, int i1, int i2, herm_matrix_timestep< T > &timestep)
Sets given matrix elements of all components at time step tstp to the components of a given herm_mat...
Class function for objects with time on real axis.
void Recv_timestep(int tstp, int root, int tag)
Recevies the herm_matrix at a given time step from a specific task.
void Send_timestep(int tstp, int dest, int tag)
Sends the herm_matrix at a given time step to a specific task.
void set_submatrix(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 all time steps to a (sub-) matrix of a given herm_mat...
Class herm_matrix for two-time contour objects with hermitian symmetry.
void write_to_hdf5(hid_t group_id)
Stores herm_matrix to a given group in HDF5 format.
void print_to_file(const char *file, int precision=16) const
Saves the herm_matrix to file in text format.
void write_to_hdf5_tavtrel(hid_t group_id, int dt)
Stores greater and lesser components of herm_matrix average-relative time representation to a given ...
void Reduce_timestep(int tstp, int root)
MPI reduce for the herm_matrix at a given time step.
void write_to_hdf5_slices(hid_t group_id, int dt)
Stores time slices of herm_matrix with a given interval to a HDF5 group handle. ...
void set_timestep(int tstp, herm_matrix &g1)
Sets all components at time step tstp to the components of a given herm_matrix. ...
void read_from_hdf5(hid_t group_id)
Reads herm_matrix from a given HDF5 group handle.
herm_matrix & operator=(const herm_matrix &g)
void Bcast_timestep(int tstp, int root)
Broadcasts the herm_matrix at a given time step to all tasks.
void resize(int nt, int ntau, int size1)
Resizes herm_matrix object with respect to the number of time points nt, points on the Matsubara bra...
void incr(herm_matrix< T > &g, cplx alpha)
Adds a herm_matrix with given weight to the herm_matrix at all time steps.
void read_from_file(const char *file)
Reads the herm_matrix from file in text format.