From a2f8d4be6a334ff5caaf782d804330d545b58b53 Mon Sep 17 00:00:00 2001 From: Benoit Jacob Date: Fri, 29 Feb 2008 14:35:14 +0000 Subject: [PATCH] Patch by Gael Guennebaud: coeff-wise binary operators. This unifies + and - and moreover this patch introduces coeff-wise * and / based on this. Also, corresponding test. --- Eigen/Core | 3 +- Eigen/src/Core/Block.h | 1 + Eigen/src/Core/CwiseBinaryOp.h | 222 +++++++++++++++++++++++++++ Eigen/src/Core/Dot.h | 6 +- Eigen/src/Core/ForwardDeclarations.h | 5 +- Eigen/src/Core/MatrixBase.h | 17 +- test/CMakeLists.txt | 1 + test/cwiseop.cpp | 84 ++++++++++ test/main.h | 1 + 9 files changed, 331 insertions(+), 9 deletions(-) create mode 100644 Eigen/src/Core/CwiseBinaryOp.h create mode 100644 test/cwiseop.cpp diff --git a/Eigen/Core b/Eigen/Core index 8cdd812ac..3ced3e88a 100644 --- a/Eigen/Core +++ b/Eigen/Core @@ -18,8 +18,7 @@ namespace Eigen { #include "src/Core/Matrix.h" #include "src/Core/Cast.h" #include "src/Core/Eval.h" -#include "src/Core/Sum.h" -#include "src/Core/Difference.h" +#include "src/Core/CwiseBinaryOp.h" #include "src/Core/Opposite.h" #include "src/Core/ScalarMultiple.h" #include "src/Core/Product.h" diff --git a/Eigen/src/Core/Block.h b/Eigen/src/Core/Block.h index 4771c94cb..9ed69d03f 100644 --- a/Eigen/src/Core/Block.h +++ b/Eigen/src/Core/Block.h @@ -1,6 +1,7 @@ // This file is part of Eigen, a lightweight C++ template library // for linear algebra. Eigen itself is part of the KDE project. // +// Copyright (C) 2008 Gael Guennebaud // Copyright (C) 2006-2008 Benoit Jacob // // Eigen is free software; you can redistribute it and/or diff --git a/Eigen/src/Core/CwiseBinaryOp.h b/Eigen/src/Core/CwiseBinaryOp.h new file mode 100644 index 000000000..8a344b540 --- /dev/null +++ b/Eigen/src/Core/CwiseBinaryOp.h @@ -0,0 +1,222 @@ +// This file is part of Eigen, a lightweight C++ template library +// for linear algebra. Eigen itself is part of the KDE project. +// +// Copyright (C) 2008 Gael Guennebaud +// Copyright (C) 2006-2008 Benoit Jacob +// +// 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 . + +#ifndef EIGEN_CWISE_BINARY_OP_H +#define EIGEN_CWISE_BINARY_OP_H + +/** \class CwiseBinaryOp + * + * \brief Generic expression of a coefficient-wise operator between two matrices or vectors + * + * \param BinaryOp template functor implementing the operator + * \param Lhs the type of the left-hand side + * \param Rhs the type of the right-hand side + * + * This class represents an expression of a generic binary operator of two matrices or vectors. + * It is the return type of the operator+, operator-, cwiseiseProduct, cwiseQuotient between matrices or vectors, and most + * of the time this is the only way it is used. + * + * \sa class CwiseProductOp, class CwiseQuotientOp + */ +template class BinaryOp, typename Lhs, typename Rhs> +class CwiseBinaryOp : NoOperatorEquals, + public MatrixBase > +{ + public: + typedef typename Lhs::Scalar Scalar; + typedef typename Lhs::Ref LhsRef; + typedef typename Rhs::Ref RhsRef; + friend class MatrixBase; + typedef MatrixBase Base; + + CwiseBinaryOp(const LhsRef& lhs, const RhsRef& rhs) + : m_lhs(lhs), m_rhs(rhs) + { + assert(lhs.rows() == rhs.rows() && lhs.cols() == rhs.cols()); + } + + private: + enum { + RowsAtCompileTime = Lhs::Traits::RowsAtCompileTime, + ColsAtCompileTime = Lhs::Traits::ColsAtCompileTime, + MaxRowsAtCompileTime = Lhs::Traits::MaxRowsAtCompileTime, + MaxColsAtCompileTime = Lhs::Traits::MaxColsAtCompileTime + }; + + const CwiseBinaryOp& _ref() const { return *this; } + int _rows() const { return m_lhs.rows(); } + int _cols() const { return m_lhs.cols(); } + + Scalar _coeff(int row, int col) const + { + return BinaryOp::op(m_lhs.coeff(row, col), m_rhs.coeff(row, col)); + } + + protected: + const LhsRef m_lhs; + const RhsRef m_rhs; +}; + +/** \brief Template functor to compute the sum of two scalars + * + * \sa class CwiseBinaryOp, MatrixBase::operator+ + */ +template struct CwiseSumOp { + static Scalar op(const Scalar& a, const Scalar& b) { return a + b; } +}; + +/** \brief Template functor to compute the difference of two scalars + * + * \sa class CwiseBinaryOp, MatrixBase::operator- + */ +template struct CwiseDifferenceOp { + static Scalar op(const Scalar& a, const Scalar& b) { return a - b; } +}; + +/** \brief Template functor to compute the product of two scalars + * + * \sa class CwiseBinaryOp, MatrixBase::cwiseProduct() + */ +template struct CwiseProductOp { + static Scalar op(const Scalar& a, const Scalar& b) { return a * b; } +}; + +/** \brief Template functor to compute the quotient of two scalars + * + * \sa class CwiseBinaryOp, MatrixBase::cwiseQuotient() + */ +template struct CwiseQuotientOp { + static Scalar op(const Scalar& a, const Scalar& b) { return a / b; } +}; + +/** \relates MatrixBase + * + * \returns an expression of the difference of \a mat1 and \a mat2 + * + * \sa class CwiseBinaryOp, MatrixBase::operator-=() + */ +template +const CwiseBinaryOp +operator-(const MatrixBase &mat1, const MatrixBase &mat2) +{ + return CwiseBinaryOp(mat1.ref(), mat2.ref()); +} + +/** replaces \c *this by \c *this - \a other. + * + * \returns a reference to \c *this + */ +template +template +Derived & +MatrixBase::operator-=(const MatrixBase &other) +{ + return *this = *this - other; +} + + +/** \relates MatrixBase + * + * \returns an expression of the sum of \a mat1 and \a mat2 + * + * \sa class CwiseBinaryOp, MatrixBase::operator+=() + */ +template +const CwiseBinaryOp +operator+(const MatrixBase &mat1, const MatrixBase &mat2) +{ + return CwiseBinaryOp(mat1.ref(), mat2.ref()); +} + +/** replaces \c *this by \c *this + \a other. + * + * \returns a reference to \c *this + */ +template +template +Derived & +MatrixBase::operator+=(const MatrixBase& other) +{ + return *this = *this + other; +} + + +/** \returns an expression of the Schur product (coefficient wise product) of *this and \a other + * + * \sa class CwiseBinaryOp + */ +template +template +const CwiseBinaryOp +MatrixBase::cwiseProduct(const MatrixBase &other) const +{ + return CwiseBinaryOp(ref(), other.ref()); +} + + +/** \returns an expression of the coefficient-wise quotient of *this and \a other + * + * \sa class CwiseBinaryOp + */ +template +template +const CwiseBinaryOp +MatrixBase::cwiseQuotient(const MatrixBase &other) const +{ + return CwiseBinaryOp(ref(), other.ref()); +} + + +/** \relates MatrixBase + * + * \returns an expression of a custom coefficient-wise operator of \a mat1 and \a mat2 + * + * \param CustomBinaryOp template functor of the custom operator + * + * \sa class CwiseBinaryOp, MatrixBase::operator+, MatrixBase::operator-, MatrixBase::cwiseProduct, MatrixBase::cwiseQuotient + */ +template class CustomBinaryOp, typename Scalar, typename Derived1, typename Derived2> +const CwiseBinaryOp +cwise(const MatrixBase &mat1, const MatrixBase &mat2) +{ + return CwiseBinaryOp(mat1.ref(), mat2.ref()); +} + +/** \returns an expression of a custom coefficient-wise operator of *this and \a other + * + * \param CustomBinaryOp template functor of the custom operator + * + * \sa class CwiseBinaryOp, MatrixBase::operator+, MatrixBase::operator-, MatrixBase::cwiseProduct, MatrixBase::cwiseQuotient + */ +template +template class CustomBinaryOp, typename OtherDerived> +const CwiseBinaryOp +MatrixBase::cwise(const MatrixBase &other) const +{ + return CwiseBinaryOp(ref(), other.ref()); +} + + +#endif // EIGEN_CWISE_BINARY_OP_H diff --git a/Eigen/src/Core/Dot.h b/Eigen/src/Core/Dot.h index b4a24506f..c3ca1e2ac 100644 --- a/Eigen/src/Core/Dot.h +++ b/Eigen/src/Core/Dot.h @@ -69,7 +69,7 @@ struct DotUnroller */ template template -Scalar MatrixBase::dot(const OtherDerived& other) const +Scalar MatrixBase::dot(const MatrixBase& other) const { assert(Traits::IsVectorAtCompileTime && OtherDerived::Traits::IsVectorAtCompileTime @@ -79,7 +79,7 @@ Scalar MatrixBase::dot(const OtherDerived& other) const && Traits::SizeAtCompileTime != Dynamic && Traits::SizeAtCompileTime <= 16) DotUnroller + Derived, MatrixBase > ::run(*static_cast(this), other, res); else { @@ -136,7 +136,7 @@ MatrixBase::normalized() const template template bool MatrixBase::isOrtho -(const OtherDerived& other, +(const MatrixBase& other, typename NumTraits::Real prec) const { return ei_abs2(dot(other)) <= prec * prec * norm2() * other.norm2(); diff --git a/Eigen/src/Core/ForwardDeclarations.h b/Eigen/src/Core/ForwardDeclarations.h index 3ddcbe0b4..59c088c19 100644 --- a/Eigen/src/Core/ForwardDeclarations.h +++ b/Eigen/src/Core/ForwardDeclarations.h @@ -35,8 +35,9 @@ template clas template class Transpose; template class Conjugate; template class Opposite; -template class Sum; -template class Difference; +template class BinaryOp, typename Lhs, typename Rhs> class CwiseBinaryOp; +template struct CwiseProductOp; +template struct CwiseQuotientOp; template class Product; template class ScalarMultiple; template class Random; diff --git a/Eigen/src/Core/MatrixBase.h b/Eigen/src/Core/MatrixBase.h index ba88e09c9..1d3e90d40 100644 --- a/Eigen/src/Core/MatrixBase.h +++ b/Eigen/src/Core/MatrixBase.h @@ -212,12 +212,13 @@ template class MatrixBase Scalar trace() const; template - Scalar dot(const OtherDerived& other) const; + Scalar dot(const MatrixBase& other) const; RealScalar norm2() const; RealScalar norm() const; const ScalarMultiple normalized() const; template - bool isOrtho(const OtherDerived& other, RealScalar prec = precision()) const; + bool isOrtho(const MatrixBase& other, + RealScalar prec = precision()) const; bool isOrtho(RealScalar prec = precision()) const; static const Eval > random(int rows, int cols); @@ -262,6 +263,18 @@ template class MatrixBase const Opposite operator-() const; + template + const CwiseBinaryOp + cwiseProduct(const MatrixBase &other) const; + + template + const CwiseBinaryOp + cwiseQuotient(const MatrixBase &other) const; + + template class CustomBinaryOp, typename OtherDerived> + const CwiseBinaryOp + cwise(const MatrixBase &other) const; + template Derived& operator+=(const MatrixBase& other); template diff --git a/test/CMakeLists.txt b/test/CMakeLists.txt index 6fff26b43..504b30611 100644 --- a/test/CMakeLists.txt +++ b/test/CMakeLists.txt @@ -5,6 +5,7 @@ FIND_PACKAGE(Qt4 REQUIRED) INCLUDE_DIRECTORIES( ${QT_INCLUDE_DIR} ) SET(test_SRCS + cwiseop.cpp main.cpp basicstuff.cpp linearstructure.cpp diff --git a/test/cwiseop.cpp b/test/cwiseop.cpp new file mode 100644 index 000000000..3e8179aa3 --- /dev/null +++ b/test/cwiseop.cpp @@ -0,0 +1,84 @@ +// 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 +// +// 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 . + +#include "main.h" +#include +#include +#include + +namespace Eigen { + +template struct AddIfNull { + static Scalar op(const Scalar a, const Scalar b) {return a<1e-3 ? b : a;} +}; + +template void cwiseops(const MatrixType& m) +{ + typedef typename MatrixType::Scalar Scalar; + typedef Matrix VectorType; + + int rows = m.rows(); + int cols = m.cols(); + + MatrixType m1 = MatrixType::random(rows, cols), + m2 = MatrixType::random(rows, cols), + m3(rows, cols), + mzero = MatrixType::zero(rows, cols), + mones = MatrixType::ones(rows, cols), + identity = Matrix + ::identity(rows, rows), + square = Matrix + ::random(rows, rows); + VectorType v1 = VectorType::random(rows), + v2 = VectorType::random(rows), + vzero = VectorType::zero(rows); + + + m2 = m2.cwise(mones); + //m2 = cwise(m2,mones); + std::cout << m2 << "\n"; + + VERIFY_IS_APPROX( mzero, m1-m1); + VERIFY_IS_APPROX( m2, m1+m2-m1); + VERIFY_IS_APPROX( mones, m2.cwiseQuotient(m2)); + VERIFY_IS_APPROX( m1.cwiseProduct(m2), m2.cwiseProduct(m1)); + //VERIFY_IS_APPROX( m1, m2.cwiseProduct(m1).cwiseQuotient(m2)); + +// VERIFY_IS_APPROX( cwiseMin(m1,m2), cwiseMin(m2,m1) ); +// VERIFY_IS_APPROX( cwiseMin(m1,m1+mones), m1 ); +// VERIFY_IS_APPROX( cwiseMin(m1,m1-mones), m1-mones ); +} + +void EigenTest::testCwiseops() +{ + for(int i = 0; i < 1/*m_repeat*/ ; i++) { +// cwiseops(Matrix()); +// cwiseops(Matrix4d()); +// cwiseops(MatrixXcf(3, 3)); + cwiseops(MatrixXi(8, 12)); +// cwiseops(MatrixXcd(20, 20)); + } +} + +} // namespace Eigen diff --git a/test/main.h b/test/main.h index cacdd4aeb..54f6039cf 100644 --- a/test/main.h +++ b/test/main.h @@ -121,6 +121,7 @@ class EigenTest : public QObject void testMiscMatrices(); void testSmallVectors(); void testMap(); + void testCwiseops(); protected: int m_repeat; };