-new: recursive costs system, useful to determine automatically

 when to evaluate arguments and when to meta-unroll.
-use it in Product to determine when to eval args. not yet used
 to determine when to unroll. for now, not used anywhere else but
 that'll follow.
-fix badness of my last commit

28 files changed, 379 insertions, 212 deletions

diff --git a/Eigen/src/Core/AssociativeFunctors.h b/Eigen/src/Core/AssociativeFunctors.h
index fd42de9ec..3b3a09ed9 100644
--- a/Eigen/src/Core/AssociativeFunctors.h
+++ b/Eigen/src/Core/AssociativeFunctors.h
@@ -30,8 +30,9 @@
   *
   * \sa class CwiseBinaryOp, MatrixBase::operator+, class PartialRedux, MatrixBase::sum()
   */
-struct ei_scalar_sum_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return a + b; }
+template<typename Scalar> struct ei_scalar_sum_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a, const Scalar& b) const { return a + b; }
+  enum { Cost = NumTraits<Scalar>::AddCost };
 };
 
 /** \internal
@@ -39,8 +40,9 @@ struct ei_scalar_sum_op EIGEN_EMPTY_STRUCT {
   *
   * \sa class CwiseBinaryOp, MatrixBase::cwiseProduct(), class PartialRedux, MatrixBase::redux()
   */
-struct ei_scalar_product_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return a * b; }
+template<typename Scalar> struct ei_scalar_product_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a, const Scalar& b) const { return a * b; }
+  enum { Cost = NumTraits<Scalar>::MulCost };
 };
 
 /** \internal
@@ -48,8 +50,9 @@ struct ei_scalar_product_op EIGEN_EMPTY_STRUCT {
   *
   * \sa class CwiseBinaryOp, MatrixBase::cwiseMin, class PartialRedux, MatrixBase::minCoeff()
   */
-struct ei_scalar_min_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return std::min(a, b); }
+template<typename Scalar> struct ei_scalar_min_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a, const Scalar& b) const { return std::min(a, b); }
+  enum { Cost = ConditionalJumpCost + NumTraits<Scalar>::AddCost };
 };
 
 /** \internal
@@ -57,8 +60,9 @@ struct ei_scalar_min_op EIGEN_EMPTY_STRUCT {
   *
   * \sa class CwiseBinaryOp, MatrixBase::cwiseMax, class PartialRedux, MatrixBase::maxCoeff()
   */
-struct ei_scalar_max_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return std::max(a, b); }
+template<typename Scalar> struct ei_scalar_max_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a, const Scalar& b) const { return std::max(a, b); }
+  enum { Cost = ConditionalJumpCost + NumTraits<Scalar>::AddCost };
 };
 
 #endif // EIGEN_ASSOCIATIVE_FUNCTORS_H
diff --git a/Eigen/src/Core/Block.h b/Eigen/src/Core/Block.h
index 4c672e4f2..34f98030e 100644
--- a/Eigen/src/Core/Block.h
+++ b/Eigen/src/Core/Block.h
@@ -69,7 +69,8 @@ struct ei_traits<Block<MatrixType, BlockRows, BlockCols> >
       : (BlockCols==Dynamic ? MatrixType::MaxColsAtCompileTime : BlockCols),
     Flags = RowsAtCompileTime == Dynamic || ColsAtCompileTime == Dynamic
             ? (unsigned int)MatrixType::Flags
-            : (unsigned int)MatrixType::Flags &~ LargeBit
+            : (unsigned int)MatrixType::Flags &~ LargeBit,
+    CoeffReadCost = MatrixType::CoeffReadCost
   };
 };
 
diff --git a/Eigen/src/Core/CommaInitializer.h b/Eigen/src/Core/CommaInitializer.h
index 656ee3002..d7a5fcb47 100644
--- a/Eigen/src/Core/CommaInitializer.h
+++ b/Eigen/src/Core/CommaInitializer.h
@@ -23,8 +23,8 @@
 // License and a copy of the GNU General Public License along with
 // Eigen. If not, see <http://www.gnu.org/licenses/>.
 
-#ifndef EIGEN_COMMA_INITIALIZER_H
-#define EIGEN_COMMA_INITIALIZER_H
+#ifndef EIGEN_COMMAINITIALIZER_H
+#define EIGEN_COMMAINITIALIZER_H
 
 /** \internal
   * Helper class to define the MatrixBase::operator<<
@@ -120,4 +120,4 @@ MatrixBase<Derived>::operator<<(const MatrixBase<OtherDerived>& other)
   return CommaInitializer(*static_cast<Derived *>(this), other);
 }
 
-#endif // EIGEN_COMMA_INITIALIZER_H
+#endif // EIGEN_COMMAINITIALIZER_H
diff --git a/Eigen/src/Core/CwiseBinaryOp.h b/Eigen/src/Core/CwiseBinaryOp.h
index 9a83f4257..12da94a6c 100644
--- a/Eigen/src/Core/CwiseBinaryOp.h
+++ b/Eigen/src/Core/CwiseBinaryOp.h
@@ -60,7 +60,8 @@ struct ei_traits<CwiseBinaryOp<BinaryOp, Lhs, Rhs> >
     ColsAtCompileTime = Lhs::ColsAtCompileTime,
     MaxRowsAtCompileTime = Lhs::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = Lhs::MaxColsAtCompileTime,
-    Flags = Lhs::Flags | Rhs::Flags
+    Flags = Lhs::Flags | Rhs::Flags,
+    CoeffReadCost = Lhs::CoeffReadCost + Rhs::CoeffReadCost + BinaryOp::Cost
   };
 };
 
@@ -99,8 +100,9 @@ class CwiseBinaryOp : ei_no_assignment_operator,
   *
   * \sa class CwiseBinaryOp, MatrixBase::operator-
   */
-struct ei_scalar_difference_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return a - b; }
+template<typename Scalar> struct ei_scalar_difference_op EIGEN_EMPTY_STRUCT {
+    const Scalar operator() (const Scalar& a, const Scalar& b) const { return a - b; }
+    enum { Cost = NumTraits<Scalar>::AddCost };
 };
 
 /** \internal
@@ -108,8 +110,9 @@ struct ei_scalar_difference_op EIGEN_EMPTY_STRUCT {
   *
   * \sa class CwiseBinaryOp, MatrixBase::cwiseQuotient()
   */
-struct ei_scalar_quotient_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a, const Scalar& b) const { return a / b; }
+template<typename Scalar> struct ei_scalar_quotient_op EIGEN_EMPTY_STRUCT {
+    const Scalar operator() (const Scalar& a, const Scalar& b) const { return a / b; }
+    enum { Cost = 2 * NumTraits<Scalar>::MulCost };
 };
 
 /**\returns an expression of the difference of \c *this and \a other
@@ -118,10 +121,12 @@ struct ei_scalar_quotient_op EIGEN_EMPTY_STRUCT {
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_difference_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_difference_op<typename ei_traits<Derived>::Scalar>,
+                    Derived, OtherDerived>
 MatrixBase<Derived>::operator-(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_difference_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_difference_op<Scalar>,
+                       Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** replaces \c *this by \c *this - \a other.
@@ -144,10 +149,10 @@ MatrixBase<Derived>::operator-=(const MatrixBase<OtherDerived> &other)
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_sum_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_sum_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
 MatrixBase<Derived>::operator+(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_sum_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_sum_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** replaces \c *this by \c *this + \a other.
@@ -168,10 +173,10 @@ MatrixBase<Derived>::operator+=(const MatrixBase<OtherDerived>& other)
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_product_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_product_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
 MatrixBase<Derived>::cwiseProduct(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_product_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_product_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** \returns an expression of the coefficient-wise quotient of *this and \a other
@@ -180,10 +185,10 @@ MatrixBase<Derived>::cwiseProduct(const MatrixBase<OtherDerived> &other) const
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_quotient_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_quotient_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
 MatrixBase<Derived>::cwiseQuotient(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_quotient_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_quotient_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** \returns an expression of the coefficient-wise min of *this and \a other
@@ -192,10 +197,10 @@ MatrixBase<Derived>::cwiseQuotient(const MatrixBase<OtherDerived> &other) const
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_min_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_min_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
 MatrixBase<Derived>::cwiseMin(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_min_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_min_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** \returns an expression of the coefficient-wise max of *this and \a other
@@ -204,10 +209,10 @@ MatrixBase<Derived>::cwiseMin(const MatrixBase<OtherDerived> &other) const
   */
 template<typename Derived>
 template<typename OtherDerived>
