NESSi  v1.0.2
The NonEquilibrium Systems Simulation Library
cntr_herm_matrix_timestep_impl.hpp
Go to the documentation of this file.
1 #ifndef CNTR_HERM_MATRIX_TIMESTEP_IMPL_H
2 #define CNTR_HERM_MATRIX_TIMESTEP_IMPL_H
3 
5 #include "cntr_elements.hpp"
6 #include "cntr_function_decl.hpp"
9 
10 namespace cntr {
11 
12 /* #######################################################################################
13 #
14 # CONSTRUCTION/DESTRUCTION
15 #
16 ########################################################################################*/
17 template <typename T>
19  data_ = 0;
20  ntau_ = 0;
21  tstp_ = 0;
22  size1_ = 0;
23  size2_ = 0;
24  element_size_ = 0;
25  total_size_ = 0;
26  sig_ = -1;
27  }
28 template <typename T>
30  delete[] data_;
31  }
32 
54 template <typename T> herm_matrix_timestep<T>::herm_matrix_timestep(int tstp,int ntau,int size1){
55  int len=((tstp+1)*2+(ntau+1))*size1*size1;
56  assert(size1>=0 && tstp>=-1 && ntau>=0);
57  if(len==0) data_=0;
58  else{
59  data_ = new cplx [len];
60  memset(data_, 0, sizeof(cplx)*len);
61  }
62  size1_=size1;
63  size2_=size1;
64  element_size_=size1*size1_;
65  tstp_=tstp;
66  ntau_=ntau;
67  total_size_=len;
68  sig_=-1;
69 }
70 
95 template <typename T> herm_matrix_timestep<T>::herm_matrix_timestep(int tstp,int ntau,int size1,int size2,int sig){
96 int len=((tstp+1)*2+(ntau+1))*size1*size2;
97 assert(size1>=0 && tstp>=-1 && ntau>=0);
98 if(len==0) data_=0;
99 else{
100  data_ = new cplx [len];
101  memset(data_, 0, sizeof(cplx)*len);
102 }
103 size1_=size1;
104 size2_=size2;
105 element_size_=size1*size2;
106 tstp_=tstp;
107 ntau_=ntau;
108 total_size_=len;
109 sig_=sig;
110 }
111 
133 template <typename T> herm_matrix_timestep<T>::herm_matrix_timestep(int tstp,int ntau,int size1,int sig){
134 int len=((tstp+1)*2+(ntau+1))*size1*size1;
135 assert(size1>=0 && tstp>=-1 && ntau>=0 && sig*sig==1);
136 if(len==0) data_=0;
137 else{
138  data_ = new cplx [len];
139  memset(data_, 0, sizeof(cplx)*len);
140 }
141 size1_=size1;
142 size2_=size1;
143 element_size_=size1*size1_;
144 tstp_=tstp;
145 ntau_=ntau;
146 total_size_=len;
147 sig_=sig;
148 }
149 
167 template <typename T>
169  tstp_ = g.tstp_;
170  ntau_ = g.ntau_;
171  size1_ = g.size1_;
172  size2_ = g.size1_;
173  element_size_ = size1_ * size1_;
174  total_size_ = g.total_size_;
175  if (total_size_ > 0) {
176  data_ = new cplx[total_size_];
177  memcpy(data_, g.data_, sizeof(cplx) * total_size_);
178  } else {
179  data_ = 0;
180  }
181  sig_ = g.sig_;
182 }
183 template <typename T>
186  if (this == &g)
187  return *this;
188  if (total_size_ != g.total_size_) {
189  delete[] data_;
190  total_size_ = g.total_size_;
191  if (total_size_ > 0)
192  data_ = new cplx[total_size_];
193  else
194  data_ = 0;
195  }
196  tstp_ = g.tstp_;
197  ntau_ = g.ntau_;
198  size1_ = g.size1_;
199  size2_ = g.size1_;
200  sig_ = g.sig_;
201  element_size_ = size1_ * size1_;
202  if (total_size_ > 0)
203  memcpy(data_, g.data_, sizeof(cplx) * total_size_);
204  return *this;
205 }
206 #if __cplusplus >= 201103L
207 template <typename T>
209  herm_matrix_timestep &&g) noexcept : data_(g.data_),
210 tstp_(g.tstp_),
211 ntau_(g.ntau_),
212 size1_(g.size1_),
213 size2_(g.size2_),
214 element_size_(g.element_size_),
215 total_size_(g.total_size_),
216 sig_(g.sig_) {
217  g.data_ = nullptr;
218  g.tstp_ = 0;
219  g.ntau_ = 0;
220  g.size1_ = 0;
221  g.size2_ = 0;
222  g.element_size_ = 0;
223  g.total_size_ = 0;
224 }
225 template <typename T>
228  if (&g == this)
229  return *this;
230 
231  data_ = g.data_;
232  tstp_ = g.tstp_;
233  ntau_ = g.ntau_;
234  size1_ = g.size1_;
235  size2_ = g.size2_;
236  element_size_ = g.element_size_;
237  total_size_ = g.total_size_;
238  sig_ = g.sig_;
239 
240  g.data_ = nullptr;
241  g.tstp_ = 0;
242  g.ntau_ = 0;
243  g.size1_ = 0;
244  g.size2_ = 0;
245  g.element_size_ = 0;
246  g.total_size_ = 0;
247 
248  return *this;
249 }
250 #endif
251 
274 template <typename T>
275 void herm_matrix_timestep<T>::resize(int tstp, int ntau, int size1) {
276  int len = ((tstp + 1) * 2 + (ntau + 1)) * size1 * size1;
277  assert(ntau >= 0 && tstp >= -1 && size1 >= 0);
278  delete[] data_;
279  if (len == 0)
280  data_ = 0;
281  else {
282  data_ = new cplx[len];
283  }
284  size1_ = size1;
285  size2_ = size1;
286  element_size_ = size1_ * size1_;
287  tstp_ = tstp;
288  ntau_ = ntau;
289  total_size_ = len;
290 }
291 
303 template <typename T>
305  if (total_size_ > 0)
306  memset(data_, 0, sizeof(cplx) * total_size_);
307 }
308 
309 
328 template <typename T>
330  assert(tstp == tstp_);
331  assert(tstp >= -1);
332  if (tstp == -1) {
333  memset(matptr(0), 0, sizeof(cplx) * (ntau_ + 1) * element_size_);
334  } else {
335  memset(retptr(0), 0, sizeof(cplx) * (tstp + 1) * element_size_);
336  memset(tvptr(0), 0, sizeof(cplx) * (ntau_ + 1) * element_size_);
337  memset(lesptr(0), 0, sizeof(cplx) * (tstp + 1) * element_size_);
338  }
339 }
340 
341 
364 template <typename T>
366  assert(tstp == tstp_);
367  assert(tstp >= -1 && tstp <= g1.nt() && "tstp >= -1 && tstp <= g1.nt()");
368  assert(g1.size1() == size1_ && "g1.size1() == size1_");
369  assert(g1.ntau() == ntau_ && "g1.ntau() == ntau_");
370  if (tstp == -1) {
371  memcpy(matptr(0), g1.matptr(0), sizeof(cplx) * (ntau_ + 1) * element_size_);
372  } else {
373  memcpy(retptr(0), g1.retptr(tstp, 0),
374  sizeof(cplx) * (tstp + 1) * element_size_);
375  memcpy(tvptr(0), g1.tvptr(tstp, 0),
376  sizeof(cplx) * (ntau_ + 1) * element_size_);
377  memcpy(lesptr(0), g1.lesptr(0, tstp),
378  sizeof(cplx) * (tstp + 1) * element_size_);
379  }
380 }
381 
404 template <typename T>
406  assert(tstp == tstp_);
407  assert(tstp == g1.tstp());
408  assert(tstp >= -1 && tstp <= g1.tstp() && "tstp >= -1 && tstp <= g1.tstp()");
409  assert(g1.size1() == size1_ && "g1.size1() == size1_");
410  assert(g1.ntau() == ntau_ && "g1.ntau() == ntau_");
411  if (tstp == -1) {
412  memcpy(matptr(0), g1.matptr(0), sizeof(cplx) * (ntau_ + 1) * element_size_);
413  } else {
414  memcpy(retptr(0), g1.retptr(0),
415  sizeof(cplx) * (tstp + 1) * element_size_);
416  memcpy(tvptr(0), g1.tvptr(0),
417  sizeof(cplx) * (ntau_ + 1) * element_size_);
418  memcpy(lesptr(0), g1.lesptr(0),
419  sizeof(cplx) * (tstp + 1) * element_size_);
420  }
421 }
422 
423 
424 
426 // the following routines are not very "efficent" but sometimes simple to
427 // implement
428 #define herm_matrix_READ_ELEMENT \
429 { \
430  int r, s, dim = size1_; \
431  M.resize(dim, dim); \
432  for (r = 0; r < dim; r++) \
433  for (s = 0; s < dim; s++) \
434  M(r, s) = x[r * dim + s]; \
435  }
436 #define herm_matrix_READ_ELEMENT_MINUS_CONJ \
437  { \
438  cplx w; \
439  int r, s, dim = size1_; \
440  M.resize(dim, dim); \
441  for (r = 0; r < dim; r++) \
442  for (s = 0; s < dim; s++) { \
443  w = x[s * dim + r]; \
444  M(r, s) = std::complex<T>(-w.real(), w.imag()); \
445  } \
446  }
447 
448 
469 template <typename T> template <class Matrix>
470  void herm_matrix_timestep<T>::get_les(int i, int j, Matrix &M){
471  assert(j == tstp_ && i <= tstp_);
472  cplx *x;
473  x = lesptr(i);
474  herm_matrix_READ_ELEMENT
475  }
476 
477 
498 template <typename T>
499 template <class Matrix>
500  void herm_matrix_timestep<T>::get_les_t_tstp(int i, Matrix &M) {
501  cplx *x;
502  x = lesptr(i);
503  herm_matrix_READ_ELEMENT
504  }
505 
526 template <typename T>
527 template <class Matrix>
528  void herm_matrix_timestep<T>::get_les_tstp_t(int i, Matrix &M) {
529  cplx *x;
530  x = lesptr(i);
531  herm_matrix_READ_ELEMENT_MINUS_CONJ
532  }
533 
534 
535 
556 template <typename T>
557 template <class Matrix>
558  void herm_matrix_timestep<T>::get_ret(int i, int j, Matrix &M) {
559  assert(i == tstp_ && j <= tstp_);
560  cplx *x;
561  x = retptr(j);
562  herm_matrix_READ_ELEMENT
563  }
564 
565 
566 
587 template <typename T>
588 template <class Matrix>
589  void herm_matrix_timestep<T>::get_ret_tstp_t(int j, Matrix &M) {
590  cplx *x;
591  x = retptr(j);
592  herm_matrix_READ_ELEMENT
593  }
594 
615 template <typename T>
616 template <class Matrix>
617  void herm_matrix_timestep<T>::get_ret_t_tstp(int i, Matrix &M) {
618  cplx *x;
619  x = retptr(i);
620  herm_matrix_READ_ELEMENT_MINUS_CONJ
621  }
622 
623 
644 template <typename T>
645 template <class Matrix>
646  void herm_matrix_timestep<T>::get_tv(int i, int j, Matrix &M) {
647  assert(i == tstp_);
648  cplx *x = tvptr(j);
649  herm_matrix_READ_ELEMENT
650  }
651 
652 
673 template <typename T>
674 template <class Matrix>
675  void herm_matrix_timestep<T>::get_tv(int j, Matrix &M) {
676  cplx *x = tvptr(j);
677  herm_matrix_READ_ELEMENT
678  }
679 
702 template <typename T>
703 template <class Matrix>
704  void herm_matrix_timestep<T>::get_vt(int i, Matrix &M, int sig) {
705  cplx *x = tvptr(ntau_ - i);
706  herm_matrix_READ_ELEMENT_MINUS_CONJ if (sig == -1) M = -M;
707  }
708 
727 template <typename T>
728 template <class Matrix>
729  void herm_matrix_timestep<T>::get_mat(int i, Matrix &M) {
730  cplx *x = matptr(i);
731  herm_matrix_READ_ELEMENT
732  }
733 
755 template <typename T>
756 template <class Matrix>
757  void herm_matrix_timestep<T>::get_matminus(int i, Matrix &M, int sig) {
758  cplx *x = matptr(ntau_ - i);
759  herm_matrix_READ_ELEMENT if (sig == -1) M = -M;
760  }
761 
762 
784 template <typename T>
785 template <class Matrix>
786  void herm_matrix_timestep<T>::get_matminus(int i, Matrix &M) {
787  cplx *x = matptr(ntau_ - i);
788  herm_matrix_READ_ELEMENT if (sig_ == -1) M = -M;
789  }
790 
809 template <typename T>
810 template <class Matrix>
811  void herm_matrix_timestep<T>::get_gtr_tstp_t(int i, Matrix &M) {
812  Matrix M1;
813  get_ret_tstp_t(i, M);
814  get_les_tstp_t(i, M1);
815  M += M1;
816  }
817 
836 template <typename T>
837 template <class Matrix>
838  void herm_matrix_timestep<T>::get_gtr_t_tstp(int i, Matrix &M) {
839  Matrix M1;
840  get_ret_t_tstp(i, M);
841  get_les_t_tstp(i, M1);
842  M += M1;
843  }
844 
846 // Note: same syntax like for herm_matrix (with dummy argument tstp).
847 // The following routines are not very "efficent" but sometimes simple to implement
849 
877 template <typename T>
878  inline void herm_matrix_timestep<T>::get_ret(int i, int j, cplx &x) {
879  assert(i == tstp_ || j == tstp_);
880 
881  if (i == tstp_) {
882  x = *retptr(j);
883  } else {
884  x = *retptr(i);
885  x = -std::conj(x);
886  }
887  }
888 
915 template <typename T>
916  inline void herm_matrix_timestep<T>::get_les(int i, int j, cplx &x) {
917  assert(i == tstp_ || j == tstp_);
918 
919  if (j == tstp_) {
920  x = *lesptr(i);
921  } else {
922  x = *lesptr(j);
923  x = -std::conj(x);
924  }
925  }
926 
953 template <typename T>
954  inline void herm_matrix_timestep<T>::get_tv(int i, int j, cplx &x) {
955  assert(i == tstp_);
956  x = *tvptr(j);
957  }
958 
985 template <typename T>
986  inline void herm_matrix_timestep<T>::get_vt(int i, int j, cplx &x) {
987  assert(j == tstp_);
988  x = *tvptr(ntau_ - i);
989  if (sig_ == -1)
990  x = std::conj(x);
991  else
992  x = -std::conj(x);
993  }
994 
1014 template <typename T>
1015  inline void herm_matrix_timestep<T>::get_mat(int i, cplx &x) {
1016  assert(tstp_ == -1);
1017  x = *matptr(i);
1018  }
1019 
1039 template <typename T>
1040  inline void herm_matrix_timestep<T>::get_matminus(int i, cplx &x) {
1041  assert(tstp_ == -1);
1042  x = *matptr(ntau_ - i);
1043  if (sig_ == -1)
1044  x = -x;
1045  }
1046 
1047 
1067 template <typename T>
1068  std::complex<T> herm_matrix_timestep<T>::density_matrix(int tstp) {
1069  assert(tstp_ == tstp);
1070  cplx x1;
1071  if (tstp_ == -1) {
1072  get_mat(ntau_, x1);
1073  return -x1;
1074  } else {
1075  get_les(tstp_, tstp_, x1);
1076  return std::complex<T>(0.0, sig_) * x1;
1077  }
1078  }
1079 
1101 template <typename T>
1102  inline void herm_matrix_timestep<T>::density_matrix(int tstp, cplx &rho) {
1103  assert(tstp_ == tstp);
1104  cplx x1;
1105  if (tstp_ == -1) {
1106  get_mat(ntau_, x1);
1107  rho = -x1;
1108  } else {
1109  get_les(tstp_, tstp_, x1);
1110  rho = std::complex<T>(0.0, sig_) * x1;
1111  }
1112  }
1113 
1114 
1133 template<typename T> template<class Matrix>
1135 
1136  if(tstp_==-1){
1137  get_mat(ntau_,M);
1138  M *= (-1.0);
1139  }else{
1140  get_les_tstp_t(tstp_,M);
1141  M *= std::complex<T>(0.0,1.0*sig_);
1142  }
1143  }
1144 
1145 
1166 template<typename T> template<class Matrix>
1167  void herm_matrix_timestep<T>::density_matrix(int tstp, Matrix &M){
1168  assert(tstp == tstp_);
1169  if(tstp_==-1){
1170  get_mat(ntau_,M);
1171  M *= (-1.0);
1172  }else{
1173  get_les_tstp_t(tstp_,M);
1174  M *= std::complex<T>(0.0,1.0*sig_);
1175  }
1176  }
1177 
1178 
1179 /* #######################################################################################
1180 #
1181 # WRITING ELEMENTS FROM ANY MATRIX TYPE
1182 # OR FROM COMPLEX NUMBERS (then only the (0,0) element is addressed for dim>0)
1183 #
1184 ########################################################################################*/
1185 #define herm_matrix_SET_ELEMENT_MATRIX \
1186  { \
1187  int r, s; \
1188  for (r = 0; r < size1_; r++) \
1189  for (s = 0; s < size2_; s++) \
1190  x[r * size2_ + s] = M(r, s); \
1191  }
1192 
1212 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_ret(int j,Matrix &M){
1213  cplx *x=retptr(j);
1214  herm_matrix_SET_ELEMENT_MATRIX
1215  }
1216 
1217 
1238 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_ret(int i, int j,Matrix &M){
1239  assert(i == tstp_);
1240  assert(j <= i);
1241  cplx *x=retptr(j);
1242  herm_matrix_SET_ELEMENT_MATRIX
1243  }
1244 
1245 
1266 template<typename T> inline void herm_matrix_timestep<T>::set_ret(int i, int j, cplx &x){
1267  assert(i == tstp_ && j <= i);
1268  *retptr(j) = x;
1269  }
1270 
1271 
1291 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_les(int j,Matrix &M){
1292  cplx *x=lesptr(j);
1293  herm_matrix_SET_ELEMENT_MATRIX
1294  }
1295 
1316 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_les(int i, int j, Matrix &M){
1317  assert(j == tstp_);
1318  assert(i <= j);
1319  cplx *x=lesptr(i);
1320  herm_matrix_SET_ELEMENT_MATRIX
1321  }
1322 
1343 template<typename T> inline void herm_matrix_timestep<T>::set_les(int i, int j, cplx &x){
1344  assert(j == tstp_ && i <= j);
1345  *lesptr(i) = x;
1346  }
1347 
1348 
1368 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_tv(int j,Matrix &M){
1369  cplx *x=tvptr(j);
1370  herm_matrix_SET_ELEMENT_MATRIX
1371 }
1372 
1373 
1394 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_tv(int i, int j, Matrix &M){
1395  assert(i == tstp_);
1396  assert(j <= ntau_);
1397  cplx *x=tvptr(j);
1398  herm_matrix_SET_ELEMENT_MATRIX
1399 }
1400 
1401 
1422 template<typename T> inline void herm_matrix_timestep<T>::set_tv(int i, int j, cplx &x){
1423  assert(i == tstp_);
1424  assert(j <= ntau_);
1425  *tvptr(j) = x;
1426 }
1427 
1428 
1447 template<typename T> template <class Matrix> void herm_matrix_timestep<T>::set_mat(int i,Matrix &M){
1448  assert(i <= ntau_ );
1449  cplx *x=matptr(i);
1450  herm_matrix_SET_ELEMENT_MATRIX
1451 }
1452 
1471 template<typename T> inline void herm_matrix_timestep<T>::set_mat(int i,cplx &x){
1472  assert(i <= ntau_ );
1473  *matptr(i) = x;
1474 }
1475 
1476 
1477 
1478 
1479 
1480 
1481 
1482 
1483 
1484 
1485 
1486 
1487 
1488 
1489 
1490 // G(t,t') ==> F(t)G(t,t') ... ft+t*element_size_ points to F(t)
1491 
1512 template <typename T>
1513 void herm_matrix_timestep<T>::left_multiply(std::complex<T> *f0,
1514  std::complex<T> *ft, T weight) {
1515  int m;
1516  cplx *x0, *xtemp, *ftemp;
1517  xtemp = new cplx[element_size_];
1518  if (tstp_ == -1) {
1519  x0 = data_;
1520  for (m = 0; m <= ntau_; m++) {
1521  element_mult<T, LARGESIZE>(size1_, xtemp, f0,
1522  x0 + m * element_size_);
1523  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1524  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1525  }
1526  } else {
1527  ftemp = ft + tstp_ * element_size_;
1528  x0 = data_;
1529  for (m = 0; m <= tstp_; m++) {
1530  element_mult<T, LARGESIZE>(size1_, xtemp, ftemp,
1531  x0 + m * element_size_);
1532  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1533  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1534  }
1535  x0 = data_ + (tstp_ + 1) * element_size_;
1536  for (m = 0; m <= ntau_; m++) {
1537  element_mult<T, LARGESIZE>(size1_, xtemp, ftemp,
1538  x0 + m * element_size_);
1539  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1540  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1541  }
1542  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
1543  for (m = 0; m <= tstp_; m++) {
1544  element_mult<T, LARGESIZE>(size1_, xtemp, ft + m * element_size_,
1545  x0 + m * element_size_);
1546  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1547  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1548  }
1549  }
1550  delete[] xtemp;
1551 }
1552 
1571 template <typename T>
1573  assert( ft.nt() >= tstp_);
1574 
1575  this->left_multiply(ft.ptr(-1), ft.ptr(0), weight);
1576 }
1577 
1598 template <typename T>
1600  assert(tstp == tstp_);
1601  assert( ft.nt() >= tstp_);
1602 
1603  this->left_multiply(ft.ptr(-1), ft.ptr(0), weight);
1604 }
1605 
1606 
1629 template <typename T>
1631  T weight) {
1632  assert(tstp == tstp_);
1633  int m;
1634  cplx *xtemp, *ftemp, *x0, *f0, *fcc;
1635  xtemp = new cplx[element_size_];
1636  fcc = new cplx[element_size_];
1637  assert(tstp >= -1 && ft.nt() >= tstp &&
1638  ft.size1() == size1_ && ft.size2() == size2_ &&
1639  "tstp >= -1 && ft.nt() >= tstp && ft.size1() == size1_ && ft.size2() == size2_");
1640 
1641  f0 = ft.ptr(-1);
1642  element_conj<T, LARGESIZE>(size1_, fcc, f0);
1643  if (tstp == -1) {
1644  x0 = data_;
1645  for (m = 0; m <= ntau_; m++) {
1646  element_mult<T, LARGESIZE>(size1_, xtemp, fcc,
1647  x0 + m * element_size_);
1648  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1649  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1650  }
1651  } else {
1652  ftemp = ft.ptr(tstp);
1653  element_conj<T, LARGESIZE>(size1_, fcc, ftemp);
1654  x0 = data_;
1655  for (m = 0; m <= tstp; m++) {
1656  element_mult<T, LARGESIZE>(size1_, xtemp, fcc,
1657  x0 + m * element_size_);
1658  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1659  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1660  }
1661  x0 = data_ + (tstp_ + 1) * element_size_;
1662  for (m = 0; m <= ntau_; m++) {
1663  element_mult<T, LARGESIZE>(size1_, xtemp, fcc,
1664  x0 + m * element_size_);
1665  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1666  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1667  }
1668  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
1669  for (m = 0; m <= tstp; m++) {
1670  element_conj<T, LARGESIZE>(size1_, fcc, ft.ptr(m));
1671  element_mult<T, LARGESIZE>(size1_, xtemp, fcc,
1672  x0 + m * element_size_);
1673  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1674  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1675  }
1676  }
1677  delete[] xtemp;
1678  delete[] fcc;
1679 }
1680 
1681 
1682 
1683 
1684 
1685 // G(t,t') ==> G(t,t')F(t')
1686 
1707 template <typename T>
1708 void herm_matrix_timestep<T>::right_multiply(std::complex<T> *f0,
1709  std::complex<T> *ft, T weight) {
1710  int m;
1711  cplx *xtemp, *ftemp, *x0;
1712  xtemp = new cplx[element_size_];
1713  if (tstp_ == -1) {
1714  x0 = data_;
1715  for (m = 0; m <= ntau_; m++) {
1716  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1717  f0);
1718  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1719  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1720  }
1721  } else {
1722  x0 = data_;
1723  for (m = 0; m <= tstp_; m++) {
1724  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1725  ft + m * element_size_);
1726  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1727  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1728  }
1729  x0 = data_ + (tstp_ + 1) * element_size_;
1730  for (m = 0; m <= ntau_; m++) {
1731  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1732  f0);
1733  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1734  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1735  }
1736  ftemp = ft + tstp_ * element_size_;
1737  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
1738  for (m = 0; m <= tstp_; m++) {
1739  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1740  ftemp);
1741  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1742  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1743  }
1744  }
1745  delete[] xtemp;
1746 }
1747 
1766 template <typename T>
1768  assert(ft.nt() >= tstp_);
1769 
1770  this->right_multiply(ft.ptr(-1), ft.ptr(0), weight);
1771 }
1772 
1773 
1794 template <typename T>
1796  assert(tstp == tstp_);
1797  assert(ft.nt() >= tstp_);
1798 
1799  this->right_multiply(ft.ptr(-1), ft.ptr(0), weight);
1800 }
1801 
1802 
1803 
1827 template <typename T>
1829  T weight) {
1830  assert(tstp == tstp_);
1831  int m;
1832  cplx *xtemp, *ftemp, *x0, *f0, *fcc;
1833  xtemp = new cplx[element_size_];
1834  fcc = new cplx[element_size_];
1835  assert(ft.size1() == size1_ && "ft.size1() == size1_");
1836  assert(ft.nt() >= tstp && "ft.nt() >= tstp");
1837  assert(tstp >= -1&& "ft.nt() >= tstp");
1838  f0 = ft.ptr(-1);
1839  element_conj<T, LARGESIZE>(size1_, fcc, f0);
1840  if (tstp == -1) {
1841  x0 = data_;
1842  for (m = 0; m <= ntau_; m++) {
1843  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1844  fcc);
1845  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1846  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1847  }
1848  } else {
1849  x0 = data_;
1850  for (m = 0; m <= tstp; m++) {
1851  element_conj<T, LARGESIZE>(size1_, fcc, ft.ptr(m));
1852  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1853  fcc);
1854  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1855  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1856  }
1857  x0 = data_ + (tstp_ + 1) * element_size_;
1858  element_conj<T, LARGESIZE>(size1_, fcc, f0);
1859  for (m = 0; m <= ntau_; m++) {
1860  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1861  fcc);
1862  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1863  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1864  }
1865  ftemp = ft.ptr(tstp);
1866  element_conj<T, LARGESIZE>(size1_, fcc, ftemp);
1867  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
1868  for (m = 0; m <= tstp; m++) {
1869  element_mult<T, LARGESIZE>(size1_, xtemp, x0 + m * element_size_,
1870  fcc);
1871  element_smul<T, LARGESIZE>(size1_, xtemp, weight);
1872  element_set<T, LARGESIZE>(size1_, x0 + m * element_size_, xtemp);
1873  }
1874  }
1875  delete[] fcc;
1876  delete[] xtemp;
1877 }
1878 
1879 
1899 template <typename T>
1901  int m;
1902  assert(g1.size1_ == size1_ && g1.ntau_ == ntau_ && g1.tstp_ == tstp_);
1903  for (m = 0; m < total_size_; m++)
1904  data_[m] += weight * g1.data_[m];
1905 }
1906 
1907 
1930 template <typename T>
1932  assert(tstp == tstp_);
1933  assert(g1.size1_ == size1_ && g1.ntau_ == ntau_ && g1.tstp_ == tstp_);
1934  for (int m = 0; m < total_size_; m++)
1935  data_[m] += weight * g1.data_[m];
1936 }
1937 
1938 
1939 
1940 #define HERM_MATRIX_INCR_TSTP \
1941 if (alpha == cplx(1.0, 0.0)) { \
1942  for (i = 0; i < len; i++) \
1943  x0[i] += x[i]; \
1944 } else { \
1945  for (i = 0; i < len; i++) \
1946  x0[i] += alpha * x[i]; \
1947 }
1948 
1969 template <typename T>
1971  int i, len;
1972  cplx *x, *x0;
1973  assert(tstp_ <= g.nt() && ntau_ == g.ntau() && size1_ == g.size1());
1974  if (tstp_ == -1) {
1975  len = (ntau_ + 1) * element_size_;
1976  x = g.matptr(0);
1977  x0 = data_;
1978  HERM_MATRIX_INCR_TSTP
1979  } else {
1980  len = (tstp_ + 1) * element_size_;
1981  x = g.retptr(tstp_, 0);
1982  x0 = data_;
1983  HERM_MATRIX_INCR_TSTP
1984  len = (ntau_ + 1) * element_size_;
1985  x = g.tvptr(tstp_, 0);
1986  x0 = data_ + (tstp_ + 1) * element_size_;
1987  HERM_MATRIX_INCR_TSTP
1988  len = (tstp_ + 1) * element_size_;
1989  x = g.lesptr(0, tstp_);
1990  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
1991  HERM_MATRIX_INCR_TSTP
1992  }
1993 }
1994 
2018 template <typename T>
2020  assert(tstp == tstp_);
2021  assert(tstp_ <= g.nt() && ntau_ == g.ntau() && size1_ == g.size1());
2022  int i, len;
2023  cplx *x, *x0;
2024  if (tstp_ == -1) {
2025  len = (ntau_ + 1) * element_size_;
2026  x = g.matptr(0);
2027  x0 = data_;
2028  HERM_MATRIX_INCR_TSTP
2029  } else {
2030  len = (tstp_ + 1) * element_size_;
2031  x = g.retptr(tstp_, 0);
2032  x0 = data_;
2033  HERM_MATRIX_INCR_TSTP
2034  len = (ntau_ + 1) * element_size_;
2035  x = g.tvptr(tstp_, 0);
2036  x0 = data_ + (tstp_ + 1) * element_size_;
2037  HERM_MATRIX_INCR_TSTP
2038  len = (tstp_ + 1) * element_size_;
2039  x = g.lesptr(0, tstp_);
2040  x0 = data_ + (tstp_ + 1 + ntau_ + 1) * element_size_;
2041  HERM_MATRIX_INCR_TSTP
2042  }
2043 }
2044 
2045 
2046 #undef HERM_MATRIX_INCR_TSTP
2047 
2067 template <typename T>
2069  herm_matrix_timestep_view<T> g_view(*this);
2070 
2071  g_view.smul(weight);
2072 
2073 }
2074 
2075 
2098 template <typename T>
2099 void herm_matrix_timestep<T>::smul(int tstp, T weight) {
2100  assert(tstp == tstp_);
2101  herm_matrix_timestep_view<T> g_view(*this);
2102 
2103  g_view.smul(weight);
2104 
2105 }
2106 
2107 
2108 
2135 template <typename T>
2137  herm_matrix_timestep<T> &g, int j1, int j2) {
2139  herm_matrix_timestep_view<T> tmp1(*this);
2140  tmp1.set_matrixelement(i1, i2, tmp, j1, j2);
2141 }
2142 
2143 
2172 template <typename T>
2173 void herm_matrix_timestep<T>::set_matrixelement(int tstp, int i1, int i2,
2174  herm_matrix_timestep<T> &g, int j1, int j2) {
2175  assert(tstp == tstp_);
2177  herm_matrix_timestep_view<T> tmp1(*this);
2178  tmp1.set_matrixelement(i1, i2, tmp, j1, j2);
2179 }
2180 
2181 
2209 template <typename T>
2211  int i1, int i2, herm_matrix_timestep_view<T> &g, int j1, int j2) {
2212  herm_matrix_timestep_view<T> tmp1(*this);
2213  tmp1.set_matrixelement(i1, i2, g, j1, j2);
2214 }
2215 
2216 
2246 template <typename T>
2248  int i1, int i2, herm_matrix_timestep_view<T> &g, int j1, int j2) {
2249  herm_matrix_timestep_view<T> tmp1(*this);
2250  tmp1.set_matrixelement(i1, i2, g, j1, j2);
2251 }
2252 
2253 
2281 template <typename T>
2283  herm_matrix<T> &g, int j1,
2284  int j2) {
2285  herm_matrix_timestep_view<T> tmp(tstp_, g);
2286  herm_matrix_timestep_view<T> tmp1(*this);
2287  tmp1.set_matrixelement(i1, i2, tmp, j1, j2);
2288 }
2289 
2290 
2320 template <typename T>
2321 void herm_matrix_timestep<T>::set_matrixelement(int tstp, int i1, int i2,
2322  herm_matrix<T> &g, int j1,
2323  int j2) {
2324  assert(tstp == tstp_);
2325  herm_matrix_timestep_view<T> tmp(tstp_, g);
2326  herm_matrix_timestep_view<T> tmp1(*this);
2327  tmp1.set_matrixelement(i1, i2, tmp, j1, j2);
2328 }
2329 
2359 template<typename T> void herm_matrix_timestep<T>::set_submatrix(int tstp, std::vector<int> &i1,
2360  std::vector<int> &i2,herm_matrix<T> &g,std::vector<int> &j1,std::vector<int> &j2){
2361  assert(tstp == tstp_);
2362  assert(tstp <= g.nt());
2363  assert(i1.size()==i2.size() && i1.size()==j1.size() && j1.size()==j2.size());
2364  assert(size1_*size2_==i1.size());
2365  for(int k1=0;k1<i1.size();k1++){
2366  this->set_matrixelement(tstp,i1[k1],i2[k1],g,j1[k1],j2[k1]);
2367  }
2368 }
2369 
2370 
2400 template<typename T> void herm_matrix_timestep<T>::set_submatrix(int tstp, std::vector<int> &i1,
2401  std::vector<int> &i2,herm_matrix_timestep<T> &g,std::vector<int> &j1,std::vector<int> &j2){
2402  assert(tstp == tstp_);
2403  assert(tstp == g.tstp());
2404  assert(i1.size()==i2.size() && i1.size()==j1.size() && j1.size()==j2.size());
2405  assert(size1_*size2_==i1.size());
2406  for(int k1=0;k1<i1.size();k1++){
2407  this->set_matrixelement(tstp,i1[k1],i2[k1],g,j1[k1],j2[k1]);
2408  }
2409 }
2410 
2411 /* #########################################################################
2412 #
2413 # MPI UTILS
2414 #
2415 ############################################################################ */
2416 
2417 
2418 #if CNTR_USE_MPI == 1
2419 
2420 
2437 template <typename T>
2439  int taskid;
2440  int len = 2 * (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2441  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2442  if (sizeof(T) == sizeof(double)) {
2443  if (taskid == root) {
2444  MPI_Reduce(MPI_IN_PLACE, (double *)this->data_, len, MPI_DOUBLE_PRECISION, MPI_SUM, root, MPI_COMM_WORLD);
2445  } else {
2446  MPI_Reduce((double *)this->data_,
2447  (double *)this->data_, len, MPI_DOUBLE_PRECISION, MPI_SUM, root, MPI_COMM_WORLD);
2448  }
2449  } else {
2450  std::cerr << "herm_matrix_timestep<T>::MPI_Reduce only for double "
2451  << std::endl;
2452  exit(0);
2453  }
2454 }
2455 
2456 
2475 template <typename T>
2477  assert(tstp == tstp_);
2478  int taskid;
2479  int len = 2 * (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2480  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2481  if (sizeof(T) == sizeof(double)) {
2482  if (taskid == root) {
2483  MPI_Reduce(MPI_IN_PLACE, (double *)this->data_, len,
2484  MPI_DOUBLE_PRECISION, MPI_SUM, root, MPI_COMM_WORLD);
2485  } else {
2486  MPI_Reduce((double *)this->data_,
2487  (double *)this->data_, len, MPI_DOUBLE_PRECISION,
2488  MPI_SUM, root, MPI_COMM_WORLD);
2489  }
2490  } else {
2491  std::cerr << "herm_matrix_timestep<T>::MPI_Reduce only for double "
2492  << std::endl;
2493  exit(0);
2494  }
2495 }
2496 
2497 
2498 
2525 template <typename T>
2526 void herm_matrix_timestep<T>::Bcast_timestep(int tstp, int ntau, int size1,
2527  int root) {
2528  int numtasks,taskid;
2529  MPI_Comm_size(MPI_COMM_WORLD, &numtasks);
2530  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2531  if (taskid != root)
2532  resize(tstp, ntau, size1);
2533  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2534 // test effective on root:
2535  assert(tstp == tstp_);
2536  assert(ntau == ntau_);
2537  assert(size1 == size1_);
2538  if (sizeof(T) == sizeof(double))
2539  MPI_Bcast(data_, len, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
2540  else
2541  MPI_Bcast(data_, len, MPI_COMPLEX, root, MPI_COMM_WORLD);
2542 }
2543 
2566 template <typename T>
2568  int numtasks,taskid;
2569  MPI_Comm_size(MPI_COMM_WORLD, &numtasks);
2570  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2571  if (taskid != root)
2572  resize(tstp, ntau_, size1_);
2573  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2574 // test effective on root:
2575  assert(tstp == tstp_);
2576  if (sizeof(T) == sizeof(double))
2577  MPI_Bcast(data_, len, MPI_DOUBLE_COMPLEX, root, MPI_COMM_WORLD);
2578  else
2579  MPI_Bcast(data_, len, MPI_COMPLEX, root, MPI_COMM_WORLD);
2580 }
2581 
2605 template <typename T>
2606 void herm_matrix_timestep<T>::Send_timestep(int tstp, int ntau, int size1,
2607  int dest, int tag) {
2608  int taskid;
2609  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2610  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2611  if (!(taskid == dest)) {
2612  assert(tstp == tstp_);
2613  assert(ntau == ntau_);
2614  assert(size1 == size1_);
2615  if (sizeof(T) == sizeof(double))
2616  MPI_Send(data_, len, MPI_DOUBLE_COMPLEX, dest, tag, MPI_COMM_WORLD);
2617  else
2618  MPI_Send(data_, len, MPI_COMPLEX, dest, tag, MPI_COMM_WORLD);
2619  } else {
2620  }
2621 }
2622 
2642 template <typename T>
2643 void herm_matrix_timestep<T>::Send_timestep(int tstp, int dest, int tag) {
2644  int taskid;
2645  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2646  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2647  if (!(taskid == dest)) {
2648  assert(tstp == tstp_);
2649  if (sizeof(T) == sizeof(double))
2650  MPI_Send(data_, len, MPI_DOUBLE_COMPLEX, dest, tag, MPI_COMM_WORLD);
2651  else
2652  MPI_Send(data_, len, MPI_COMPLEX, dest, tag, MPI_COMM_WORLD);
2653  } else {
2654  }
2655 }
2656 
2680 template <typename T>
2681 void herm_matrix_timestep<T>::Recv_timestep(int tstp, int ntau, int size1,
2682  int root, int tag) {
2683  int taskid;
2684  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2685  if (!(taskid == root)) {
2686  resize(tstp, ntau, size1);
2687  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2688  if (sizeof(T) == sizeof(double))
2689  MPI_Recv(data_, len, MPI_DOUBLE_COMPLEX, root, tag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
2690  else
2691  MPI_Recv(data_, len, MPI_COMPLEX, root, tag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
2692  }
2693 }
2694 
2695 
2715 template <typename T>
2716 void herm_matrix_timestep<T>::Recv_timestep(int tstp, int root, int tag) {
2717  int taskid;
2718  MPI_Comm_rank(MPI_COMM_WORLD, &taskid);
2719  if (!(taskid == root)) {
2720  resize(tstp, ntau_, size1_);
2721  int len = (2 * (tstp_ + 1) + ntau_ + 1) * element_size_;
2722  if (sizeof(T) == sizeof(double))
2723  MPI_Recv(data_, len, MPI_DOUBLE_COMPLEX, root, tag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
2724  else
2725  MPI_Recv(data_, len, MPI_COMPLEX, root, tag, MPI_COMM_WORLD, MPI_STATUS_IGNORE);
2726  }
2727 }
2728 
2729 
2730 #endif
2731 
2732 #undef herm_matrix_READ_ELEMENT
2733 #undef herm_matrix_READ_ELEMENT_MINUS_CONJ
2734 
2735 /* ##########################################################################
2736 #
2737 # HDF5 I/O (using implementation in herm_matrix_timestep_view)
2738 #
2739 #############################################################################*/
2740 
2741 #if CNTR_USE_HDF5 == 1
2742 
2761 template <typename T>
2763  herm_matrix_timestep_view<T> tmp(*this);
2764  tmp.write_to_hdf5(group_id);
2765 }
2766 
2787 template <typename T>
2789  const char *groupname) {
2790  herm_matrix_timestep_view<T> tmp(*this);
2791  tmp.write_to_hdf5(group_id, groupname);
2792 }
2793 
2814 template <typename T>
2815 void herm_matrix_timestep<T>::write_to_hdf5(const char *filename,
2816  const char *groupname) {
2817  herm_matrix_timestep_view<T> tmp(*this);
2818  tmp.write_to_hdf5(filename, groupname);
2819 }
2820 
2839 template <typename T>
2841 // -- implementation is different from herm_matrix_timestep_view, because G can be resized:
2842  int tstp = read_primitive_type<int>(group_id, "tstp");
2843  int ntau = read_primitive_type<int>(group_id, "ntau");
2844  int sig = read_primitive_type<int>(group_id, "sig");
2845  int size1 = read_primitive_type<int>(group_id, "size1");
2846 // RESIZE G
2847  this->resize(tstp, ntau, size1);
2848  sig_ = sig;
2849  herm_matrix_timestep_view<T> tmp(*this);
2850  tmp.read_from_hdf5(group_id);
2851 }
2852 
2872 template <typename T>
2874  const char *groupname) {
2875  hid_t sub_group_id = open_group(group_id, groupname);
2876  this->read_from_hdf5(sub_group_id);
2877  close_group(sub_group_id);
2878 }
2879 
2899 template <typename T>
2901  const char *groupname) {
2902  hid_t file_id = read_hdf5_file(filename);
2903  this->read_from_hdf5(file_id, groupname);
2904  close_hdf5_file(file_id);
2905 }
2906 #endif
2907 
2908 } // namespace cntr
2909 
2910 #endif // CNTR_HERM_MATRIX_TIMESTEP_IMPL_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 get_mat(const int m, std::complex< T > &G_mat, herm_matrix< T > &G, herm_matrix< T > &Gcc)
void density_matrix(int tstp, Matrix &M)
Returns the density matrix at given time step.
int size1(void) const
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 smul(T alpha)
Multiply herm_matrix_timestep_view with scalar weight.
void set_timestep_zero(int tstp)
Sets all components at time step tstp to zero.
int size2(void) const
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...
std::complex< double > cplx
Definition: fourier.cpp:11
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
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 read_from_hdf5(hid_t group_id, const char *groupname)
Read herm_matrix_timestep_view from hdf5 group given by group ID group_id and group name groupname...
int nt(void) const
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 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 write_to_hdf5(hid_t group_id)
Write herm_matrix_timestep_view to hdf5 group given by group_id.
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 get_les(const int i, const int j, std::complex< T > &G_les, herm_matrix< T > &G, herm_matrix< T > &Gcc)
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