NESSi  v1.0.2
The NonEquilibrium Systems Simulation Library
cntr_function_decl.hpp
Go to the documentation of this file.
1 #ifndef CNTR_FUNCTION_DECL_H
2 #define CNTR_FUNCTION_DECL_H
3 
5 
6 namespace cntr {
7 
8 template <typename T>
26 class function {
27  public:
28  // this is like a vector with pointer access
29  typedef std::complex<T> cplx;
30  function();
31  ~function();
32  function(int nt, int size1 = 1);
33  function(int nt,int size1,int size2);
34  function(const function &f);
35  function &operator=(const function &f);
36 #if __cplusplus >= 201103L
37  function(function &&f) noexcept;
38  function &operator=(function &&f) noexcept;
39 #endif
40  int element_size(void) const { return element_size_; }
42  int size1(void) const { return size1_; }
43  int size2(void) const { return size2_; }
44  int nt(void) const { return nt_; }
45  inline cplx *ptr(int t) { return data_ + (t + 1) * element_size_; }
46  inline const cplx *ptr(int t) const {
47  return data_ + (t + 1) * element_size_;
48  }
49  void resize(int nt, int size1);
50  void set_zero(void);
51  void set_constant(cplx *f0); // f0 must be size*size
52  template<class EigenMatrix>
53  void set_constant(EigenMatrix &M);
54  template <class EigenMatrix>
55  void set_value(int tstp, EigenMatrix &M);
56  template <class EigenMatrix>
57  void set_value(int tstp, cplx x);
58  template <class EigenMatrix>
59  void get_value(int tstp, EigenMatrix &M) const;
60  void smul(T weight);
61  template<class EigenMatrix>
62  void set_matrixelement(int tstp,int i1,int i2,EigenMatrix &M,int j1,int j2);
63  void set_matrixelement(int i1,int i2,function &g,int j1,int j2);
64  void left_multiply(cplx *f, T weight = 1.0);
65  void right_multiply(cplx *f, T weight = 1.0);
66  void left_multiply(function<T> &f, T weight = 1.0);
67  void right_multiply(function<T> &f, T weight = 1.0);
68  void get_matrixelement(int i1, int i2, function<T> &g);
69  void incr(function<T> &g, T weight = 1.0);
70 
71  cplx &operator[](int t) { return *ptr(t); } // useful only for size=1
72  const cplx &operator[](int t) const {
73  return *ptr(t);
74  } // useful only for size=1
75  void print_to_file(const char *file, int precision = 16) const;
76  void read_from_file(const char *file);
77  // READ UP TO NT1, NO RESIZE
78  void read_from_file(int nt1, const char *file);
79 #if CNTR_USE_HDF5 == 1
80  // HDF5 I/O
81  void write_to_hdf5(hid_t group_id) const;
82  void write_to_hdf5(hid_t group_id, const char *groupname) const;
83  void write_to_hdf5(const char *filename, const char *groupname) const;
84  void read_from_hdf5(hid_t group_id);
85  void read_from_hdf5(hid_t group_id, const char *groupname);
86  void read_from_hdf5(const char *filename, const char *groupname);
87  void read_from_hdf5(int nt1, hid_t group_id);
88  void read_from_hdf5(int nt1, hid_t group_id, const char *groupname);
89  void read_from_hdf5(int nt1, const char *filename, const char *groupname);
90 #endif
91 #if CNTR_USE_MPI==1
92  void Bcast_timestep(int tstp,int root);
93 #endif
94 
97  int nt_;
99  int size1_;
101  int size2_;
103 
104  int element_size_;
107 };
108 
109 } // namespace cntr
110 
111 #endif // CNTR_FUNCTION_DECL_H
int size2_
Number of the rows in the Matrix form.
int size1(void) const
void Bcast_timestep(int tstp, int root)
Broadcast a time step of the function using MPI
const cplx * ptr(int t) const
void resize(int nt, int size1)
Resize the time length and/or matrix dimension of a square-matrix function
int size2(void) const
function & operator=(const function &f)
Overloaded operator &#39;=&#39;. It copies data from an existing function object
void set_constant(cplx *f0)
Set all data to a constant(scalar or matrix) for the function class
int total_size_
Size of the data stored for the function on the real-time axis including ; * size1 * size2 ...
void set_zero(void)
Set all data to zero for the function class
cplx * data_
Pointer to the function in the Matrix form on the real-time axis ( ) ; &#39;data_+ element_size&#39; corresp...
int size1_
Number of the colums in the Matrix form.
Class function for objects with time on real axis.
void get_value(int tstp, EigenMatrix &M) const
Get matrix value of this function object at a specific time point
void left_multiply(cplx *f, T weight=1.0)
Multiply function values with a matrix-valued function from the left
void read_from_hdf5(hid_t group_id)
Read the function data from a hdf5 file
cplx & operator[](int t)
void incr(function< T > &g, T weight=1.0)
Increase the function object by another function.
void right_multiply(cplx *f, T weight=1.0)
Multiply the function with another matrix-valued function from the right
void set_matrixelement(int tstp, int i1, int i2, EigenMatrix &M, int j1, int j2)
Set a matrix element from an Eigen matrix
int nt(void) const
void get_matrixelement(int i1, int i2, function< T > &g)
Get a matrix element of the function object
std::complex< T > cplx
int nt_
Maximum number of the time steps.
void read_from_file(const char *file)
Read the function data from a txt file
void write_to_hdf5(hid_t group_id) const
Write the function data to a hdf5 file
void print_to_file(const char *file, int precision=16) const
Store the function data to a txt file
void set_value(int tstp, EigenMatrix &M)
Set matrix value at a specific time point
void smul(T weight)
Multiply the function with a scalar
const cplx & operator[](int t) const