-const CwiseBinaryOp<ei_scalar_max_op, Derived, OtherDerived>
+const CwiseBinaryOp<ei_scalar_max_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
 MatrixBase<Derived>::cwiseMax(const MatrixBase<OtherDerived> &other) const
 {
-  return CwiseBinaryOp<ei_scalar_max_op, Derived, OtherDerived>(derived(), other.derived());
+  return CwiseBinaryOp<ei_scalar_max_op<Scalar>, Derived, OtherDerived>(derived(), other.derived());
 }
 
 /** \returns an expression of a custom coefficient-wise operator \a func of *this and \a other
diff --git a/Eigen/src/Core/CwiseUnaryOp.h b/Eigen/src/Core/CwiseUnaryOp.h
index d3cc6ad30..548d500a1 100644
--- a/Eigen/src/Core/CwiseUnaryOp.h
+++ b/Eigen/src/Core/CwiseUnaryOp.h
@@ -50,7 +50,8 @@ struct ei_traits<CwiseUnaryOp<UnaryOp, MatrixType> >
     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = MatrixType::CoeffReadCost + UnaryOp::Cost
   };
 };
 
@@ -103,17 +104,18 @@ MatrixBase<Derived>::cwise(const CustomUnaryOp& func) const
   *
   * \sa class CwiseUnaryOp, MatrixBase::operator-
   */
-struct ei_scalar_opposite_op EIGEN_EMPTY_STRUCT {
-  template<typename Scalar> Scalar operator() (const Scalar& a) const { return -a; }
+template<typename Scalar> struct ei_scalar_opposite_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return -a; }
+  enum { Cost = NumTraits<Scalar>::AddCost };
 };
 
 /** \returns an expression of the opposite of \c *this
   */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_opposite_op,Derived>
+const CwiseUnaryOp<ei_scalar_opposite_op<typename ei_traits<Derived>::Scalar>,Derived>
 MatrixBase<Derived>::operator-() const
 {
-  return CwiseUnaryOp<ei_scalar_opposite_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_opposite_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -121,17 +123,18 @@ MatrixBase<Derived>::operator-() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseAbs
   */
-struct ei_scalar_abs_op EIGEN_EMPTY_STRUCT {
-  template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_abs(a); }
+template<typename Scalar> struct ei_scalar_abs_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_abs(a); }
+  enum { Cost = NumTraits<Scalar>::AddCost };
 };
 
 /** \returns an expression of the coefficient-wise absolute value of \c *this
   */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_abs_op,Derived>
+const CwiseUnaryOp<ei_scalar_abs_op<typename ei_traits<Derived>::Scalar>,Derived>
 MatrixBase<Derived>::cwiseAbs() const
 {
-  return CwiseUnaryOp<ei_scalar_abs_op,Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_abs_op<Scalar>,Derived>(derived());
 }
 
 /** \internal
@@ -139,17 +142,18 @@ MatrixBase<Derived>::cwiseAbs() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseAbs2
   */
-struct ei_scalar_abs2_op EIGEN_EMPTY_STRUCT {
-  template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_abs2(a); }
+template<typename Scalar> struct ei_scalar_abs2_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_abs2(a); }
+  enum { Cost = NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise squared absolute value of \c *this
   */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_abs2_op,Derived>
+const CwiseUnaryOp<ei_scalar_abs2_op<typename ei_traits<Derived>::Scalar>,Derived>
 MatrixBase<Derived>::cwiseAbs2() const
 {
-  return CwiseUnaryOp<ei_scalar_abs2_op,Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_abs2_op<Scalar>,Derived>(derived());
 }
 
 /** \internal
@@ -157,18 +161,19 @@ MatrixBase<Derived>::cwiseAbs2() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::conjugate()
   */
-struct ei_scalar_conjugate_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_conj(a); }
+template<typename Scalar> struct ei_scalar_conjugate_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_conj(a); }
+  enum { Cost = NumTraits<Scalar>::IsComplex ? NumTraits<Scalar>::AddCost : 0 };
 };
 
 /** \returns an expression of the complex conjugate of *this.
   *
   * \sa adjoint() */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_conjugate_op, Derived>
