3
0
mirror of https://github.com/triqs/dft_tools synced 2024-12-25 05:43:40 +01:00

gfs: Fix bug in imfreq with <0 freq.

- there was a confusion in gf imfreq, in the new case
where freq can be <0 (non real gf, or for product gf).
- index: is the matsubara n, as in the struct matsubara_freq
  index can be >0 or <0
- linear_index : is the shift from the 0. It is always >0.
  Fixed function to compute it.

- Also changed the construction of mesh_point in the generic iterator.
  Before, was constructed with a mesh point of index 0
  Now, added a new constructor on mesh_point_t, just taking the mesh
  which construct the *first* mesh_point.
  Fixed linear, discrete, product accordingly.
  Added to the documentation of the concepts of gf.
This commit is contained in:
Olivier Parcollet 2014-01-05 20:10:33 +01:00
parent b8e29192da
commit 1ecec0b933
6 changed files with 34 additions and 29 deletions

View File

@ -146,30 +146,32 @@ MeshPoint
* **Definition** :
+------------------------------------------------+-----------------------------------------------------------------------------+
| Elements | Comment |
+================================================+=============================================================================+
| mesh_t | Type of the mesh |
+------------------------------------------------+-----------------------------------------------------------------------------+
| mesh_t const * m | A pointer to the mesh to which the point belongs. |
+------------------------------------------------+-----------------------------------------------------------------------------+
| mesh_t::index_t index | The index of the point |
+------------------------------------------------+-----------------------------------------------------------------------------+
| mesh_point_t( mesh_t const &, index_t const &) | Constructor |
+------------------------------------------------+-----------------------------------------------------------------------------+
| mesh_t::index_t [const &|] index() const | The index corresponding to the point |
+------------------------------------------------+-----------------------------------------------------------------------------+
| size_t linear_index() const | The linear index of the point (same as m->index_to_linear(index()) |
+------------------------------------------------+-----------------------------------------------------------------------------+
| void advance() | Advance to the next point on the mesh (used by iterators). |
+------------------------------------------------+-----------------------------------------------------------------------------+
| void at_end() | Is the point at the end of the grid |
+------------------------------------------------+-----------------------------------------------------------------------------+
| void reset() | Reset the mesh point to the first point |
+------------------------------------------------+-----------------------------------------------------------------------------+
| cast_t | == mesh_t::domain_t::point_t |
| operator cast_t() const | *implicit* cast to the corresponding domain point |
+------------------------------------------------+-----------------------------------------------------------------------------+
+------------------------------------------------+--------------------------------------------------------------------+
| Elements | Comment |
+================================================+====================================================================+
| mesh_t | Type of the mesh |
+------------------------------------------------+--------------------------------------------------------------------+
| mesh_t const * m | A pointer to the mesh to which the point belongs. |
+------------------------------------------------+--------------------------------------------------------------------+
| mesh_t::index_t index | The index of the point |
+------------------------------------------------+--------------------------------------------------------------------+
| mesh_point_t( mesh_t const &, index_t const &) | Constructor : a mesh point at the given index |
+------------------------------------------------+--------------------------------------------------------------------+
| mesh_point_t( mesh_t const &) | Constructor : the first mesh point |
+------------------------------------------------+--------------------------------------------------------------------+
| mesh_t::index_t [const &,] index() const | The index corresponding to the point |
+------------------------------------------------+--------------------------------------------------------------------+
| size_t linear_index() const | The linear index of the point (same as m->index_to_linear(index()) |
+------------------------------------------------+--------------------------------------------------------------------+
| void advance() | Advance to the next point on the mesh (used by iterators). |
+------------------------------------------------+--------------------------------------------------------------------+
| void at_end() | Is the point at the end of the grid |
+------------------------------------------------+--------------------------------------------------------------------+
| void reset() | Reset the mesh point to the first point |
+------------------------------------------------+--------------------------------------------------------------------+
| cast_t operator cast_t() const | == mesh_t::domain_t::point_t *implicit* cast to the corresponding |
| | domain point |
+------------------------------------------------+--------------------------------------------------------------------+
For one dimensional mesh, we also require that the MeshPoint implement the basic arithmetic operations
using the cast.

View File

@ -47,6 +47,7 @@ namespace gfs {
public:
mesh_point_t() = default;
mesh_point_t(discrete_mesh const &mesh, index_t const &index_) : m(&mesh), _index(index_) {}
mesh_point_t(discrete_mesh const &mesh) : mesh_point_t(mesh, 0){}
void advance() { ++_index; }
using cast_t = long;
operator cast_t() const { return m->index_to_point(_index); }

View File

@ -84,6 +84,7 @@ namespace gfs {
public:
mesh_point_t() : m(nullptr) {}
mesh_point_t(linear_mesh const &mesh, index_t const &index_) : m(&mesh), _index(index_) {}
mesh_point_t(linear_mesh const &mesh) : mesh_point_t(mesh,0) {}
void advance() { ++_index; }
using cast_t = domain_pt_t;
operator cast_t() const { return m->index_to_point(_index); }

View File

@ -67,7 +67,7 @@ namespace gfs {
domain_pt_t index_to_point(index_t ind) const { return 1_j * M_PI * (2 * ind + (_dom.statistic == Fermion)) / _dom.beta; }
/// Flatten the index in the positive linear index for memory storage (almost trivial here).
long index_to_linear(index_t ind) const { return ind + index_start(); }
long index_to_linear(index_t ind) const { return ind - index_start(); }
/**
* The mesh point
@ -80,8 +80,9 @@ namespace gfs {
: matsubara_freq(index_, mesh.domain().beta, mesh.domain().statistic),
index_start(mesh.index_start()),
index_stop(mesh.index_start() + mesh.size() - 1) {}
mesh_point_t(matsubara_freq_mesh const &mesh) : mesh_point_t(mesh, mesh.index_start()) {}
void advance() { ++n; }
long linear_index() const { return n; }
long linear_index() const { return n - index_start; }
long index() const { return n; }
bool at_end() const { return (n == index_stop); }
void reset() { n = index_start; }

View File

@ -40,7 +40,7 @@ namespace gfs {
//bool equal(mesh_pt_generator const & other) const { return ((mesh == other.mesh) && (other.u==u) );}
public:
mesh_pt_generator(): mesh(nullptr), u(0) {}
mesh_pt_generator( MeshType const * m, bool atEnd = false): mesh(m), u(atEnd ? m->size(): 0), pt((*m)[typename MeshType::index_t()]) {}
mesh_pt_generator( MeshType const * m, bool atEnd = false): mesh(m), u(atEnd ? m->size(): 0), pt(*m) {}
void increment() { ++u; pt.advance(); }
bool at_end() const { return (u>=mesh->size());}
typename MeshType::domain_t::point_t to_point() const { return pt;}

View File

@ -128,14 +128,14 @@ namespace gfs {
template <typename M> typename M::mesh_point_t operator()(M const &m, typename M::index_t const &i) const { return m[i]; }
};
struct F1 {
template <typename M> typename M::mesh_point_t operator()(M const &m) const { return m[typename M::index_t()]; }
template <typename M> typename M::mesh_point_t operator()(M const &m) const { return {m}; }
};
public:
mesh_point_t() = default;
mesh_point_t(mesh_product const &m_, index_t index_)
: m(&m_), _c(triqs::tuple::apply_on_zip(F2(), m_.m_tuple, index_)), _atend(false) {}
mesh_point_t(mesh_product const &m_) : m(&m_), _c(triqs::tuple::apply(F1(), m_.m_tuple)), _atend(false) {}
mesh_point_t(mesh_product const &m_) : m(&m_), _c(triqs::tuple::apply_on_tuple(F1(), m_.m_tuple)), _atend(false) {}
m_pt_tuple_t const &components_tuple() const { return _c; }
size_t linear_index() const { return m->mp_to_linear(_c); }
const mesh_product *mesh() const { return m; }