The following code block compiles fine, but throws a std::bad_alloc

Code: Select all

namespace {
  const unsigned int szx = 10, szy = 3;

Eigen::CwiseBinaryOp<Eigen::internal::scalar_sum_op<double>,
                     const Eigen::MatrixXd, const Eigen::MatrixXd>
get_x() {
  Eigen::MatrixXd A = Eigen::MatrixXd::Constant(szx, szy, 1);
  Eigen::MatrixXd B = Eigen::MatrixXd::Constant(szx, szy, 2);
  Eigen::CwiseBinaryOp<Eigen::internal::scalar_sum_op<double>,
                       const Eigen::MatrixXd, const Eigen::MatrixXd> op = A + B;
  return op;
}
}

int main(int, char * []) {
  Eigen::MatrixXd v = get_x();
}


I am writing expression templates that use Eigen as an underlying mechanism, so I would need to be able to pass these expressions around, and sometimes cast them to a matrix, like you would normally do in Eigen::MatrixXd(A+B), which works.

any insight much appreciated!

On GCC 4.4, it runs without throwing exceptions, but valgrind reports invalid reads. I do not know how these cause bad_alloc exceptions on your machine, but I think I can explain the invalid reads. The CwiseBinaryOperator you construct in get_x has references to A and B. After get_x returns to main, A and B disappear, so CwiseBinaryOperator now contains dangling references. In main, the CwiseBinaryOperator is cast to a MatrixXd, so it follows the dangling references and performs invalid reads.

This case can indeed be solved by storing the temporaries and, for example, passing them as const ref to the function. However, the following more elaborate example still fails although afaik everything is stored somewhere.

include file

Code: Select all

#ifndef PB_SIMPLE_H
#define PB_SIMPLE_H

#include <Eigen/Dense>

template <class EigenType>
class Holder {
   public:
    typedef EigenType result_type;

   private:
    result_type in_;

   public:
    Holder(const EigenType& in) : in_(in) {}
    result_type get() const { return in_; }
};

template <class HoldLeft, class HoldRight>
class Summer {
   public:
    typedef const typename Eigen::CwiseBinaryOp<
        Eigen::internal::scalar_sum_op<double>,
        const typename HoldLeft::result_type,
        const typename HoldRight::result_type> result_type;
    // typedef Eigen::MatrixXd result_type;
   private:
    HoldLeft left_;
    HoldRight right_;

   public:
    Summer(const HoldLeft& left, const HoldRight& right)
        : left_(left), right_(right) {}

    result_type get() const { return left_.get() + right_.get(); }
};

typedef Holder<Eigen::MatrixXd> MatrixHolder;
typedef Summer<MatrixHolder, MatrixHolder> MatrixSummer;

#endif /* PB_SIMPLE_H */

simple test

Code: Select all

#include "internal/PbSimple.h"

#include <Eigen/Dense>

int main(int, char * []) {
  const unsigned int szx=10,szy=3;
  Eigen::MatrixXd x(Eigen::MatrixXd::Constant(szx,szy,1));
  MatrixHolder vx(x);
  Eigen::MatrixXd y(Eigen::MatrixXd::Constant(szx,szy,2));
  MatrixHolder vy(y);
  MatrixSummer vsum(vx,vy);
  auto expr = vsum.get();
  MatrixHolder vz(expr); //force evaluation of sum into new matrix, fails here
  return 0;
}

In the include file, if you use the commented out typedef instead, it works fine.

I think the issue here is that when you call MatrixSummer::get(), it calls MatrixHolder::get() which returns a MatrixXd. This MatrixXd is a temporary; your definition of MatrixSummer::get() is equivalent (I think) to

Code: Select all

auto tmp1 = left_.get(); // in this case, tmp1 has type MatrixXd
auto tmp2 = right_.get();
return tmp1 + tmp2;


The CwiseBinaryOperator constructed has references to the temporary variables tmp1 and tmp2, which are deleted once the function exits. One possible solution would be to have MatrixHolder::get() return a constant reference to its member variable; that would make the program work (though it may well fail in other use cases you have in mind).

I find expression templates in C++ quite hard to understand and use, but it does teach you a lot about C++.

Thank you so much for this answer! I cross-posted this question on stackoverflow, if you want, you could post this answer and I'll mark it there as well.
http://stackoverflow.com/questions/2116 ... s-segfault

In addition, could you expand a little on the cases where you feel this would fail?