+const CwiseUnaryOp<ei_scalar_conjugate_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::conjugate() const
 {
-  return CwiseUnaryOp<ei_scalar_conjugate_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -176,10 +181,11 @@ MatrixBase<Derived>::conjugate() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cast()
   */
-template<typename NewType>
+template<typename Scalar, typename NewType>
 struct ei_scalar_cast_op EIGEN_EMPTY_STRUCT {
-    typedef NewType result_type;
-    template<typename Scalar> NewType operator() (const Scalar& a) const { return static_cast<NewType>(a); }
+  typedef NewType result_type;
+  const NewType operator() (const Scalar& a) const { return static_cast<NewType>(a); }
+  enum { Cost = ei_is_same_type<Scalar, NewType>::ret ? 0 : NumTraits<NewType>::AddCost };
 };
 
 /** \returns an expression of *this with the \a Scalar type casted to
@@ -191,10 +197,10 @@ struct ei_scalar_cast_op EIGEN_EMPTY_STRUCT {
   */
 template<typename Derived>
 template<typename NewType>
-const CwiseUnaryOp<ei_scalar_cast_op<NewType>, Derived>
+const CwiseUnaryOp<ei_scalar_cast_op<typename ei_traits<Derived>::Scalar, NewType>, Derived>
 MatrixBase<Derived>::cast() const
 {
-  return CwiseUnaryOp<ei_scalar_cast_op<NewType>, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_cast_op<Scalar, NewType>, Derived>(derived());
 }
 
 /** \internal
@@ -204,23 +210,26 @@ MatrixBase<Derived>::cast() const
   */
 template<typename Scalar>
 struct ei_scalar_multiple_op {
-    ei_scalar_multiple_op(const Scalar& other) : m_other(other) {}
-    Scalar operator() (const Scalar& a) const { return a * m_other; }
-    const Scalar m_other;
+  ei_scalar_multiple_op(const Scalar& other) : m_other(other) {}
+  Scalar operator() (const Scalar& a) const { return a * m_other; }
+  const Scalar m_other;
+  enum { Cost = NumTraits<Scalar>::MulCost };
 };
 
 template<typename Scalar, bool HasFloatingPoint>
 struct ei_scalar_quotient1_impl {
-    ei_scalar_quotient1_impl(const Scalar& other) : m_other(static_cast<Scalar>(1) / other) {}
-    Scalar operator() (const Scalar& a) const { return a * m_other; }
-    const Scalar m_other;
+  ei_scalar_quotient1_impl(const Scalar& other) : m_other(static_cast<Scalar>(1) / other) {}
+  Scalar operator() (const Scalar& a) const { return a * m_other; }
+  const Scalar m_other;
+  enum { Cost = NumTraits<Scalar>::MulCost };
 };
 
 template<typename Scalar>
 struct ei_scalar_quotient1_impl<Scalar,false> {
-    ei_scalar_quotient1_impl(const Scalar& other) : m_other(other) {}
-    Scalar operator() (const Scalar& a) const { return a / m_other; }
-    const Scalar m_other;
+  ei_scalar_quotient1_impl(const Scalar& other) : m_other(other) {}
+  Scalar operator() (const Scalar& a) const { return a / m_other; }
+  const Scalar m_other;
+  enum { Cost = 2 * NumTraits<Scalar>::MulCost };
 };
 
 /** \internal
@@ -274,16 +283,17 @@ MatrixBase<Derived>::operator/=(const Scalar& other)
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseSqrt()
   */
-struct ei_scalar_sqrt_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_sqrt(a); }
+template<typename Scalar> struct ei_scalar_sqrt_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_sqrt(a); }
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise square root of *this. */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_sqrt_op, Derived>
+const CwiseUnaryOp<ei_scalar_sqrt_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::cwiseSqrt() const
 {
-  return CwiseUnaryOp<ei_scalar_sqrt_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_sqrt_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -291,16 +301,17 @@ MatrixBase<Derived>::cwiseSqrt() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseExp()
   */
-struct ei_scalar_exp_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_exp(a); }
+template<typename Scalar> struct ei_scalar_exp_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_exp(a); }
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise exponential of *this. */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_exp_op, Derived>
+const CwiseUnaryOp<ei_scalar_exp_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::cwiseExp() const
 {
-  return CwiseUnaryOp<ei_scalar_exp_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_exp_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -308,16 +319,17 @@ MatrixBase<Derived>::cwiseExp() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseLog()
   */
-struct ei_scalar_log_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_log(a); }
+template<typename Scalar> struct ei_scalar_log_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_log(a); }
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise logarithm of *this. */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_log_op, Derived>
+const CwiseUnaryOp<ei_scalar_log_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::cwiseLog() const
 {
-  return CwiseUnaryOp<ei_scalar_log_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_log_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -325,16 +337,17 @@ MatrixBase<Derived>::cwiseLog() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseCos()
   */
-struct ei_scalar_cos_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_cos(a); }
+template<typename Scalar> struct ei_scalar_cos_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_cos(a); }
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise cosine of *this. */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_cos_op, Derived>
+const CwiseUnaryOp<ei_scalar_cos_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::cwiseCos() const
 {
-  return CwiseUnaryOp<ei_scalar_cos_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_cos_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -342,16 +355,17 @@ MatrixBase<Derived>::cwiseCos() const
   *
   * \sa class CwiseUnaryOp, MatrixBase::cwiseSin()
   */
-struct ei_scalar_sin_op EIGEN_EMPTY_STRUCT {
-    template<typename Scalar> Scalar operator() (const Scalar& a) const { return ei_sin(a); }
+template<typename Scalar> struct ei_scalar_sin_op EIGEN_EMPTY_STRUCT {
+  const Scalar operator() (const Scalar& a) const { return ei_sin(a); }
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \returns an expression of the coefficient-wise sine of *this. */
 template<typename Derived>
-const CwiseUnaryOp<ei_scalar_sin_op, Derived>
+const CwiseUnaryOp<ei_scalar_sin_op<typename ei_traits<Derived>::Scalar>, Derived>
 MatrixBase<Derived>::cwiseSin() const
 {
-  return CwiseUnaryOp<ei_scalar_sin_op, Derived>(derived());
+  return CwiseUnaryOp<ei_scalar_sin_op<Scalar>, Derived>(derived());
 }
 
 /** \internal
@@ -361,9 +375,10 @@ MatrixBase<Derived>::cwiseSin() const
   */
 template<typename Scalar>
 struct ei_scalar_pow_op {
-    ei_scalar_pow_op(const Scalar& exponent) : m_exponent(exponent) {}
-    Scalar operator() (const Scalar& a) const { return ei_pow(a, m_exponent); }
-    const Scalar m_exponent;
+  ei_scalar_pow_op(const Scalar& exponent) : m_exponent(exponent) {}
+  Scalar operator() (const Scalar& a) const { return ei_pow(a, m_exponent); }
+  const Scalar m_exponent;
+  enum { Cost = 5 * NumTraits<Scalar>::MulCost };
 };
 
 /** \relates MatrixBase */
diff --git a/Eigen/src/Core/DiagonalCoeffs.h b/Eigen/src/Core/DiagonalCoeffs.h
index 417daa8ce..030de5cf0 100644
--- a/Eigen/src/Core/DiagonalCoeffs.h
+++ b/Eigen/src/Core/DiagonalCoeffs.h
@@ -52,9 +52,10 @@ struct ei_traits<DiagonalCoeffs<MatrixType> >
                             : EIGEN_ENUM_MIN(MatrixType::MaxRowsAtCompileTime,
                                              MatrixType::MaxColsAtCompileTime),
     MaxColsAtCompileTime = 1,
-    Flags = RowsAtCompileTime == Dynamic || ColsAtCompileTime == Dynamic
+    Flags = RowsAtCompileTime == Dynamic && ColsAtCompileTime == Dynamic
             ? (unsigned int)MatrixType::Flags
-            : (unsigned int)MatrixType::Flags &~ LargeBit
+            : (unsigned int)MatrixType::Flags &~ LargeBit,
+    CoeffReadCost = MatrixType::CoeffReadCost
   };
 };
 
diff --git a/Eigen/src/Core/DiagonalMatrix.h b/Eigen/src/Core/DiagonalMatrix.h
index f46e8ad0d..b7fffab72 100644
--- a/Eigen/src/Core/DiagonalMatrix.h
+++ b/Eigen/src/Core/DiagonalMatrix.h
@@ -47,7 +47,8 @@ struct ei_traits<DiagonalMatrix<CoeffsVectorType> >
     ColsAtCompileTime = CoeffsVectorType::SizeAtCompileTime,
     MaxRowsAtCompileTime = CoeffsVectorType::MaxSizeAtCompileTime,
     MaxColsAtCompileTime = CoeffsVectorType::MaxSizeAtCompileTime,
-    Flags = CoeffsVectorType::Flags
+    Flags = CoeffsVectorType::Flags,
+    CoeffReadCost = CoeffsVectorType::CoeffReadCost
   };
 };
 
diff --git a/Eigen/src/Core/ForwardDeclarations.h b/Eigen/src/Core/ForwardDeclarations.h
index 58e1bf353..50cffe637 100644
--- a/Eigen/src/Core/ForwardDeclarations.h
+++ b/Eigen/src/Core/ForwardDeclarations.h
@@ -36,7 +36,7 @@ template<typename MatrixType> class Transpose;
 template<typename MatrixType> class Conjugate;
 template<typename BinaryOp, typename Lhs, typename Rhs> class CwiseBinaryOp;
 template<typename UnaryOp, typename MatrixType> class CwiseUnaryOp;
-template<typename Lhs, typename Rhs, int EvalMode=ei_product_eval_mode<Lhs,Rhs>::EvalMode > class Product;
+template<typename Lhs, typename Rhs, int EvalMode=ei_product_eval_mode<Lhs,Rhs>::value> class Product;
 template<typename MatrixType> class Random;
 template<typename MatrixType> class Zero;
 template<typename MatrixType> class Ones;
@@ -48,55 +48,53 @@ template<typename Derived> class Eval;
 template<typename Derived> class EvalOMP;
 template<int Direction, typename UnaryOp, typename MatrixType> class PartialRedux;
 
-struct ei_scalar_sum_op;
-struct ei_scalar_difference_op;
-struct ei_scalar_product_op;
-struct ei_scalar_quotient_op;
-struct ei_scalar_opposite_op;
-struct ei_scalar_conjugate_op;
-struct ei_scalar_abs_op;
-struct ei_scalar_abs2_op;
-struct ei_scalar_sqrt_op;
-struct ei_scalar_exp_op;
-struct ei_scalar_log_op;
-struct ei_scalar_cos_op;
-struct ei_scalar_sin_op;
-template<typename Scalar>  struct ei_scalar_pow_op;
-template<typename NewType> struct ei_scalar_cast_op;
-template<typename Scalar>  struct ei_scalar_multiple_op;
-template<typename Scalar>  struct ei_scalar_quotient1_op;
-struct ei_scalar_min_op;
-struct ei_scalar_max_op;
+template<typename Scalar> struct ei_scalar_sum_op;
+template<typename Scalar> struct ei_scalar_difference_op;
+template<typename Scalar> struct ei_scalar_product_op;
+template<typename Scalar> struct ei_scalar_quotient_op;
+template<typename Scalar> struct ei_scalar_opposite_op;
+template<typename Scalar> struct ei_scalar_conjugate_op;
+template<typename Scalar> struct ei_scalar_abs_op;
+template<typename Scalar> struct ei_scalar_abs2_op;
+template<typename Scalar> struct ei_scalar_sqrt_op;
+template<typename Scalar> struct ei_scalar_exp_op;
+template<typename Scalar> struct ei_scalar_log_op;
+template<typename Scalar> struct ei_scalar_cos_op;
+template<typename Scalar> struct ei_scalar_sin_op;
+template<typename Scalar> struct ei_scalar_pow_op;
+template<typename Scalar, typename NewType> struct ei_scalar_cast_op;
+template<typename Scalar> struct ei_scalar_multiple_op;
+template<typename Scalar> struct ei_scalar_quotient1_op;
+template<typename Scalar> struct ei_scalar_min_op;
+template<typename Scalar> struct ei_scalar_max_op;
 
 template<typename T> struct ei_xpr_copy
 {
-  typedef T Type;
+  typedef T type;
 };
 
 template<typename _Scalar, int _Rows, int _Cols, unsigned int _Flags, int _MaxRows, int _MaxCols>
 struct ei_xpr_copy<Matrix<_Scalar, _Rows, _Cols, _Flags, _MaxRows, _MaxCols> >
 {
-  typedef const Matrix<_Scalar, _Rows, _Cols, _Flags, _MaxRows, _MaxCols> & Type;
+  typedef const Matrix<_Scalar, _Rows, _Cols, _Flags, _MaxRows, _MaxCols> & type;
 };
 
