add overloads of lazyAssign to detect common aliasing issue with

transpose and adjoint

Eigen/Core

@@ -150,6 +150,7 @@ namespace Eigen {
   #include "src/Core/Assign.h"
 #endif
 
+#include "src/Core/util/BlasUtil.h"
 #include "src/Core/MatrixStorage.h"
 #include "src/Core/NestByValue.h"
 #include "src/Core/ReturnByValue.h"
@@ -178,7 +179,6 @@ namespace Eigen {
 #include "src/Core/IO.h"
 #include "src/Core/Swap.h"
 #include "src/Core/CommaInitializer.h"
-#include "src/Core/util/BlasUtil.h"
 #include "src/Core/ProductBase.h"
 #include "src/Core/Product.h"
 #include "src/Core/TriangularMatrix.h"

Eigen/src/Core/CwiseBinaryOp.h

@@ -85,6 +85,8 @@ class CwiseBinaryOp : ei_no_assignment_operator,
     EIGEN_GENERIC_PUBLIC_INTERFACE(CwiseBinaryOp)
     typedef typename ei_traits<CwiseBinaryOp>::LhsNested LhsNested;
     typedef typename ei_traits<CwiseBinaryOp>::RhsNested RhsNested;
+    typedef typename ei_traits<CwiseBinaryOp>::_LhsNested _LhsNested;
+    typedef typename ei_traits<CwiseBinaryOp>::_RhsNested _RhsNested;
 
     EIGEN_STRONG_INLINE CwiseBinaryOp(const Lhs& lhs, const Rhs& rhs, const BinaryOp& func = BinaryOp())
       : m_lhs(lhs), m_rhs(rhs), m_functor(func)
@@ -130,6 +132,10 @@ class CwiseBinaryOp : ei_no_assignment_operator,
       return m_functor.packetOp(m_lhs.template packet<LoadMode>(index), m_rhs.template packet<LoadMode>(index));
     }
 
+    const _LhsNested& lhs() const { return m_lhs; }
+    const _RhsNested& rhs() const { return m_rhs; }
+    const BinaryOp& functor() const { return m_functor; }
+
   protected:
     const LhsNested m_lhs;
     const RhsNested m_rhs;

Eigen/src/Core/MatrixBase.h

@@ -457,6 +457,21 @@ template<typename Derived> class MatrixBase
     void transposeInPlace();
     const AdjointReturnType adjoint() const;
     void adjointInPlace();
+    #ifndef EIGEN_NO_DEBUG
+    template<typename OtherDerived>
+    Derived& lazyAssign(const Transpose<OtherDerived>& other);
+    template<typename DerivedA, typename DerivedB>
+    Derived& lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,Transpose<DerivedA>,DerivedB>& other);
+    template<typename DerivedA, typename DerivedB>
+    Derived& lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,Transpose<DerivedB> >& other);
+
+    template<typename OtherDerived>
+    Derived& lazyAssign(const CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<OtherDerived> > >& other);
+    template<typename DerivedA, typename DerivedB>
+    Derived& lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedA> > >,DerivedB>& other);
+    template<typename DerivedA, typename DerivedB>
+    Derived& lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedB> > > >& other);
+    #endif
 
     RowXpr row(int i);
     const RowXpr row(int i) const;

Eigen/src/Core/Transpose.h

@@ -267,4 +267,92 @@ inline void MatrixBase<Derived>::adjointInPlace()
   derived() = adjoint().eval();
 }
 
