diff options
author | 2009-08-15 22:19:29 +0200 | |
---|---|---|
committer | 2009-08-15 22:19:29 +0200 | |
commit | 239ada95b7680c75f793086eaa35fe7ec1047204 (patch) | |
tree | 88a7f107a3866e1c7f284feb818fbbc2cd22ffa3 | |
parent | a3e6047c25a4cbc2153974e04fe124c5776a23c0 (diff) |
add overloads of lazyAssign to detect common aliasing issue with
transpose and adjoint
-rw-r--r-- | Eigen/Core | 2 | ||||
-rw-r--r-- | Eigen/src/Core/CwiseBinaryOp.h | 6 | ||||
-rw-r--r-- | Eigen/src/Core/MatrixBase.h | 15 | ||||
-rw-r--r-- | Eigen/src/Core/Transpose.h | 88 | ||||
-rw-r--r-- | test/adjoint.cpp | 14 |
5 files changed, 122 insertions, 3 deletions
diff --git a/Eigen/Core b/Eigen/Core index ff56c6ca2..854f737d6 100644 --- a/Eigen/Core +++ b/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" diff --git a/Eigen/src/Core/CwiseBinaryOp.h b/Eigen/src/Core/CwiseBinaryOp.h index bcd4ae25f..78de9fbce 100644 --- a/Eigen/src/Core/CwiseBinaryOp.h +++ b/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; diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h index a554de8e3..6ed96e286 100644 --- a/Eigen/src/Core/MatrixBase.h +++ b/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; diff --git a/Eigen/src/Core/Transpose.h b/Eigen/src/Core/Transpose.h index 0787daa61..8527edc2a 100644 --- a/Eigen/src/Core/Transpose.h +++ b/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 diff --git a/test/adjoint.cpp b/test/adjoint.cpp index 300650338..964658c65 100644 --- a/test/adjoint.cpp +++ b/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()); + } } |