-template<typename T, bool value> struct ei_conditional_eval
-{
-  typedef T Type;
-};
-
-template<typename T> struct ei_conditional_eval<T, true>
+template<typename T> struct ei_eval
 {
   typedef Matrix<typename ei_traits<T>::Scalar,
                  ei_traits<T>::RowsAtCompileTime,
                  ei_traits<T>::ColsAtCompileTime,
-                 ei_traits<T>::Flags,
+                 ei_traits<T>::Flags & ~LazyBit, // unset lazy bit after evaluation
                  ei_traits<T>::MaxRowsAtCompileTime,
-                 ei_traits<T>::MaxColsAtCompileTime> Type;
+                 ei_traits<T>::MaxColsAtCompileTime> type;
 };
 
 template<typename T> struct ei_eval_unless_lazy
 {
-  typedef typename ei_conditional_eval<T, !(ei_traits<T>::Flags & LazyBit)>::Type Type;
+  typedef typename ei_meta_if<ei_traits<T>::Flags & LazyBit,
+                              T,
+                              typename ei_eval<T>::type
+                             >::ret type;
 };
 
 #endif // EIGEN_FORWARDDECLARATIONS_H
diff --git a/Eigen/src/Core/Identity.h b/Eigen/src/Core/Identity.h
index 532167b9c..104a06e2f 100644
--- a/Eigen/src/Core/Identity.h
+++ b/Eigen/src/Core/Identity.h
@@ -40,7 +40,8 @@ struct ei_traits<Identity<MatrixType> >
     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = NumTraits<Scalar>::ReadCost
   };
 };
 
diff --git a/Eigen/src/Core/Lazy.h b/Eigen/src/Core/Lazy.h
new file mode 100644
index 000000000..0968e254f
--- /dev/null
+++ b/Eigen/src/Core/Lazy.h
@@ -0,0 +1,91 @@
+// This file is part of Eigen, a lightweight C++ template library
+// for linear algebra. Eigen itself is part of the KDE project.
+//
+// Copyright (C) 2006-2008 Benoit Jacob <jacob@math.jussieu.fr>
+//
+// Eigen is free software; you can redistribute it and/or
+// modify it under the terms of the GNU Lesser General Public
+// License as published by the Free Software Foundation; either
+// version 3 of the License, or (at your option) any later version.
+//
+// Alternatively, you can redistribute it and/or
+// modify it under the terms of the GNU General Public License as
+// published by the Free Software Foundation; either version 2 of
+// the License, or (at your option) any later version.
+//
+// Eigen is distributed in the hope that it will be useful, but WITHOUT ANY
+// WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
+// FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License or the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU Lesser General Public
+// License and a copy of the GNU General Public License along with
+// Eigen. If not, see <http://www.gnu.org/licenses/>.
+
+#ifndef EIGEN_LAZY_H
+#define EIGEN_LAZY_H
+
+/** \class Lazy
+  *
+  * \brief Expression with the lazy flag set
+  *
+  * \param ExpressionType the type of the object of which we are taking the lazy version
+  *
+  * This class represents the lazy version of an expression.
+  * It is the return type of MatrixBase::lazy()
+  * and most of the time this is the only way it is used.
+  *
+  * \sa MatrixBase::lazy()
+  */
+template<typename ExpressionType>
+struct ei_traits<Lazy<ExpressionType> >
+{
+  typedef typename ExpressionType::Scalar Scalar;
+  enum {
+    RowsAtCompileTime = ExpressionType::RowsAtCompileTime,
+    ColsAtCompileTime = ExpressionType::ColsAtCompileTime,
+    MaxRowsAtCompileTime = ExpressionType::MaxRowsAtCompileTime,
+    MaxColsAtCompileTime = ExpressionType::MaxColsAtCompileTime,
+    Flags = ExpressionType::Flags | LazyBit,
+    CoeffReadCost = ExpressionType::CoeffReadCost
+  };
+};
+
+template<typename ExpressionType> class Lazy
+  : public MatrixBase<Lazy<ExpressionType> >
+{
+  public:
+
+    EIGEN_GENERIC_PUBLIC_INTERFACE(Lazy)
+
+    Lazy(const ExpressionType& matrix) : m_expression(matrix) {}
+
+    EIGEN_INHERIT_ASSIGNMENT_OPERATORS(Lazy)
+
+  private:
+
+    int _rows() const { return m_expression.rows(); }
+    int _cols() const { return m_expression.cols(); }
+
+    const Scalar _coeff(int row, int col) const
+    {
+      return m_expression.coeff(row, col);
+    }
+
+  protected:
+    const typename ExpressionType::XprCopy m_expression;
+};
+
+/** \returns an expression of the lazy version of *this.
+  *
+  * Example: \include MatrixBase_lazy.cpp
+  * Output: \verbinclude MatrixBase_lazy.out
+  */
+template<typename Derived>
+const Lazy<Derived>
+MatrixBase<Derived>::lazy() const
+{
+  return Lazy<Derived>(derived());
+}
+
+#endif // EIGEN_LAZY_H
diff --git a/Eigen/src/Core/Map.h b/Eigen/src/Core/Map.h
index 4872b3f60..cbb1633ad 100644
--- a/Eigen/src/Core/Map.h
+++ b/Eigen/src/Core/Map.h
@@ -47,7 +47,8 @@ struct ei_traits<Map<MatrixType> >
     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = NumTraits<Scalar>::ReadCost
   };
 };
 
diff --git a/Eigen/src/Core/Matrix.h b/Eigen/src/Core/Matrix.h
index 9ccad4e7e..fa40364c7 100644
--- a/Eigen/src/Core/Matrix.h
+++ b/Eigen/src/Core/Matrix.h
@@ -79,7 +79,8 @@ struct ei_traits<Matrix<_Scalar, _Rows, _Cols, _Flags, _MaxRows, _MaxCols> >
     ColsAtCompileTime = _Cols,
     MaxRowsAtCompileTime = _MaxRows,
     MaxColsAtCompileTime = _MaxCols,
-    Flags = _Flags
+    Flags = _Flags,
+    CoeffReadCost = NumTraits<Scalar>::ReadCost
   };
 };
 
diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h
index 7f8283e9f..77ecb2906 100644
--- a/Eigen/src/Core/MatrixBase.h
+++ b/Eigen/src/Core/MatrixBase.h
@@ -126,7 +126,7 @@ template<typename Derived> class MatrixBase
           * we are dealing with a column-vector (if there is only one column) or with
           * a row-vector (if there is only one row). */
 
-      Flags = ei_traits<Derived>::Flags
+      Flags = ei_traits<Derived>::Flags,
         /**< This stores expression metadata which typically is inherited by new expressions
           * constructed from this one. The available flags are:
           * \li \c RowMajorBit: if this bit is set, the preferred storage order for an evaluation
@@ -136,6 +136,8 @@ template<typename Derived> class MatrixBase
           * \li \c LargeBit: if this bit is set, optimization will be tuned for large matrices (typically,
           *               at least 32x32).
           */
+
+      CoeffReadCost = ei_traits<Derived>::CoeffReadCost
     };
 
     /** This is the "real scalar" type; if the \a Scalar type is already real numbers
@@ -219,14 +221,14 @@ template<typename Derived> class MatrixBase
       * sum, scalar multiple, ...
       */
     //@{
