LCPQ/dft_tools
3
0
Fork 0
mirror of https://github.com/triqs/dft_tools synced 2024-12-26 06:14:14 +01:00
Code Issues Releases Wiki Activity

Fixed bug out of domain for mesh with positive and negative matsubara frequencies

Browse Source
This commit is contained in:
tayral 2014-01-06 15:43:16 +00:00 committed by Olivier Parcollet
parent 53b05ad899
commit 1b3a66af24
1 changed files with 31 additions and 11 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

24
triqs/gfs/imfreq.hpp
View File

@ -72,33 +72,53 @@ namespace gfs {
// dispatch for 2x2 cases : matrix/scalar and tail/no_tail ( true means no_tail) // dispatch for 2x2 cases : matrix/scalar and tail/no_tail ( true means no_tail)
template <typename G> template <typename G>
std::complex<double> _call_impl(G const *g, matsubara_freq const &f, scalar_valued, std::false_type) const { std::complex<double> _call_impl(G const *g, matsubara_freq const &f, scalar_valued, std::false_type) const {
if (g->mesh().index_start()==0){//only positive Matsubara frequencies
if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n]; if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n];
if ((f.n < 0) && (-f.n < g->mesh().size())) return conj((*g)[-f.n]); if ((f.n < 0) && (-f.n < g->mesh().size())) return conj((*g)[-f.n]);
}
else{
if ((f.n >= g->mesh().index_start()) && (f.n < g->mesh().size()+g->mesh().index_start())) return (*g)[f.n];
}
return g->singularity().evaluate(f)(0, 0); return g->singularity().evaluate(f)(0, 0);
} }
template <typename G> template <typename G>
std::complex<double> _call_impl(G const *g, matsubara_freq const &f, scalar_valued, std::true_type) const { std::complex<double> _call_impl(G const *g, matsubara_freq const &f, scalar_valued, std::true_type) const {
if (g->mesh().index_start()==0){//only positive Matsubara frequencies
if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n]; if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n];
if ((f.n < 0) && (-f.n < g->mesh().size())) return conj((*g)[-f.n]); if ((f.n < 0) && (-f.n < g->mesh().size())) return conj((*g)[-f.n]);
}
else{
if ((f.n >= g->mesh().index_start()) && (f.n < g->mesh().size()+g->mesh().index_start())) return (*g)[f.n];
}
return 0; return 0;
} }
template <typename G> template <typename G>
arrays::matrix_const_view<std::complex<double>> _call_impl(G const *g, matsubara_freq const &f, matrix_valued, arrays::matrix_const_view<std::complex<double>> _call_impl(G const *g, matsubara_freq const &f, matrix_valued,
std::false_type) const { std::false_type) const {
if (g->mesh().index_start()==0){//only positive Matsubara frequencies
if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n](); if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n]();
if ((f.n < 0) && (-f.n < g->mesh().size())) if ((f.n < 0) && (-f.n < g->mesh().size()))
return arrays::matrix<std::complex<double>>{conj((*g)[-f.n]())}; return arrays::matrix<std::complex<double>>{conj((*g)[-f.n]())};
else }
else{
if ((f.n >= g->mesh().index_start()) && (f.n < g->mesh().size()+g->mesh().index_start())) return (*g)[f.n];
}
return g->singularity().evaluate(f); return g->singularity().evaluate(f);
} }
template <typename G> template <typename G>
arrays::matrix_const_view<std::complex<double>> _call_impl(G const *g, matsubara_freq const &f, matrix_valued, arrays::matrix_const_view<std::complex<double>> _call_impl(G const *g, matsubara_freq const &f, matrix_valued,
std::true_type) const { std::true_type) const {
if (g->mesh().index_start()==0){//only positive Matsubara frequencies
if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n](); if ((f.n >= 0) && (f.n < g->mesh().size())) return (*g)[f.n]();
if ((f.n < 0) && (-f.n < g->mesh().size())) return arrays::matrix<std::complex<double>>{conj((*g)[-f.n]())}; if ((f.n < 0) && (-f.n < g->mesh().size()))
return arrays::matrix<std::complex<double>>{conj((*g)[-f.n]())};
}
else{
if ((f.n >= g->mesh().index_start()) && (f.n < g->mesh().size()+g->mesh().index_start())) return (*g)[f.n];
}
auto r = arrays::matrix<std::complex<double>>{get_target_shape(*g)}; auto r = arrays::matrix<std::complex<double>>{get_target_shape(*g)};
r() = 0; r() = 0;
return r; return r;