+#ifndef EIGEN_NO_DEBUG
+
+// The following is to detect aliasing problems in the following common cases:
+// a = a.transpose()
+// a = a.transpose() + X
+// a = X + a.transpose()
+// a = a.adjoint()
+// a = a.adjoint() + X
+// a = X + a.adjoint()
+
+template<typename T, int Access=ei_blas_traits<T>::ActualAccess>
+struct ei_extract_data_selector {
+  static typename T::Scalar* run(const T& m)
+  {
+    return &ei_blas_traits<T>::extract(m).const_cast_derived().coeffRef(0,0);
+  }
+};
+
+template<typename T>
+struct ei_extract_data_selector<T,NoDirectAccess> {
+  static typename T::Scalar* run(const T&) { return 0; }
+};
+
+template<typename T> typename T::Scalar* ei_extract_data(const T& m)
+{
+  return ei_extract_data_selector<T>::run(m);
+}
+
+template<typename Derived>
+template<typename OtherDerived>
+Derived& MatrixBase<Derived>::lazyAssign(const Transpose<OtherDerived>& other)
+{
+  ei_assert(ei_extract_data(other) != ei_extract_data(derived())
+            && "aliasing detected during tranposition, please use transposeInPlace()");
+  return lazyAssign(static_cast<const MatrixBase<Transpose<OtherDerived> >& >(other));
+}
+
+template<typename Derived>
+template<typename DerivedA, typename DerivedB>
+Derived& MatrixBase<Derived>::
+lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,Transpose<DerivedA>,DerivedB>& other)
+{
+  ei_assert(ei_extract_data(derived()) != ei_extract_data(other.lhs())
+            && "aliasing detected during tranposition, please evaluate your expression");
+  return lazyAssign(static_cast<const MatrixBase<CwiseBinaryOp<ei_scalar_sum_op<Scalar>,Transpose<DerivedA>,DerivedB> >& >(other));
+}
+
+template<typename Derived>
+template<typename DerivedA, typename DerivedB>
+Derived& MatrixBase<Derived>::
+lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,Transpose<DerivedB> >& other)
+{
+  ei_assert(ei_extract_data(derived()) != ei_extract_data(other.rhs())
+            && "aliasing detected during tranposition, please evaluate your expression");
+  return lazyAssign(static_cast<const MatrixBase<CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,Transpose<DerivedB> > >& >(other));
+}
+
+template<typename Derived>
+template<typename OtherDerived> Derived&
+MatrixBase<Derived>::
+lazyAssign(const CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<OtherDerived> > >& other)
+{
+  ei_assert(ei_extract_data(other) != ei_extract_data(derived())
+            && "aliasing detected during tranposition, please use adjointInPlace()");
+  return lazyAssign(static_cast<const MatrixBase<CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<OtherDerived> > > >& >(other));
+}
+
+template<typename Derived>
+template<typename DerivedA, typename DerivedB>
+Derived& MatrixBase<Derived>::
+lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedA> > >,DerivedB>& other)
+{
+  ei_assert(ei_extract_data(derived()) != ei_extract_data(other.lhs())
+            && "aliasing detected during tranposition, please evaluate your expression");
+  return lazyAssign(static_cast<const MatrixBase<CwiseBinaryOp<ei_scalar_sum_op<Scalar>,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedA> > >,DerivedB> >& >(other));
+}
+
+template<typename Derived>
+template<typename DerivedA, typename DerivedB>
+Derived& MatrixBase<Derived>::
+lazyAssign(const CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedB> > > >& other)
+{
+  ei_assert(ei_extract_data(derived()) != ei_extract_data(other.rhs())
+            && "aliasing detected during tranposition, please evaluate your expression");
+  return lazyAssign(static_cast<const MatrixBase<CwiseBinaryOp<ei_scalar_sum_op<Scalar>,DerivedA,CwiseUnaryOp<ei_scalar_conjugate_op<Scalar>, NestByValue<Eigen::Transpose<DerivedB> > > > >& >(other));
+}
+#endif
+    
 #endif // EIGEN_TRANSPOSE_H

test/adjoint.cpp

@@ -21,7 +21,7 @@
 // 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/>.
-#define EIGEN_NO_ASSERTION_CHECKING
+
 #include "main.h"
 
 template<typename MatrixType> void adjoint(const MatrixType& m)
@@ -110,7 +110,6 @@ template<typename MatrixType> void adjoint(const MatrixType& m)
   m3.transposeInPlace();
   VERIFY_IS_APPROX(m3,m1.conjugate());
 
-
 }
 
 void test_adjoint()
@@ -125,5 +124,16 @@ void test_adjoint()
   }
   // test a large matrix only once
   CALL_SUBTEST( adjoint(Matrix<float, 100, 100>()) );
+  
+  {
+    MatrixXcf a(10,10), b(10,10);
+    VERIFY_RAISES_ASSERT(a = a.transpose());
+    VERIFY_RAISES_ASSERT(a = a.transpose() + b);
+    VERIFY_RAISES_ASSERT(a = b + a.transpose());
+    VERIFY_RAISES_ASSERT(a = a.conjugate().transpose());
+    VERIFY_RAISES_ASSERT(a = a.adjoint());
+    VERIFY_RAISES_ASSERT(a = a.adjoint() + b);
+    VERIFY_RAISES_ASSERT(a = b + a.adjoint());
+  }
 }