-    const CwiseUnaryOp<ei_scalar_opposite_op,Derived> operator-() const;
+    const CwiseUnaryOp<ei_scalar_opposite_op<typename ei_traits<Derived>::Scalar>,Derived> operator-() const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_sum_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_sum_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     operator+(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_difference_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_difference_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     operator-(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
@@ -237,10 +239,10 @@ template<typename Derived> class MatrixBase
     Derived& operator*=(const Scalar& other);
     Derived& operator/=(const Scalar& other);
 
-    const CwiseUnaryOp<ei_scalar_multiple_op<Scalar>, Derived> operator*(const Scalar& scalar) const;
-    const CwiseUnaryOp<ei_scalar_quotient1_op<Scalar>, Derived> operator/(const Scalar& scalar) const;
+    const CwiseUnaryOp<ei_scalar_multiple_op<typename ei_traits<Derived>::Scalar>, Derived> operator*(const Scalar& scalar) const;
+    const CwiseUnaryOp<ei_scalar_quotient1_op<typename ei_traits<Derived>::Scalar>, Derived> operator/(const Scalar& scalar) const;
 
-    friend const CwiseUnaryOp<ei_scalar_multiple_op<Scalar>, Derived>
+    friend const CwiseUnaryOp<ei_scalar_multiple_op<typename ei_traits<Derived>::Scalar>, Derived>
     operator*(const Scalar& scalar, const MatrixBase& matrix)
     { return matrix*scalar; }
     //@}
@@ -250,7 +252,7 @@ template<typename Derived> class MatrixBase
       */
     //@{
     template<typename OtherDerived>
-    const typename ei_eval_unless_lazy<Product<Derived, OtherDerived> >::Type
+    const typename ei_eval_unless_lazy<Product<Derived, OtherDerived> >::type
     operator*(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
@@ -265,11 +267,11 @@ template<typename Derived> class MatrixBase
     Scalar dot(const MatrixBase<OtherDerived>& other) const;
     RealScalar norm2() const;
     RealScalar norm()  const;
-    const CwiseUnaryOp<ei_scalar_multiple_op<Scalar>, Derived> normalized() const;
+    const CwiseUnaryOp<ei_scalar_multiple_op<typename ei_traits<Derived>::Scalar>, Derived> normalized() const;
 
     Transpose<Derived> transpose();
     const Transpose<Derived> transpose() const;
-    const Transpose<CwiseUnaryOp<ei_scalar_conjugate_op, Derived> > adjoint() const;
+    const Transpose<CwiseUnaryOp<ei_scalar_conjugate_op<typename ei_traits<Derived>::Scalar>, Derived> > adjoint() const;
     //@}
 
     /// \name Sub-matrices
@@ -310,9 +312,9 @@ template<typename Derived> class MatrixBase
 
     /// \name Generating special matrices
     //@{
-    static const typename ei_eval_unless_lazy<Random<Derived> >::Type random(int rows, int cols);
-    static const typename ei_eval_unless_lazy<Random<Derived> >::Type random(int size);
-    static const typename ei_eval_unless_lazy<Random<Derived> >::Type random();
+    static const typename ei_eval_unless_lazy<Random<Derived> >::type random(int rows, int cols);
+    static const typename ei_eval_unless_lazy<Random<Derived> >::type random(int size);
+    static const typename ei_eval_unless_lazy<Random<Derived> >::type random();
     static const Zero<Derived> zero(int rows, int cols);
     static const Zero<Derived> zero(int size);
     static const Zero<Derived> zero();
@@ -354,11 +356,11 @@ template<typename Derived> class MatrixBase
     /// \name Special functions
     //@{
     template<typename NewType>
-    const CwiseUnaryOp<ei_scalar_cast_op<NewType>, Derived> cast() const;
+    const CwiseUnaryOp<ei_scalar_cast_op<typename ei_traits<Derived>::Scalar, NewType>, Derived> cast() const;
 
-    const typename ei_eval_unless_lazy<Derived>::Type eval() const EIGEN_ALWAYS_INLINE
+    const typename ei_eval_unless_lazy<Derived>::type eval() const EIGEN_ALWAYS_INLINE
     {
-      return typename ei_eval_unless_lazy<Derived>::Type(derived());
+      return typename ei_eval_unless_lazy<Derived>::type(derived());
     }
 
     template<typename OtherDerived>
@@ -369,31 +371,31 @@ template<typename Derived> class MatrixBase
 
     /// \name Coefficient-wise operations
     //@{
-    const CwiseUnaryOp<ei_scalar_conjugate_op, Derived> conjugate() const;
+    const CwiseUnaryOp<ei_scalar_conjugate_op<typename ei_traits<Derived>::Scalar>, Derived> conjugate() const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_product_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_product_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     cwiseProduct(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_quotient_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_quotient_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     cwiseQuotient(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_min_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_min_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     cwiseMin(const MatrixBase<OtherDerived> &other) const;
 
     template<typename OtherDerived>
-    const CwiseBinaryOp<ei_scalar_max_op, Derived, OtherDerived>
+    const CwiseBinaryOp<ei_scalar_max_op<typename ei_traits<Derived>::Scalar>, Derived, OtherDerived>
     cwiseMax(const MatrixBase<OtherDerived> &other) const;
 
-    const CwiseUnaryOp<ei_scalar_abs_op, Derived> cwiseAbs() const;
-    const CwiseUnaryOp<ei_scalar_abs2_op, Derived> cwiseAbs2() const;
-    const CwiseUnaryOp<ei_scalar_sqrt_op, Derived> cwiseSqrt() const;
-    const CwiseUnaryOp<ei_scalar_exp_op, Derived> cwiseExp() const;
-    const CwiseUnaryOp<ei_scalar_log_op, Derived> cwiseLog() const;
-    const CwiseUnaryOp<ei_scalar_cos_op, Derived> cwiseCos() const;
-    const CwiseUnaryOp<ei_scalar_sin_op, Derived> cwiseSin() const;
+    const CwiseUnaryOp<ei_scalar_abs_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseAbs() const;
+    const CwiseUnaryOp<ei_scalar_abs2_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseAbs2() const;
+    const CwiseUnaryOp<ei_scalar_sqrt_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseSqrt() const;
+    const CwiseUnaryOp<ei_scalar_exp_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseExp() const;
+    const CwiseUnaryOp<ei_scalar_log_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseLog() const;
+    const CwiseUnaryOp<ei_scalar_cos_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseCos() const;
+    const CwiseUnaryOp<ei_scalar_sin_op<typename ei_traits<Derived>::Scalar>, Derived> cwiseSin() const;
     const CwiseUnaryOp<ei_scalar_pow_op<typename ei_traits<Derived>::Scalar>, Derived>
     cwisePow(const Scalar& exponent) const;
 
diff --git a/Eigen/src/Core/Minor.h b/Eigen/src/Core/Minor.h
index 62941d764..911ac2151 100644
--- a/Eigen/src/Core/Minor.h
+++ b/Eigen/src/Core/Minor.h
@@ -50,7 +50,8 @@ struct ei_traits<Minor<MatrixType> >
                                 MatrixType::MaxRowsAtCompileTime - 1 : Dynamic,
     MaxColsAtCompileTime = (MatrixType::MaxColsAtCompileTime != Dynamic) ?
                                 MatrixType::MaxColsAtCompileTime - 1 : Dynamic,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = MatrixType::CoeffReadCost
   };
 };
 
diff --git a/Eigen/src/Core/NumTraits.h b/Eigen/src/Core/NumTraits.h
index cc541af8f..cb601b160 100644
--- a/Eigen/src/Core/NumTraits.h
+++ b/Eigen/src/Core/NumTraits.h
@@ -53,7 +53,10 @@ template<> struct NumTraits<int>
   typedef double FloatingPoint;
   enum {
     IsComplex = 0,
-    HasFloatingPoint = 0
+    HasFloatingPoint = 0,
+    ReadCost = 1,
+    AddCost = 1,
+    MulCost = 3
   };
 };
 
@@ -63,7 +66,10 @@ template<> struct NumTraits<float>
   typedef float FloatingPoint;
   enum {
     IsComplex = 0,
-    HasFloatingPoint = 1
+    HasFloatingPoint = 1,
+    ReadCost = 1,
+    AddCost = 2,
+    MulCost = 6
   };
 };
 
@@ -73,7 +79,10 @@ template<> struct NumTraits<double>
   typedef double FloatingPoint;
   enum {
     IsComplex = 0,
-    HasFloatingPoint = 1
+    HasFloatingPoint = 1,
+    ReadCost = 1,
+    AddCost = 2,
+    MulCost = 6
   };
 };
 
