Hello, I'm using the gcc: *Version: 3.2, as well as 3.3 *System type: linux *options given when GCC was configured/built: Configured with: ../configure --prefix=/usr --libdir=/usr/lib --with-slibdir=/lib --mandir=/u sr/share/man --infodir=/usr/share/info --enable-shared --enable-threads=posi x --disable-checking --enable-long-long --enable-__cxa_atexit --enable-langu ages=c,c++,ada,f77,objc,java --host=i586-mandrake-linux-gnu --with-system-zl ib Thread model: posix gcc version 3.2 (Mandrake Linux 9.0 3.2-1mdk) * compilerflags were g++ -O2 g++ -O2 -W -Wall -Winline -finline-limit=5000 -ftemplate-depth-200 -fstrict-alia sing -malign-double I've got no error or warnings compiling the attached file. The program's output is: K = Matrix<d, 2, 2> = [ [6.25, -5001.25], [-3.75, 2506.25] ] Where K has the wrong result. It should be: K = Matrix<d, 2, 2> = [ [6.25, -3.75], [-3.75, 6.25] ] The problem is related to the temporary temp_lhs used by the prod function. As used here in the example, using RVO results in wrong results. Parts of the result are correct, other not. The algorithm/code used is correct (checked with a lot of regression tests). A simple product of matrizes (using prod(const Matrix<>&, const Matrix<>&) not shown here) got the expected result. The attached file uses (depends on define) dynamic allocated memory or static arrays. valgrind shows some "Invalid read of size 8" (see comments on the lines). The interesting behavior is by defined VMET_RVO_BUG_WO_STATIC_TEMP, where static temps are used for evaluating temporaries inside the prod() function - there are no errors found by valgrind! These looks for me as wrong code produced by gcc. Unfortunally I'm not able to understand the asm code produced. It doesn't seems to be a problem with the FPU - the same with long types too. Regards Olaf ----8<---- #include <iostream> #include <iomanip> #define TVMET_DYNAMIC_DATA #define TVMET_RVO_BUG_WO_STATIC_TEMP #ifndef _tvmet_restrict #define _tvmet_restrict __restrict__ #endif namespace tvmet { template<class T, std::size_t Rows, std::size_t Cols> class Matrix; template<class T1, class T2> struct PromoteTraits { }; template<> struct PromoteTraits<double, double> { typedef double value_type; }; class XprNull { XprNull& operator=(const XprNull&); public: explicit XprNull() { } }; template< class T > static inline T operator+(const T& lhs, XprNull) { return lhs; } // valgrind: Invalid read of size 8 namespace meta { template<std::size_t Rows, std::size_t Cols, std::size_t RowStride=0, std::size_t ColStride=0> class Matrix { private: enum { doRows = (RowStride < Rows - 1) ? 1 : 0, doCols = (ColStride < Cols - 1) ? 1 : 0 }; public: template<class Mtrx, class E, class Fcnl> static inline void assign2(Mtrx& mat, const E& expr, const Fcnl& fn) { fn.applyOn(mat(RowStride, ColStride), expr(RowStride, ColStride)); Matrix<Rows * doCols, Cols * doCols, RowStride * doCols, (ColStride+1) * doCols>::assign2(mat, expr, fn); } template<class Mtrx, class E, class Fcnl> static inline void assign(Mtrx& mat, const E& expr, const Fcnl& fn) { Matrix<Rows, Cols, RowStride, 0>::assign2(mat, expr, fn); Matrix<Rows * doRows, Cols * doRows, (RowStride+1) * doRows, 0>::assign(mat, expr, fn); } }; template<> class Matrix<0, 0, 0, 0> { public: template<class Mtrx, class E, class Fcnl> static inline void assign2(Mtrx&, const E&, const Fcnl&) { } template<class Mtrx, class E, class Fcnl> static inline void assign(Mtrx&, const E&, const Fcnl&) { } }; template<std::size_t Rows1, std::size_t Cols1, std::size_t Cols2, std::size_t RowStride1, std::size_t ColStride1, std::size_t RowStride2, std::size_t ColStride2, std::size_t K> class gemm { private: enum { doIt = (K != Cols1 - 1) }; public: template<class T1, class T2> static inline typename PromoteTraits<T1, T2>::value_type prod(const T1* _tvmet_restrict lhs, const T2* _tvmet_restrict rhs, std::size_t i, std::size_t j) { // valgrind: Invalid read of size 8 return lhs[i * RowStride1 + K * ColStride1] * rhs[K * RowStride2 + j * ColStride2] + gemm<Rows1 * doIt, Cols1 * doIt, Cols2 * doIt, RowStride1 * doIt, ColStride1 * doIt, RowStride2 * doIt, ColStride2 * doIt, (K+1) * doIt>::prod(lhs, rhs, i, j); } }; template<> class gemm<0,0,0,0,0,0,0,0> { public: static inline XprNull prod(const void*, const void*, std::size_t, std::size_t) { return XprNull(); } }; } // namespace meta template <class T1, class T2> struct fcnl_Assign { typedef void return_type; static inline return_type applyOn(T1& _tvmet_restrict lhs, T2 rhs) { lhs = static_cast<T1>(rhs); } }; template<class T1, std::size_t Rows1, std::size_t Cols1, class T2, std::size_t Cols2, std::size_t RowStride1, std::size_t ColStride1, std::size_t RowStride2, std::size_t ColStride2> class XprMMProduct { public: typedef typename PromoteTraits<T1, T2>::value_type value_type; public: explicit XprMMProduct(const T1* _tvmet_restrict lhs, const T2* _tvmet_restrict rhs) : m_lhs(lhs), m_rhs(rhs) { } value_type operator()(std::size_t i, std::size_t j) const { return meta::gemm<Rows1, Cols1, Cols2, RowStride1, ColStride1, RowStride2, ColStride2, 0>::prod(m_lhs, m_rhs, i, j); } private: const T1* _tvmet_restrict m_lhs; const T2* _tvmet_restrict m_rhs; }; template<class E> class XprMatrixTranspose { public: typedef E expr_type; typedef typename expr_type::value_type value_type; public: explicit XprMatrixTranspose(const expr_type& e) : m_expr(e) { } value_type operator()(std::size_t i, std::size_t j) const { return m_expr(j, i); } private: const expr_type& _tvmet_restrict m_expr; }; template<class E, std::size_t Rows, std::size_t Cols> class XprMatrix { public: typedef E expr_type; typedef typename expr_type::value_type value_type; public: explicit XprMatrix(const expr_type& e) : m_expr(e) { } value_type operator()(std::size_t i, std::size_t j) const { return m_expr(i, j); } private: const expr_type& _tvmet_restrict m_expr; }; template<class T, std::size_t Rows, std::size_t Cols, std::size_t RowStride=Cols, std::size_t ColStride=1> class MatrixConstReference { public: typedef T value_type; public: explicit MatrixConstReference(const Matrix<T, Rows, Cols>& rhs) : m_data(rhs.data()) { } value_type operator()(std::size_t i, std::size_t j) const { return m_data[i * RowStride + j * ColStride]; } private: const value_type* _tvmet_restrict m_data; }; template<class T, std::size_t Rows, std::size_t Cols> class Matrix { public: typedef T value_type; typedef Matrix<T, Rows, Cols> this_type; public: explicit Matrix() { #if defined(TVMET_DYNAMIC_DATA) m_data = new double [Rows*Cols]; #endif } template<class E> explicit Matrix(const XprMatrix<E, Rows, Cols>& expr) { #if defined(TVMET_DYNAMIC_DATA) m_data = new double [Rows*Cols]; #endif this->assign(expr, fcnl_Assign<value_type, typename E::value_type>()); } ~Matrix() { #if defined(TVMET_DYNAMIC_DATA) delete [] m_data; #endif } public: value_type* _tvmet_restrict data() { return m_data; } const value_type* _tvmet_restrict data() const { return m_data; } value_type& _tvmet_restrict operator()(std::size_t i, std::size_t j) { return m_data[i * Cols + j]; } value_type operator()(std::size_t i, std::size_t j) const { return m_data[i * Cols + j]; } public: typedef MatrixConstReference<T, Rows, Cols> ConstReference; ConstReference const_ref() const { return ConstReference(*this); } private: template<class E, class Fcnl> void assign(const E& expr, const Fcnl& fn) { meta::Matrix<Rows, Cols>::assign(*this, expr, fn); } public: template <class E> this_type& operator=(const XprMatrix<E, Rows, Cols>& rhs) { this->assign(rhs, fcnl_Assign<value_type, typename E::value_type>()); return *this; } public: std::ostream& print_on(std::ostream& os) const; private: #if defined(TVMET_DYNAMIC_DATA) value_type* m_data; #else value_type m_data[Rows*Cols]; #endif }; template<class T, std::size_t Rows, std::size_t Cols> inline std::ostream& Matrix<T, Rows, Cols>::print_on(std::ostream& os) const { std::streamsize w = os.width(); os << std::setw(0) << "Matrix<" << typeid(T).name() << ", " << Rows << ", " << Cols << "> = [\n"; for(std::size_t i = 0; i < Rows; ++i) { os << " ["; for(std::size_t j = 0; j < (Cols - 1); ++j) { os << std::setw(w) << this->operator()(i, j) << ", "; } os << std::setw(w) << this->operator()(i, Cols - 1) << (i != (Rows-1) ? "],\n" : "]\n"); } os << "]"; return os; } template<class T, std::size_t Rows, std::size_t Cols> inline std::ostream& operator<<(std::ostream& os, const Matrix<T, Rows, Cols>& rhs) { return rhs.print_on(os); } template<class E1, std::size_t Rows1, std::size_t Cols1, class T2, std::size_t Cols2> inline XprMatrix< XprMMProduct< typename E1::value_type, Rows1, Cols1, T2, Cols2, Cols1, 1, Cols2, 1 >, Rows1, Cols2 > prod(const XprMatrix<E1, Rows1, Cols1>& lhs, const Matrix<T2, Cols1, Cols2>& rhs) { typedef Matrix<typename E1::value_type, Rows1, Cols1> temp_matrix_type; typedef XprMMProduct< typename E1::value_type, Rows1, Cols1, T2, Cols2, Cols1, 1, Cols2, 1 > expr_type; #if defined(TVMET_RVO_BUG_WO_STATIC_TEMP) static temp_matrix_type temp_lhs(lhs); return XprMatrix<expr_type, Rows1, Cols2>(expr_type(temp_lhs.data(), rhs.data ())); #else return XprMatrix<expr_type, Rows1, Cols2>(expr_type(temp_matrix_type(lhs).data (), rhs.data())); #endif } template<class T, std::size_t Rows, std::size_t Cols> inline XprMatrix< XprMatrixTranspose< MatrixConstReference<T, Rows, Cols> >, Cols, Rows > trans(const Matrix<T, Rows, Cols>& rhs) { typedef XprMatrixTranspose< MatrixConstReference<T, Rows, Cols> > expr_type; return XprMatrix<expr_type, Cols, Rows>(expr_type(rhs.const_ref())); } } // namespace tvmet /** * Test driver */ using namespace std; int main() { tvmet::Matrix<double,3,2> B; tvmet::Matrix<double,3,3> D; tvmet::Matrix<double,2,2> K; B(0,0) = -0.05; B(0,1) = 0; B(1,0) = 0; B(1,1) = 0.05; B(2,0) = 0.05; B(2,1) = -0.05; D(0,0) = 2000; D(0,1) = 1000; D(0,2) = 0; D(1,0) = 1000; D(1,1) = 2000; D(1,2) = 0; D(2,0) = 0; D(2,1) = 0; D(2,2) = 500; K = prod(prod(trans(B), D), B); cout << "K = " << K << endl; } ---->8----
Olaf, thanks for your report. However, I think it will be some time until someone comes around and looks at it -- that code is just to complicated if you are not familiar with what it does. Can you try to reduce it somehow? Some suggestions are here: http://gcc.gnu.org/bugs/minimize.html In your case, leave in the includes, as they should be fine, but it would be tremendously helpful if you replaced classes by structs and removed accessor functions, reduced the number of template arguments by hardwiring their types or values, etc. Everything that makes the testcase shorter is ok. Basically, the shorter it is, the higher chances are that people are able to see what's supposed to happen and what is really going on. I think, looking at the code, that you are quite close to getting it down to a page or two, and that it should be possible in a reasonable amount of time if you know what the code does. Thanks Wolfgang
I cannot reproduce this on the mainline (20030920), 3.3.1 (20030707), 3.2.3, or 3.2.2.
But I can reproduce with the example you gave in <http://gcc.gnu.org/ml/gcc/2003-09/ msg00637.html>.
Hi, I've minimized the problematic code. The problem of minimizing is the use of templates and operator calls to evaluate templated expressions at compile time. A short introduction of technique used can be found at http://osl.iu.edu/~tveldhui/papers/Expression-Templates/exprtmpl.html and http://osl.iu.edu/~tveldhui/papers/Template-Metaprograms/meta-art.html . Valgrind reports invalid reads and the result isn't what expected. I haven't checked it with 3.3.1 (20030707), 3.2.3, or 3.2.2 but, 3.3.0 suffers too (A work arround is by defined TVMET_RVO_BUG_WO_STATIC_TEMP, where valgrind to errors found and the results are correct!) Well, the short version: ---8<--- extern "C" int printf(const char*, ...); //#define TVMET_RVO_BUG_WO_STATIC_TEMP #ifndef restrict #define restrict __restrict__ #endif template<unsigned Rows, unsigned Cols> class Matrix; struct XprNull { explicit XprNull() { } }; static inline double operator+(const double& lhs, XprNull) { return lhs; } // valgrind: Invalid read of size 8 struct fcnl_Assign { static inline void applyOn(double& restrict lhs, double rhs) { lhs = rhs; } }; template<unsigned Rows, unsigned Cols, unsigned RowStride, unsigned ColStride> struct MetaMatrix { enum { doRows = (RowStride < Rows - 1) ? 1 : 0, doCols = (ColStride < Cols - 1) ? 1 : 0 }; template<class Mtrx, class E, class Fcnl> static inline void assign2(Mtrx& mat, const E& expr, const Fcnl& fn) { fn.applyOn(mat(RowStride, ColStride), expr(RowStride, ColStride)); MetaMatrix<Rows * doCols, Cols * doCols, RowStride * doCols, (ColStride+1) * doCols>::assign2(mat, expr, fn); } template<class Mtrx, class E, class Fcnl> static inline void assign(Mtrx& mat, const E& expr, const Fcnl& fn) { MetaMatrix<Rows, Cols, RowStride, 0>::assign2(mat, expr, fn); MetaMatrix<Rows * doRows, Cols * doRows, (RowStride+1) * doRows, 0>::assign (mat, expr, fn); } }; template<> struct MetaMatrix<0, 0, 0, 0> { template<class Mtrx, class E, class Fcnl> static inline void assign2(Mtrx&, const E&, const Fcnl&) { } template<class Mtrx, class E, class Fcnl> static inline void assign(Mtrx&, const E&, const Fcnl&) { } }; template<unsigned Rows1, unsigned Cols1, unsigned Cols2, unsigned RowStride1, unsigned ColStride1, unsigned RowStride2, unsigned ColStride2, unsigned K> struct MetaGemm { enum { doIt = (K != Cols1 - 1) }; static inline double prod(const double* restrict lhs, const double* restrict rhs, unsigned i, unsigned j) { // valgrind: Invalid read of size 8 return lhs[i * RowStride1 + K * ColStride1] * rhs[K * RowStride2 + j * ColStride2] + MetaGemm<Rows1 * doIt, Cols1 * doIt, Cols2 * doIt, RowStride1 * doIt, ColStride1 * doIt, RowStride2 * doIt, ColStride2 * doIt, (K+1) * doIt>::prod(lhs, rhs, i, j); } }; template<> struct MetaGemm<0,0,0,0,0,0,0,0> { static inline XprNull prod(const void*, const void*, unsigned, unsigned) { return XprNull(); } }; template<unsigned Rows1, unsigned Cols1, unsigned Cols2, unsigned RowStride1, unsigned ColStride1, unsigned RowStride2, unsigned ColStride2> struct XprMMProduct { explicit XprMMProduct(const double* restrict lhs, const double* restrict rhs) : m_lhs(lhs), m_rhs(rhs) { } double operator()(unsigned i, unsigned j) const { return MetaGemm<Rows1, Cols1, Cols2, RowStride1, ColStride1, RowStride2, ColStride2, 0>::prod(m_lhs, m_rhs, i, j); } private: const double* restrict m_lhs; const double* restrict m_rhs; }; template<class E> struct XprMatrixTranspose { explicit XprMatrixTranspose(const E& e) : m_expr(e) { } double operator()(unsigned i, unsigned j) const { return m_expr(j, i); } private: const E& restrict m_expr; }; template<class E, unsigned Rows, unsigned Cols> struct XprMatrix { explicit XprMatrix(const E& e) : m_expr(e) { } double operator()(unsigned i, unsigned j) const { return m_expr(i, j); } private: const E& restrict m_expr; }; template<unsigned Rows, unsigned Cols, unsigned RowStride, unsigned ColStride> struct MatrixConstReference { explicit MatrixConstReference(const Matrix<Rows, Cols>& rhs) : m_data (rhs.m_data) { } double operator()(unsigned i, unsigned j) const { return m_data[i * RowStride + j * ColStride]; } private: const double* restrict m_data; }; template<unsigned Rows, unsigned Cols> struct Matrix { explicit Matrix() { m_data = new double [Rows*Cols]; } template<class E> explicit Matrix(const XprMatrix<E, Rows, Cols>& rhs) { m_data = new double [Rows*Cols]; MetaMatrix<Rows, Cols, 0, 0>::assign(*this, rhs, fcnl_Assign()); } ~Matrix() { delete [] m_data; } double& restrict operator()(unsigned i, unsigned j) { return m_data[i * Cols + j]; } double operator()(unsigned i, unsigned j) const { return m_data[i * Cols + j]; } MatrixConstReference<Rows,Cols,Cols,1> const_ref() const { return MatrixConstReference<Rows,Cols,Cols,1>(*this); } template <class E> Matrix& operator=(const XprMatrix<E, Rows, Cols>& rhs) { MetaMatrix<Rows, Cols, 0, 0>::assign(*this, rhs, fcnl_Assign()); return *this; } void print() const { printf("[\n"); for(unsigned i = 0; i != Rows; ++i) { printf("\t["); for(unsigned j = 0; j != Cols; ++j) printf("\t%+4.2f", this->operator()(i, j)); printf("]\n"); } printf("]\n"); } // private: double* m_data; }; template<class E1, unsigned Rows1, unsigned Cols1, unsigned Cols2> inline XprMatrix< XprMMProduct< Rows1, Cols1, Cols2, Cols1, 1, Cols2, 1 >, Rows1, Cols2 > prod(const XprMatrix<E1, Rows1, Cols1>& lhs, const Matrix<Cols1, Cols2>& rhs) { typedef XprMMProduct< Rows1, Cols1, Cols2, Cols1, 1, Cols2, 1 > expr_type; #if defined(TVMET_RVO_BUG_WO_STATIC_TEMP) static Matrix<Rows1, Cols1> temp_lhs(lhs); return XprMatrix<expr_type, Rows1, Cols2>(expr_type(temp_lhs.m_data, rhs.m_data)); #else return XprMatrix<expr_type, Rows1, Cols2>(expr_type(Matrix<Rows1, Cols1> (lhs).m_data, rhs.m_data)); #endif } template<unsigned Rows, unsigned Cols> inline XprMatrix< XprMatrixTranspose< MatrixConstReference<Rows, Cols, Cols, 1> >, Cols, Rows > trans(const Matrix<Rows, Cols>& rhs) { typedef XprMatrixTranspose< MatrixConstReference<Rows, Cols, Cols, 1> > expr_type; return XprMatrix<expr_type, Cols, Rows>(expr_type(rhs.const_ref())); } /** * Test driver */ using namespace std; int main() { Matrix<3,2> B; Matrix<3,3> D; Matrix<2,2> K; B(0,0) = -0.05; B(0,1) = 0; B(1,0) = 0; B(1,1) = 0.05; B(2,0) = 0.05; B(2,1) = -0.05; D(0,0) = 2000; D(0,1) = 1000; D(0,2) = 0; D(1,0) = 1000; D(1,1) = 2000; D(1,2) = 0; D(2,0) = 0; D(2,1) = 0; D(2,2) = 500; K = prod(prod(trans(B), D), B); printf("K = "); K.print(); // wrong result, should be symetric } --->8---
I believe the code to be invalid. Can someone confirm this please? Explanation: We have in function prod: template<...> inline XprMatrix<...> prod(const XprMatrix<...>& lhs, const Matrix<...>& rhs) { typedef XprMMProduct<...> expr_type; return XprMatrix<expr_type, ...> ( expr_type (...) ); } So this constructs a temporary object of type epxr_type (actually XprMMProduct<...>) and uses this temporary object to construct another temporary object of type XprMatrix. The second temporary is then returned by value. Now XprMatrix looks like this: template<class E, unsigned Rows, unsigned Cols> struct XprMatrix { explicit XprMatrix(const E& e) : m_expr(e) { } [ ... ] const E& restrict m_expr; }; This means that m_expr is a direct reference to the temprory of type XprMMProduct<...> created in function prod. However, the life time of that temporary ends at exit of function prod. Thus prod returns an object containing a reference to a tempory object that is dead after prod returns. This invokes undefined behaviour. regards Christian
yes, the lifetime ends at the end of the return statement
Hello, how can I solve the problem? A simple change of XprMatrix's private data to: template<class E, unsigned Rows, unsigned Cols> struct XprMatrix { explicit XprMatrix(const E& e) : m_expr(e) { } [ ... ] E m_expr; }; doesn't solve the problem. The prod function uses XprMatrix's copy constructor without a reference now - or I'm wrong? Thanks Olaf