diff --git a/test/triqs/clef_examples/sum3.cpp b/test/triqs/clef_examples/sum3.cpp new file mode 100644 index 00000000..70d21f0b --- /dev/null +++ b/test/triqs/clef_examples/sum3.cpp @@ -0,0 +1,21 @@ +#include +/// ---- test ------------ + +using triqs::arrays::range; + +int main() { + triqs::clef::placeholder<0> a_; + triqs::clef::placeholder<1> i_; + triqs::clef::placeholder<2> j_; + + auto expr = i_ + 2* j_; + + auto xa = sum(a_*i_ , i_=range(0,5)); + auto x = eval(xa, a_=2); + + std::cout << x << std::endl; + + auto y = sum(2*i_ + j_ , i_=range(0,3), j_= range (0,2)); + std::cout << y << std::endl; +} + diff --git a/test/triqs/gfs/g_k_om2.cpp b/test/triqs/gfs/g_k_om2.cpp index 2139b3ff..66d27001 100644 --- a/test/triqs/gfs/g_k_om2.cpp +++ b/test/triqs/gfs/g_k_om2.cpp @@ -11,8 +11,8 @@ using namespace triqs::arrays; using namespace triqs::lattice; template // requires ( is_function_on_mesh()) -auto sum(Function const &f, Mesh const &m) ->decltype(make_matrix(0*f(*(m.begin())))) { - //auto sum(Function const &f, Mesh const &m) { +auto sum_gf(Function const &f, Mesh const &m) ->decltype(make_matrix(0*f(*(m.begin())))) { + //auto sum_gf(Function const &f, Mesh const &m) { //auto res = typename triqs::regular_type_if_exists_else_type::type (f(m.begin())); auto res = make_matrix(0*f(*(m.begin()))); for (auto const &x : m) res = res + f(x); @@ -20,7 +20,7 @@ auto sum(Function const &f, Mesh const &m) ->decltype(make_matrix(0*f(*(m.begin } namespace triqs { namespace clef { - TRIQS_CLEF_MAKE_FNT_LAZY(sum); + TRIQS_CLEF_MAKE_FNT_LAZY(sum_gf); TRIQS_CLEF_MAKE_FNT_LAZY(conj); } } @@ -49,8 +49,8 @@ int main() { auto r = G_k_iom(k_t{0, 0}, matsubara_freq{0, beta, Fermion}); - auto r5 = sum(k_ >> G_k_iom(k_,0), g_eps.mesh()); - G_loc(w_) << sum(k_ >> G_k_iom(k_,w_), g_eps.mesh()); + auto r5 = sum_gf(k_ >> G_k_iom(k_,0), g_eps.mesh()); + G_loc(w_) << sum_gf(k_ >> G_k_iom(k_,w_), g_eps.mesh()); TEST(G_loc(0)); diff --git a/triqs/arrays/indexmaps/range.hpp b/triqs/arrays/indexmaps/range.hpp index e129829f..43ef283f 100644 --- a/triqs/arrays/indexmaps/range.hpp +++ b/triqs/arrays/indexmaps/range.hpp @@ -76,6 +76,12 @@ namespace arrays { const_iterator cend() const { return const_iterator(this, true); } }; + // foreach + template void foreach(range const& r, F const& f) { + std::ptrdiff_t i = r.first(), last = r.last(), step = r.step(); + for (; i < last; i += step) f(i); + } + /** */ class ellipsis : public range { diff --git a/triqs/clef/clef.hpp b/triqs/clef/clef.hpp index 257aafc1..16fb0cde 100644 --- a/triqs/clef/clef.hpp +++ b/triqs/clef/clef.hpp @@ -585,5 +585,31 @@ namespace triqs { namespace clef { \ template \ auto operator()(Args&&... args) && DECL_AND_RETURN(make_expr_call(std::move(*this), std::forward(args)...)); + +/* -------------------------------------------------------------------------------------------------- + * sum of expressions + * --------------------------------------------------------------------------------------------------- */ + + // sum a function f on a domain D, using a simple foreach + template + auto sum_f_domain_impl(F const& f, D const& d) + -> std::c14::enable_if_t::value, decltype(f(*(d.begin())))> { + auto res = decltype(f(*(d.begin()))) {}; + using p_t = typename std::decay::type; + foreach(d, [&res, &f](p_t const& x) { res = res + f(x); }); + return res; + } + + TRIQS_CLEF_MAKE_FNT_LAZY(sum_f_domain_impl); + + // sum( expression, i = domain) + template auto sum(Expr const& f, clef::pair const& d) + DECL_AND_RETURN(sum_f_domain_impl(make_function(f, clef::placeholder()), d.rhs)); + + // two or more indices : sum recursively + template + auto sum(Expr const& f, clef::pair const& d0, clef::pair const& d1, clef::pair const&... d) + DECL_AND_RETURN(sum(sum(f, d0), d1, d...)); + }} // namespace triqs::clef #endif