From 7aedaef945a13c3436090ad20a4bc2437bf4e4e2 Mon Sep 17 00:00:00 2001 From: Olivier Parcollet Date: Wed, 12 Feb 2014 21:50:58 +0100 Subject: [PATCH] arrays : fix matrix = scalar for complex types Pb : M() = rhs; // rhs of type RHS Currenlty does : M(i,j) = (i==j ? rhs : RHS{}) Changed to M(i,j) = (i==j ? rhs : RHS{0*rhs}) If RHS is a double, int ... Same result. If RHS is a matrix, gf, currently the offdiag elements are default constructed (i.e. of 0 size !). Which can break operations later (matrix>) After change : all elements have the same size ! --- triqs/arrays/impl/assignment.hpp | 6 +++++- 1 file changed, 5 insertions(+), 1 deletion(-) diff --git a/triqs/arrays/impl/assignment.hpp b/triqs/arrays/impl/assignment.hpp index 3b2901dd..6303ea22 100644 --- a/triqs/arrays/impl/assignment.hpp +++ b/triqs/arrays/impl/assignment.hpp @@ -139,8 +139,12 @@ namespace triqs { namespace arrays { typedef typename LHS::value_type value_type; LHS & lhs; const RHS & rhs; impl(LHS & lhs_, const RHS & rhs_): lhs(lhs_), rhs(rhs_){} //, p(*(lhs_.data_start())) {} + // we MUST make off_diag like this, if value_type is a complicated type (i.e. gf, matrix) with a size + // off diagonal element is 0*rhs, i.e. a 0, but with the SAME SIZE as the diagonal part. + // otherwise further operation may fail later. + // TO DO : look at performance issue ?? (we can remote the multiplication by 0 using an auxiliary function) template - void operator()(Args const & ... args) const {_ops_::invoke(lhs(args...), (kronecker(args...) ? rhs : RHS()));} + void operator()(Args const & ... args) const {_ops_::invoke(lhs(args...), (kronecker(args...) ? rhs : RHS{0*rhs}));} void invoke() { foreach(lhs,*this); } };