@@ -83,7 +92,10 @@ template<typename _Real> struct NumTraits<std::complex<_Real> >
   typedef std::complex<_Real> FloatingPoint;
   enum {
     IsComplex = 1,
-    HasFloatingPoint = 1 // anyway we don't allow std::complex<int>
+    HasFloatingPoint = NumTraits<Real>::HasFloatingPoint,
+    ReadCost = 2,
+    AddCost = 2 * NumTraits<Real>::AddCost,
+    MulCost = 4 * NumTraits<Real>::MulCost + 2 * NumTraits<Real>::AddCost
   };
 };
 
@@ -93,7 +105,10 @@ template<> struct NumTraits<long long int>
   typedef long double FloatingPoint;
   enum {
     IsComplex = 0,
-    HasFloatingPoint = 0
+    HasFloatingPoint = 0,
+    ReadCost = 1,
+    AddCost = 2,
+    MulCost = 6
   };
 };
 
@@ -103,7 +118,10 @@ template<> struct NumTraits<long double>
   typedef long double FloatingPoint;
   enum {
     IsComplex = 0,
-    HasFloatingPoint = 1
+    HasFloatingPoint = 1,
+    ReadCost = 1,
+    AddCost = 2,
+    MulCost = 6
   };
 };
 
diff --git a/Eigen/src/Core/Ones.h b/Eigen/src/Core/Ones.h
index 05d9b4270..4cb4bc348 100644
--- a/Eigen/src/Core/Ones.h
+++ b/Eigen/src/Core/Ones.h
@@ -41,7 +41,8 @@ struct ei_traits<Ones<MatrixType> >
     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = NumTraits<Scalar>::ReadCost
   };
 };
 
diff --git a/Eigen/src/Core/Product.h b/Eigen/src/Core/Product.h
index 62e1b6317..c71637779 100644
--- a/Eigen/src/Core/Product.h
+++ b/Eigen/src/Core/Product.h
@@ -78,6 +78,23 @@ template<typename Lhs, typename Rhs, int EvalMode>
 struct ei_traits<Product<Lhs, Rhs, EvalMode> >
 {
   typedef typename Lhs::Scalar Scalar;
+  typedef typename ei_meta_if<
+                     (int)NumTraits<Scalar>::ReadCost < (int)Lhs::CoeffReadCost,
+                     typename Lhs::Eval,
+                     Lhs>::ret ActualLhs;
+  typedef typename ei_meta_if<
+                     (int)NumTraits<Scalar>::ReadCost < (int)Lhs::CoeffReadCost,
+                     typename Lhs::Eval,
+                     typename Lhs::XprCopy>::ret ActualLhsXprCopy;
+
+  typedef typename ei_meta_if<
+                     (int)NumTraits<Scalar>::ReadCost < (int)Rhs::CoeffReadCost,
+                     typename Rhs::Eval,
+                     Rhs>::ret ActualRhs;
+  typedef typename ei_meta_if<
+                     (int)NumTraits<Scalar>::ReadCost < (int)Rhs::CoeffReadCost,
+                     typename Rhs::Eval,
+                     typename Rhs::XprCopy>::ret ActualRhsXprCopy;
   enum {
     RowsAtCompileTime = Lhs::RowsAtCompileTime,
     ColsAtCompileTime = Rhs::ColsAtCompileTime,
@@ -85,20 +102,20 @@ struct ei_traits<Product<Lhs, Rhs, EvalMode> >
     MaxColsAtCompileTime = Rhs::MaxColsAtCompileTime,
     Flags = (RowsAtCompileTime == Dynamic || ColsAtCompileTime == Dynamic)
           ? (unsigned int)(Lhs::Flags | Rhs::Flags)
-          : (unsigned int)(Lhs::Flags | Rhs::Flags) & ~LargeBit
+          : (unsigned int)(Lhs::Flags | Rhs::Flags) & ~LargeBit,
+    CoeffReadCost
+      = Lhs::ColsAtCompileTime == Dynamic
+      ? Dynamic
+      : Lhs::ColsAtCompileTime
+        * (NumTraits<Scalar>::MulCost + ActualLhs::CoeffReadCost + ActualRhs::CoeffReadCost)
+        + (Lhs::ColsAtCompileTime - 1) * NumTraits<Scalar>::AddCost
   };
 };
 
-template<typename Lhs, typename Rhs>
-struct ei_product_eval_mode
+template<typename Lhs, typename Rhs> struct ei_product_eval_mode
 {
-  enum {
-    SizeAtCompileTime = MatrixBase<Product<Lhs,Rhs,UnrolledDotProduct> >::SizeAtCompileTime,
-    MaxSizeAtCompileTime = MatrixBase<Product<Lhs,Rhs,UnrolledDotProduct> >::MaxSizeAtCompileTime,
-    EvalMode = ( EIGEN_UNROLLED_LOOPS
-              && MaxSizeAtCompileTime != Dynamic
-              && SizeAtCompileTime <= EIGEN_UNROLLING_LIMIT) ? UnrolledDotProduct : CacheOptimal,
-  };
+  enum{ value = Lhs::MaxRowsAtCompileTime == Dynamic || Rhs::MaxColsAtCompileTime == Dynamic
+                  ? CacheOptimal : UnrolledDotProduct };
 };
 
 template<typename Lhs, typename Rhs, int EvalMode> class Product : ei_no_assignment_operator,
@@ -107,6 +124,10 @@ template<typename Lhs, typename Rhs, int EvalMode> class Product : ei_no_assignm
   public:
 
     EIGEN_GENERIC_PUBLIC_INTERFACE(Product)
+    typedef typename ei_traits<Product>::ActualLhs ActualLhs;
+    typedef typename ei_traits<Product>::ActualRhs ActualRhs;
+    typedef typename ei_traits<Product>::ActualLhsXprCopy ActualLhsXprCopy;
+    typedef typename ei_traits<Product>::ActualRhsXprCopy ActualRhsXprCopy;
 
     Product(const Lhs& lhs, const Rhs& rhs)
       : m_lhs(lhs), m_rhs(rhs)
@@ -132,7 +153,7 @@ template<typename Lhs, typename Rhs, int EvalMode> class Product : ei_no_assignm
         ei_product_unroller<Lhs::ColsAtCompileTime-1,
                             Lhs::ColsAtCompileTime <= EIGEN_UNROLLING_LIMIT
                               ? Lhs::ColsAtCompileTime : Dynamic,
-                            Lhs, Rhs>
+                            ActualLhs, ActualRhs>
           ::run(row, col, m_lhs, m_rhs, res);
       else
       {
@@ -144,8 +165,8 @@ template<typename Lhs, typename Rhs, int EvalMode> class Product : ei_no_assignm
     }
 
   protected:
-    const typename Lhs::XprCopy m_lhs;
-    const typename Rhs::XprCopy m_rhs;
+    const ActualLhsXprCopy m_lhs;
+    const ActualRhsXprCopy m_rhs;
 };
 
 /** \returns the matrix product of \c *this and \a other.
@@ -157,7 +178,7 @@ template<typename Lhs, typename Rhs, int EvalMode> class Product : ei_no_assignm
   */
 template<typename Derived>
 template<typename OtherDerived>
-const typename ei_eval_unless_lazy<Product<Derived, OtherDerived> >::Type
+const typename ei_eval_unless_lazy<Product<Derived, OtherDerived> >::type
 MatrixBase<Derived>::operator*(const MatrixBase<OtherDerived> &other) const
 {
   return Product<Derived, OtherDerived>(derived(), other.derived()).eval();
@@ -199,7 +220,8 @@ void Product<Lhs,Rhs,EvalMode>::_cacheOptimalEval(DestDerived& res) const
       const Scalar tmp2 = m_rhs.coeff(j+2,k);
       const Scalar tmp3 = m_rhs.coeff(j+3,k);
       for (int i=0; i<m_lhs.rows(); ++i)
-        res.coeffRef(i,k) += tmp0 * m_lhs.coeff(i,j) + tmp1 * m_lhs.coeff(i,j+1) + tmp2 * m_lhs.coeff(i,j+2) + tmp3 * m_lhs.coeff(i,j+3);
+        res.coeffRef(i,k) += tmp0 * m_lhs.coeff(i,j) + tmp1 * m_lhs.coeff(i,j+1)
+                             + tmp2 * m_lhs.coeff(i,j+2) + tmp3 * m_lhs.coeff(i,j+3);
     }
     for (; j<m_lhs.cols(); ++j)
     {
diff --git a/Eigen/src/Core/Random.h b/Eigen/src/Core/Random.h
index 03f08ef06..4bd482958 100644
--- a/Eigen/src/Core/Random.h
+++ b/Eigen/src/Core/Random.h
@@ -41,7 +41,8 @@ struct ei_traits<Random<MatrixType> >
     ColsAtCompileTime = ei_traits<MatrixType>::ColsAtCompileTime,
     MaxRowsAtCompileTime = ei_traits<MatrixType>::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = ei_traits<MatrixType>::MaxColsAtCompileTime,
-    Flags = ei_traits<MatrixType>::Flags
+    Flags = ei_traits<MatrixType>::Flags,
+    CoeffReadCost = 2 * NumTraits<Scalar>::MulCost // FIXME: arbitrary value
   };
 };
 
@@ -92,7 +93,7 @@ template<typename MatrixType> class Random : ei_no_assignment_operator,
   * \sa ei_random(), ei_random(int)
   */
 template<typename Derived>
-const typename ei_eval_unless_lazy<Random<Derived> >::Type
+const typename ei_eval_unless_lazy<Random<Derived> >::type
 MatrixBase<Derived>::random(int rows, int cols)
 {
   return Random<Derived>(rows, cols).eval();
@@ -115,7 +116,7 @@ MatrixBase<Derived>::random(int rows, int cols)
   * \sa ei_random(), ei_random(int,int)
   */
 template<typename Derived>
-const typename ei_eval_unless_lazy<Random<Derived> >::Type
+const typename ei_eval_unless_lazy<Random<Derived> >::type
 MatrixBase<Derived>::random(int size)
 {
   ei_assert(IsVectorAtCompileTime);
@@ -135,7 +136,7 @@ MatrixBase<Derived>::random(int size)
   * \sa ei_random(int), ei_random(int,int)
   */
 template<typename Derived>
-const typename ei_eval_unless_lazy<Random<Derived> >::Type
+const typename ei_eval_unless_lazy<Random<Derived> >::type
 MatrixBase<Derived>::random()
 {
   return Random<Derived>(RowsAtCompileTime, ColsAtCompileTime).eval();
diff --git a/Eigen/src/Core/Redux.h b/Eigen/src/Core/Redux.h
index 34eefca7c..05d73afc2 100644
--- a/Eigen/src/Core/Redux.h
+++ b/Eigen/src/Core/Redux.h
@@ -26,7 +26,6 @@
 #ifndef EIGEN_REDUX_H
 #define EIGEN_REDUX_H
 
-
 template<typename BinaryOp, typename Derived, int Start, int Length>
 struct ei_redux_unroller
 {
@@ -95,7 +94,8 @@ struct ei_traits<PartialRedux<Direction, BinaryOp, MatrixType> >
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
     Flags = (RowsAtCompileTime == Dynamic || ColsAtCompileTime == Dynamic)
           ? (unsigned int)MatrixType::Flags
-          : (unsigned int)MatrixType::Flags & ~LargeBit
+          : (unsigned int)MatrixType::Flags & ~LargeBit,
+    CoeffReadCost = 1 //FIXME -- unimplemented!
   };
 };
 
@@ -198,7 +198,7 @@ template<typename Derived>
 typename ei_traits<Derived>::Scalar
 MatrixBase<Derived>::sum() const
 {
-  return this->redux(Eigen::ei_scalar_sum_op());
+  return this->redux(Eigen::ei_scalar_sum_op<Scalar>());
 }
 
 /** \returns the trace of \c *this, i.e. the sum of the coefficients on the main diagonal.
@@ -220,7 +220,7 @@ template<typename Derived>
 typename ei_traits<Derived>::Scalar
 MatrixBase<Derived>::minCoeff() const
 {
-  return this->redux(Eigen::ei_scalar_min_op());
+  return this->redux(Eigen::ei_scalar_min_op<Scalar>());
 }
 
 /** \returns the maximum of all coefficients of *this
@@ -229,7 +229,7 @@ template<typename Derived>
 typename ei_traits<Derived>::Scalar
 MatrixBase<Derived>::maxCoeff() const
 {
-  return this->redux(Eigen::ei_scalar_max_op());
+  return this->redux(Eigen::ei_scalar_max_op<Scalar>());
 }
 
 #endif // EIGEN_REDUX_H
diff --git a/Eigen/src/Core/Swap.h b/Eigen/src/Core/Swap.h
index db096e6f2..cffdb1cc6 100644
--- a/Eigen/src/Core/Swap.h
+++ b/Eigen/src/Core/Swap.h
@@ -63,7 +63,7 @@ void MatrixBase<Derived>::swap(const MatrixBase<OtherDerived>& other)
   }
   else // SizeAtCompileTime != Dynamic
   {
-    typename Eval<Derived>::MatrixType buf(*this);
+    typename Derived::Eval buf(*this);
     *this = other;
     *_other = buf;
   }
diff --git a/Eigen/src/Core/Transpose.h b/Eigen/src/Core/Transpose.h
index 51799b33e..3eb5d1f71 100644
--- a/Eigen/src/Core/Transpose.h
+++ b/Eigen/src/Core/Transpose.h
@@ -46,7 +46,8 @@ struct ei_traits<Transpose<MatrixType> >
     ColsAtCompileTime = MatrixType::RowsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxColsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
-    Flags = MatrixType::Flags ^ RowMajorBit
+    Flags = MatrixType::Flags ^ RowMajorBit,
+    CoeffReadCost = MatrixType::CoeffReadCost
   };
 };
 
@@ -108,7 +109,7 @@ MatrixBase<Derived>::transpose() const
   *
   * \sa transpose(), conjugate(), class Transpose, class ei_scalar_conjugate_op */
 template<typename Derived>
-const Transpose<CwiseUnaryOp<ei_scalar_conjugate_op, Derived> >
+const Transpose<CwiseUnaryOp<ei_scalar_conjugate_op<typename ei_traits<Derived>::Scalar>, Derived> >
 MatrixBase<Derived>::adjoint() const
 {
   return conjugate().transpose();
diff --git a/Eigen/src/Core/Util.h b/Eigen/src/Core/Util.h
index 9aaac9ae2..d82ec574f 100644
--- a/Eigen/src/Core/Util.h
+++ b/Eigen/src/Core/Util.h
@@ -109,7 +109,8 @@ EIGEN_INHERIT_SCALAR_ASSIGNMENT_OPERATOR(Derived, /=)
 #define _EIGEN_GENERIC_PUBLIC_INTERFACE(Derived, BaseClass) \
 typedef BaseClass Base; \
 typedef typename Eigen::ei_traits<Derived>::Scalar Scalar; \
-typedef typename Eigen::ei_xpr_copy<Derived>::Type XprCopy; \
+typedef typename Eigen::ei_xpr_copy<Derived>::type XprCopy; \
+typedef typename Eigen::ei_eval<Derived>::type Eval; \
 enum { RowsAtCompileTime = Base::RowsAtCompileTime, \
        ColsAtCompileTime = Base::ColsAtCompileTime, \
        MaxRowsAtCompileTime = Base::MaxRowsAtCompileTime, \
@@ -117,13 +118,8 @@ enum { RowsAtCompileTime = Base::RowsAtCompileTime, \
        SizeAtCompileTime = Base::SizeAtCompileTime, \
        MaxSizeAtCompileTime = Base::MaxSizeAtCompileTime, \
        IsVectorAtCompileTime = Base::IsVectorAtCompileTime, \
-       Flags = Base::Flags }; \
-typedef Matrix<Scalar, \
-               RowsAtCompileTime, \
-               ColsAtCompileTime, \
-               Flags, \
-               MaxRowsAtCompileTime, \
-               MaxColsAtCompileTime> Eval;
+       Flags = Base::Flags, \
+       CoeffReadCost = Base::CoeffReadCost };
 
 #define EIGEN_GENERIC_PUBLIC_INTERFACE(Derived) \
 _EIGEN_GENERIC_PUBLIC_INTERFACE(Derived, Eigen::MatrixBase<Derived>) \
@@ -138,6 +134,8 @@ const unsigned int RowMajorBit = 0x1;
 const unsigned int LazyBit = 0x2;
 const unsigned int LargeBit = 0x4;
 
+enum { ConditionalJumpCost = 5 };
+
 enum CornerType { TopLeft, TopRight, BottomLeft, BottomRight };
 
 enum DirectionType { Vertical, Horizontal };
@@ -185,6 +183,9 @@ struct ei_meta_if { typedef Then ret; };
 template <class Then, class Else>
 struct ei_meta_if <false, Then, Else> { typedef Else ret; };
 
+template<typename T, typename U> struct ei_is_same_type { enum { ret = 0 }; };
+template<typename T> struct ei_is_same_type<T,T> { enum { ret = 1 }; };
+
 
 /** \internal
   * Convenient struct to get the result type of a unary or binary functor.
diff --git a/Eigen/src/Core/Zero.h b/Eigen/src/Core/Zero.h
index 106742d32..15108b794 100644
--- a/Eigen/src/Core/Zero.h
+++ b/Eigen/src/Core/Zero.h
@@ -41,7 +41,8 @@ struct ei_traits<Zero<MatrixType> >
     ColsAtCompileTime = MatrixType::ColsAtCompileTime,
     MaxRowsAtCompileTime = MatrixType::MaxRowsAtCompileTime,
     MaxColsAtCompileTime = MatrixType::MaxColsAtCompileTime,
-    Flags = MatrixType::Flags
+    Flags = MatrixType::Flags,
+    CoeffReadCost = NumTraits<Scalar>::ReadCost
   };
 };
 
diff --git a/disabled/Eval.h b/disabled/Eval.h
index 23e35a96a..074a12df9 100644
--- a/disabled/Eval.h
+++ b/disabled/Eval.h
@@ -101,9 +101,9 @@ template<typename ExpressionType> class Eval : ei_no_assignment_operator,
   *
   * \sa class Eval */
 template<typename Derived>
-const typename ei_eval_unless_lazy<Derived>::Type MatrixBase<Derived>::eval() const
+const typename ei_eval_unless_lazy<Derived>::type MatrixBase<Derived>::eval() const
 {
-  return typename ei_eval_unless_lazy<Derived>::Type(derived());
+  return typename ei_eval_unless_lazy<Derived>::type(derived());
 }
 
 #endif // EIGEN_EVAL_H
diff --git a/doc/echelon.cpp b/doc/echelon.cpp
index d28726b0f..e305eb238 100644
--- a/doc/echelon.cpp
+++ b/doc/echelon.cpp
@@ -16,13 +16,8 @@ void echelon(MatrixBase<Derived>& m)
      .maxCoeff(&rowOfBiggest, &colOfBiggest);
     m.row(k).swap(m.row(k+rowOfBiggest));
     m.col(k).swap(m.col(k+colOfBiggest));
-    // important performance tip:
-    // in a complex expression such as below it can be very important to fine-tune
-    // exactly where evaluation occurs. The parentheses and .eval() below ensure
-    // that the quotient is computed only once, and that the evaluation caused
-    // by operator* occurs last.
     m.corner(BottomRight, cornerRows-1, cornerCols)
-      -= m.col(k).end(cornerRows-1) * (m.row(k).end(cornerCols) / m(k,k)).eval();
+      -= m.col(k).end(cornerRows-1) * (m.row(k).end(cornerCols) / m(k,k));
   }
 }
 
diff --git a/doc/examples/class_CwiseBinaryOp.cpp b/doc/examples/class_CwiseBinaryOp.cpp
index 8e4c5361e..b91abd1ef 100644
--- a/doc/examples/class_CwiseBinaryOp.cpp
+++ b/doc/examples/class_CwiseBinaryOp.cpp
@@ -1,26 +1,29 @@
+// FIXME - this example is not too good as that functionality is provided in the Eigen API
+// additionally it's quite heavy. the CwiseUnaryOp example is better.
+
 #include <Eigen/Core>
 USING_PART_OF_NAMESPACE_EIGEN
 using namespace std;
 
 // define a custom template binary functor
-struct CwiseMinOp EIGEN_EMPTY_STRUCT {
-    template<typename Scalar>
+template<typename Scalar> struct CwiseMinOp EIGEN_EMPTY_STRUCT {
     Scalar operator()(const Scalar& a, const Scalar& b) const { return std::min(a,b); }
+    enum { Cost = Eigen::ConditionalJumpCost + Eigen::NumTraits<Scalar>::AddCost };
 };
 
 // define a custom binary operator between two matrices
 template<typename Derived1, typename Derived2>
-const Eigen::CwiseBinaryOp<CwiseMinOp, Derived1, Derived2>
+const Eigen::CwiseBinaryOp<CwiseMinOp<typename Derived1::Scalar>, Derived1, Derived2>
 cwiseMin(const MatrixBase<Derived1> &mat1, const MatrixBase<Derived2> &mat2)
 {
-  return Eigen::CwiseBinaryOp<CwiseMinOp, Derived1, Derived2>(mat1, mat2);
+  return Eigen::CwiseBinaryOp<CwiseMinOp<typename Derived1::Scalar>, Derived1, Derived2>(mat1, mat2);
 }
 
 int main(int, char**)
 {
   Matrix4d m1 = Matrix4d::random(), m2 = Matrix4d::random();
   cout << cwiseMin(m1,m2) << endl;            // use our new global operator
-  cout << m1.cwise<CwiseMinOp>(m2) << endl;   // directly use the generic expression member
-  cout << m1.cwise(m2, CwiseMinOp()) << endl; // directly use the generic expression member (variant)
+  cout << m1.cwise<CwiseMinOp<double> >(m2) << endl;   // directly use the generic expression member
+  cout << m1.cwise(m2, CwiseMinOp<double>()) << endl; // directly use the generic expression member (variant)
   return 0;
 }
diff --git a/doc/examples/class_CwiseUnaryOp.cpp b/doc/examples/class_CwiseUnaryOp.cpp
index 042e051e9..ea27876eb 100644
--- a/doc/examples/class_CwiseUnaryOp.cpp
+++ b/doc/examples/class_CwiseUnaryOp.cpp
@@ -5,9 +5,10 @@ using namespace std;
 // define a custom template binary functor
 template<typename Scalar>
 struct CwiseClampOp EIGEN_EMPTY_STRUCT {
-    CwiseClampOp(const Scalar& inf, const Scalar& sup) : m_inf(inf), m_sup(sup) {}
-    Scalar operator()(const Scalar& x) const { return x<m_inf ? m_inf : (x>m_sup ? m_sup : x); }
-    Scalar m_inf, m_sup;
+  CwiseClampOp(const Scalar& inf, const Scalar& sup) : m_inf(inf), m_sup(sup) {}
+  const Scalar operator()(const Scalar& x) const { return x<m_inf ? m_inf : (x>m_sup ? m_sup : x); }
+  Scalar m_inf, m_sup;
+  enum { Cost = Eigen::ConditionalJumpCost + Eigen::NumTraits<Scalar>::AddCost };
 };
 
 int main(int, char**)
diff --git a/test/cwiseop.cpp b/test/cwiseop.cpp
index 5ca51b74e..ee9f068a4 100644
--- a/test/cwiseop.cpp
+++ b/test/cwiseop.cpp
@@ -30,8 +30,9 @@
 
 namespace Eigen {
 
-struct AddIfNull {
-    template<typename Scalar> Scalar operator() (const Scalar a, const Scalar b) const {return a<=1e-3 ? b : a;}
+template<typename Scalar> struct AddIfNull {
+    const Scalar operator() (const Scalar a, const Scalar b) const {return a<=1e-3 ? b : a;}
+    enum { Cost = NumTraits<Scalar>::AddCost };
 };
 
 template<typename MatrixType> void cwiseops(const MatrixType& m)
@@ -55,7 +56,7 @@ template<typename MatrixType> void cwiseops(const MatrixType& m)
              v2 = VectorType::random(rows),
              vzero = VectorType::zero(rows);
 
-  m2 = m2.template cwise<AddIfNull>(mones);
+  m2 = m2.template cwise<AddIfNull<Scalar> >(mones);
   
   VERIFY_IS_APPROX(               mzero,    m1-m1);
   VERIFY_IS_APPROX(               m2,       m1+m2-m1);