1 // This file is part of Eigen, a lightweight C++ template library 2 // for linear algebra. 3 // 4 // Copyright (C) 2006-2009 Benoit Jacob <jacob.benoit.1 (at) gmail.com> 5 // Copyright (C) 2008 Gael Guennebaud <gael.guennebaud (at) inria.fr> 6 // 7 // This Source Code Form is subject to the terms of the Mozilla 8 // Public License v. 2.0. If a copy of the MPL was not distributed 9 // with this file, You can obtain one at http://mozilla.org/MPL/2.0/. 10 11 #ifndef EIGEN_MATRIXBASE_H 12 #define EIGEN_MATRIXBASE_H 13 14 namespace Eigen { 15 16 /** \class MatrixBase 17 * \ingroup Core_Module 18 * 19 * \brief Base class for all dense matrices, vectors, and expressions 20 * 21 * This class is the base that is inherited by all matrix, vector, and related expression 22 * types. Most of the Eigen API is contained in this class, and its base classes. Other important 23 * classes for the Eigen API are Matrix, and VectorwiseOp. 24 * 25 * Note that some methods are defined in other modules such as the \ref LU_Module LU module 26 * for all functions related to matrix inversions. 27 * 28 * \tparam Derived is the derived type, e.g. a matrix type, or an expression, etc. 29 * 30 * When writing a function taking Eigen objects as argument, if you want your function 31 * to take as argument any matrix, vector, or expression, just let it take a 32 * MatrixBase argument. As an example, here is a function printFirstRow which, given 33 * a matrix, vector, or expression \a x, prints the first row of \a x. 34 * 35 * \code 36 template<typename Derived> 37 void printFirstRow(const Eigen::MatrixBase<Derived>& x) 38 { 39 cout << x.row(0) << endl; 40 } 41 * \endcode 42 * 43 * This class can be extended with the help of the plugin mechanism described on the page 44 * \ref TopicCustomizingEigen by defining the preprocessor symbol \c EIGEN_MATRIXBASE_PLUGIN. 45 * 46 * \sa \ref TopicClassHierarchy 47 */ 48 template<typename Derived> class MatrixBase 49 : public DenseBase<Derived> 50 { 51 public: 52 #ifndef EIGEN_PARSED_BY_DOXYGEN 53 typedef MatrixBase StorageBaseType; 54 typedef typename internal::traits<Derived>::StorageKind StorageKind; 55 typedef typename internal::traits<Derived>::Index Index; 56 typedef typename internal::traits<Derived>::Scalar Scalar; 57 typedef typename internal::packet_traits<Scalar>::type PacketScalar; 58 typedef typename NumTraits<Scalar>::Real RealScalar; 59 60 typedef DenseBase<Derived> Base; 61 using Base::RowsAtCompileTime; 62 using Base::ColsAtCompileTime; 63 using Base::SizeAtCompileTime; 64 using Base::MaxRowsAtCompileTime; 65 using Base::MaxColsAtCompileTime; 66 using Base::MaxSizeAtCompileTime; 67 using Base::IsVectorAtCompileTime; 68 using Base::Flags; 69 using Base::CoeffReadCost; 70 71 using Base::derived; 72 using Base::const_cast_derived; 73 using Base::rows; 74 using Base::cols; 75 using Base::size; 76 using Base::coeff; 77 using Base::coeffRef; 78 using Base::lazyAssign; 79 using Base::eval; 80 using Base::operator+=; 81 using Base::operator-=; 82 using Base::operator*=; 83 using Base::operator/=; 84 85 typedef typename Base::CoeffReturnType CoeffReturnType; 86 typedef typename Base::ConstTransposeReturnType ConstTransposeReturnType; 87 typedef typename Base::RowXpr RowXpr; 88 typedef typename Base::ColXpr ColXpr; 89 #endif // not EIGEN_PARSED_BY_DOXYGEN 90 91 92 93 #ifndef EIGEN_PARSED_BY_DOXYGEN 94 /** type of the equivalent square matrix */ 95 typedef Matrix<Scalar,EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime), 96 EIGEN_SIZE_MAX(RowsAtCompileTime,ColsAtCompileTime)> SquareMatrixType; 97 #endif // not EIGEN_PARSED_BY_DOXYGEN 98 99 /** \returns the size of the main diagonal, which is min(rows(),cols()). 100 * \sa rows(), cols(), SizeAtCompileTime. */ 101 inline Index diagonalSize() const { return (std::min)(rows(),cols()); } 102 103 /** \brief The plain matrix type corresponding to this expression. 104 * 105 * This is not necessarily exactly the return type of eval(). In the case of plain matrices, 106 * the return type of eval() is a const reference to a matrix, not a matrix! It is however guaranteed 107 * that the return type of eval() is either PlainObject or const PlainObject&. 108 */ 109 typedef Matrix<typename internal::traits<Derived>::Scalar, 110 internal::traits<Derived>::RowsAtCompileTime, 111 internal::traits<Derived>::ColsAtCompileTime, 112 AutoAlign | (internal::traits<Derived>::Flags&RowMajorBit ? RowMajor : ColMajor), 113 internal::traits<Derived>::MaxRowsAtCompileTime, 114 internal::traits<Derived>::MaxColsAtCompileTime 115 > PlainObject; 116 117 #ifndef EIGEN_PARSED_BY_DOXYGEN 118 /** \internal Represents a matrix with all coefficients equal to one another*/ 119 typedef CwiseNullaryOp<internal::scalar_constant_op<Scalar>,Derived> ConstantReturnType; 120 /** \internal the return type of MatrixBase::adjoint() */ 121 typedef typename internal::conditional<NumTraits<Scalar>::IsComplex, 122 CwiseUnaryOp<internal::scalar_conjugate_op<Scalar>, ConstTransposeReturnType>, 123 ConstTransposeReturnType 124 >::type AdjointReturnType; 125 /** \internal Return type of eigenvalues() */ 126 typedef Matrix<std::complex<RealScalar>, internal::traits<Derived>::ColsAtCompileTime, 1, ColMajor> EigenvaluesReturnType; 127 /** \internal the return type of identity */ 128 typedef CwiseNullaryOp<internal::scalar_identity_op<Scalar>,Derived> IdentityReturnType; 129 /** \internal the return type of unit vectors */ 130 typedef Block<const CwiseNullaryOp<internal::scalar_identity_op<Scalar>, SquareMatrixType>, 131 internal::traits<Derived>::RowsAtCompileTime, 132 internal::traits<Derived>::ColsAtCompileTime> BasisReturnType; 133 #endif // not EIGEN_PARSED_BY_DOXYGEN 134 135 #define EIGEN_CURRENT_STORAGE_BASE_CLASS Eigen::MatrixBase 136 # include "../plugins/CommonCwiseUnaryOps.h" 137 # include "../plugins/CommonCwiseBinaryOps.h" 138 # include "../plugins/MatrixCwiseUnaryOps.h" 139 # include "../plugins/MatrixCwiseBinaryOps.h" 140 # ifdef EIGEN_MATRIXBASE_PLUGIN 141 # include EIGEN_MATRIXBASE_PLUGIN 142 # endif 143 #undef EIGEN_CURRENT_STORAGE_BASE_CLASS 144 145 /** Special case of the template operator=, in order to prevent the compiler 146 * from generating a default operator= (issue hit with g++ 4.1) 147 */ 148 Derived& operator=(const MatrixBase& other); 149 150 // We cannot inherit here via Base::operator= since it is causing 151 // trouble with MSVC. 152 153 template <typename OtherDerived> 154 Derived& operator=(const DenseBase<OtherDerived>& other); 155 156 template <typename OtherDerived> 157 Derived& operator=(const EigenBase<OtherDerived>& other); 158 159 template<typename OtherDerived> 160 Derived& operator=(const ReturnByValue<OtherDerived>& other); 161 162 #ifndef EIGEN_PARSED_BY_DOXYGEN 163 template<typename ProductDerived, typename Lhs, typename Rhs> 164 Derived& lazyAssign(const ProductBase<ProductDerived, Lhs,Rhs>& other); 165 166 template<typename MatrixPower, typename Lhs, typename Rhs> 167 Derived& lazyAssign(const MatrixPowerProduct<MatrixPower, Lhs,Rhs>& other); 168 #endif // not EIGEN_PARSED_BY_DOXYGEN 169 170 template<typename OtherDerived> 171 Derived& operator+=(const MatrixBase<OtherDerived>& other); 172 template<typename OtherDerived> 173 Derived& operator-=(const MatrixBase<OtherDerived>& other); 174 175 template<typename OtherDerived> 176 const typename ProductReturnType<Derived,OtherDerived>::Type 177 operator*(const MatrixBase<OtherDerived> &other) const; 178 179 template<typename OtherDerived> 180 const typename LazyProductReturnType<Derived,OtherDerived>::Type 181 lazyProduct(const MatrixBase<OtherDerived> &other) const; 182 183 template<typename OtherDerived> 184 Derived& operator*=(const EigenBase<OtherDerived>& other); 185 186 template<typename OtherDerived> 187 void applyOnTheLeft(const EigenBase<OtherDerived>& other); 188 189 template<typename OtherDerived> 190 void applyOnTheRight(const EigenBase<OtherDerived>& other); 191 192 template<typename DiagonalDerived> 193 const DiagonalProduct<Derived, DiagonalDerived, OnTheRight> 194 operator*(const DiagonalBase<DiagonalDerived> &diagonal) const; 195 196 template<typename OtherDerived> 197 typename internal::scalar_product_traits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType 198 dot(const MatrixBase<OtherDerived>& other) const; 199 200 #ifdef EIGEN2_SUPPORT 201 template<typename OtherDerived> 202 Scalar eigen2_dot(const MatrixBase<OtherDerived>& other) const; 203 #endif 204 205 RealScalar squaredNorm() const; 206 RealScalar norm() const; 207 RealScalar stableNorm() const; 208 RealScalar blueNorm() const; 209 RealScalar hypotNorm() const; 210 const PlainObject normalized() const; 211 void normalize(); 212 213 const AdjointReturnType adjoint() const; 214 void adjointInPlace(); 215 216 typedef Diagonal<Derived> DiagonalReturnType; 217 DiagonalReturnType diagonal(); 218 typedef typename internal::add_const<Diagonal<const Derived> >::type ConstDiagonalReturnType; 219 ConstDiagonalReturnType diagonal() const; 220 221 template<int Index> struct DiagonalIndexReturnType { typedef Diagonal<Derived,Index> Type; }; 222 template<int Index> struct ConstDiagonalIndexReturnType { typedef const Diagonal<const Derived,Index> Type; }; 223 224 template<int Index> typename DiagonalIndexReturnType<Index>::Type diagonal(); 225 template<int Index> typename ConstDiagonalIndexReturnType<Index>::Type diagonal() const; 226 227 // Note: The "MatrixBase::" prefixes are added to help MSVC9 to match these declarations with the later implementations. 228 // On the other hand they confuse MSVC8... 229 #if (defined _MSC_VER) && (_MSC_VER >= 1500) // 2008 or later 230 typename MatrixBase::template DiagonalIndexReturnType<DynamicIndex>::Type diagonal(Index index); 231 typename MatrixBase::template ConstDiagonalIndexReturnType<DynamicIndex>::Type diagonal(Index index) const; 232 #else 233 typename DiagonalIndexReturnType<DynamicIndex>::Type diagonal(Index index); 234 typename ConstDiagonalIndexReturnType<DynamicIndex>::Type diagonal(Index index) const; 235 #endif 236 237 #ifdef EIGEN2_SUPPORT 238 template<unsigned int Mode> typename internal::eigen2_part_return_type<Derived, Mode>::type part(); 239 template<unsigned int Mode> const typename internal::eigen2_part_return_type<Derived, Mode>::type part() const; 240 241 // huuuge hack. make Eigen2's matrix.part<Diagonal>() work in eigen3. Problem: Diagonal is now a class template instead 242 // of an integer constant. Solution: overload the part() method template wrt template parameters list. 243 template<template<typename T, int N> class U> 244 const DiagonalWrapper<ConstDiagonalReturnType> part() const 245 { return diagonal().asDiagonal(); } 246 #endif // EIGEN2_SUPPORT 247 248 template<unsigned int Mode> struct TriangularViewReturnType { typedef TriangularView<Derived, Mode> Type; }; 249 template<unsigned int Mode> struct ConstTriangularViewReturnType { typedef const TriangularView<const Derived, Mode> Type; }; 250 251 template<unsigned int Mode> typename TriangularViewReturnType<Mode>::Type triangularView(); 252 template<unsigned int Mode> typename ConstTriangularViewReturnType<Mode>::Type triangularView() const; 253 254 template<unsigned int UpLo> struct SelfAdjointViewReturnType { typedef SelfAdjointView<Derived, UpLo> Type; }; 255 template<unsigned int UpLo> struct ConstSelfAdjointViewReturnType { typedef const SelfAdjointView<const Derived, UpLo> Type; }; 256 257 template<unsigned int UpLo> typename SelfAdjointViewReturnType<UpLo>::Type selfadjointView(); 258 template<unsigned int UpLo> typename ConstSelfAdjointViewReturnType<UpLo>::Type selfadjointView() const; 259 260 const SparseView<Derived> sparseView(const Scalar& m_reference = Scalar(0), 261 const typename NumTraits<Scalar>::Real& m_epsilon = NumTraits<Scalar>::dummy_precision()) const; 262 static const IdentityReturnType Identity(); 263 static const IdentityReturnType Identity(Index rows, Index cols); 264 static const BasisReturnType Unit(Index size, Index i); 265 static const BasisReturnType Unit(Index i); 266 static const BasisReturnType UnitX(); 267 static const BasisReturnType UnitY(); 268 static const BasisReturnType UnitZ(); 269 static const BasisReturnType UnitW(); 270 271 const DiagonalWrapper<const Derived> asDiagonal() const; 272 const PermutationWrapper<const Derived> asPermutation() const; 273 274 Derived& setIdentity(); 275 Derived& setIdentity(Index rows, Index cols); 276 277 bool isIdentity(const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 278 bool isDiagonal(const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 279 280 bool isUpperTriangular(const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 281 bool isLowerTriangular(const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 282 283 template<typename OtherDerived> 284 bool isOrthogonal(const MatrixBase<OtherDerived>& other, 285 const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 286 bool isUnitary(const RealScalar& prec = NumTraits<Scalar>::dummy_precision()) const; 287 288 /** \returns true if each coefficients of \c *this and \a other are all exactly equal. 289 * \warning When using floating point scalar values you probably should rather use a 290 * fuzzy comparison such as isApprox() 291 * \sa isApprox(), operator!= */ 292 template<typename OtherDerived> 293 inline bool operator==(const MatrixBase<OtherDerived>& other) const 294 { return cwiseEqual(other).all(); } 295 296 /** \returns true if at least one pair of coefficients of \c *this and \a other are not exactly equal to each other. 297 * \warning When using floating point scalar values you probably should rather use a 298 * fuzzy comparison such as isApprox() 299 * \sa isApprox(), operator== */ 300 template<typename OtherDerived> 301 inline bool operator!=(const MatrixBase<OtherDerived>& other) const 302 { return cwiseNotEqual(other).any(); } 303 304 NoAlias<Derived,Eigen::MatrixBase > noalias(); 305 306 inline const ForceAlignedAccess<Derived> forceAlignedAccess() const; 307 inline ForceAlignedAccess<Derived> forceAlignedAccess(); 308 template<bool Enable> inline typename internal::add_const_on_value_type<typename internal::conditional<Enable,ForceAlignedAccess<Derived>,Derived&>::type>::type forceAlignedAccessIf() const; 309 template<bool Enable> inline typename internal::conditional<Enable,ForceAlignedAccess<Derived>,Derived&>::type forceAlignedAccessIf(); 310 311 Scalar trace() const; 312 313 /////////// Array module /////////// 314 315 template<int p> RealScalar lpNorm() const; 316 317 MatrixBase<Derived>& matrix() { return *this; } 318 const MatrixBase<Derived>& matrix() const { return *this; } 319 320 /** \returns an \link Eigen::ArrayBase Array \endlink expression of this matrix 321 * \sa ArrayBase::matrix() */ 322 ArrayWrapper<Derived> array() { return derived(); } 323 const ArrayWrapper<const Derived> array() const { return derived(); } 324 325 /////////// LU module /////////// 326 327 const FullPivLU<PlainObject> fullPivLu() const; 328 const PartialPivLU<PlainObject> partialPivLu() const; 329 330 #if EIGEN2_SUPPORT_STAGE < STAGE20_RESOLVE_API_CONFLICTS 331 const LU<PlainObject> lu() const; 332 #endif 333 334 #ifdef EIGEN2_SUPPORT 335 const LU<PlainObject> eigen2_lu() const; 336 #endif 337 338 #if EIGEN2_SUPPORT_STAGE > STAGE20_RESOLVE_API_CONFLICTS 339 const PartialPivLU<PlainObject> lu() const; 340 #endif 341 342 #ifdef EIGEN2_SUPPORT 343 template<typename ResultType> 344 void computeInverse(MatrixBase<ResultType> *result) const { 345 *result = this->inverse(); 346 } 347 #endif 348 349 const internal::inverse_impl<Derived> inverse() const; 350 template<typename ResultType> 351 void computeInverseAndDetWithCheck( 352 ResultType& inverse, 353 typename ResultType::Scalar& determinant, 354 bool& invertible, 355 const RealScalar& absDeterminantThreshold = NumTraits<Scalar>::dummy_precision() 356 ) const; 357 template<typename ResultType> 358 void computeInverseWithCheck( 359 ResultType& inverse, 360 bool& invertible, 361 const RealScalar& absDeterminantThreshold = NumTraits<Scalar>::dummy_precision() 362 ) const; 363 Scalar determinant() const; 364 365 /////////// Cholesky module /////////// 366 367 const LLT<PlainObject> llt() const; 368 const LDLT<PlainObject> ldlt() const; 369 370 /////////// QR module /////////// 371 372 const HouseholderQR<PlainObject> householderQr() const; 373 const ColPivHouseholderQR<PlainObject> colPivHouseholderQr() const; 374 const FullPivHouseholderQR<PlainObject> fullPivHouseholderQr() const; 375 376 #ifdef EIGEN2_SUPPORT 377 const QR<PlainObject> qr() const; 378 #endif 379 380 EigenvaluesReturnType eigenvalues() const; 381 RealScalar operatorNorm() const; 382 383 /////////// SVD module /////////// 384 385 JacobiSVD<PlainObject> jacobiSvd(unsigned int computationOptions = 0) const; 386 387 #ifdef EIGEN2_SUPPORT 388 SVD<PlainObject> svd() const; 389 #endif 390 391 /////////// Geometry module /////////// 392 393 #ifndef EIGEN_PARSED_BY_DOXYGEN 394 /// \internal helper struct to form the return type of the cross product 395 template<typename OtherDerived> struct cross_product_return_type { 396 typedef typename internal::scalar_product_traits<typename internal::traits<Derived>::Scalar,typename internal::traits<OtherDerived>::Scalar>::ReturnType Scalar; 397 typedef Matrix<Scalar,MatrixBase::RowsAtCompileTime,MatrixBase::ColsAtCompileTime> type; 398 }; 399 #endif // EIGEN_PARSED_BY_DOXYGEN 400 template<typename OtherDerived> 401 typename cross_product_return_type<OtherDerived>::type 402 cross(const MatrixBase<OtherDerived>& other) const; 403 template<typename OtherDerived> 404 PlainObject cross3(const MatrixBase<OtherDerived>& other) const; 405 PlainObject unitOrthogonal(void) const; 406 Matrix<Scalar,3,1> eulerAngles(Index a0, Index a1, Index a2) const; 407 408 #if EIGEN2_SUPPORT_STAGE > STAGE20_RESOLVE_API_CONFLICTS 409 ScalarMultipleReturnType operator*(const UniformScaling<Scalar>& s) const; 410 // put this as separate enum value to work around possible GCC 4.3 bug (?) 411 enum { HomogeneousReturnTypeDirection = ColsAtCompileTime==1?Vertical:Horizontal }; 412 typedef Homogeneous<Derived, HomogeneousReturnTypeDirection> HomogeneousReturnType; 413 HomogeneousReturnType homogeneous() const; 414 #endif 415 416 enum { 417 SizeMinusOne = SizeAtCompileTime==Dynamic ? Dynamic : SizeAtCompileTime-1 418 }; 419 typedef Block<const Derived, 420 internal::traits<Derived>::ColsAtCompileTime==1 ? SizeMinusOne : 1, 421 internal::traits<Derived>::ColsAtCompileTime==1 ? 1 : SizeMinusOne> ConstStartMinusOne; 422 typedef CwiseUnaryOp<internal::scalar_quotient1_op<typename internal::traits<Derived>::Scalar>, 423 const ConstStartMinusOne > HNormalizedReturnType; 424 425 const HNormalizedReturnType hnormalized() const; 426 427 ////////// Householder module /////////// 428 429 void makeHouseholderInPlace(Scalar& tau, RealScalar& beta); 430 template<typename EssentialPart> 431 void makeHouseholder(EssentialPart& essential, 432 Scalar& tau, RealScalar& beta) const; 433 template<typename EssentialPart> 434 void applyHouseholderOnTheLeft(const EssentialPart& essential, 435 const Scalar& tau, 436 Scalar* workspace); 437 template<typename EssentialPart> 438 void applyHouseholderOnTheRight(const EssentialPart& essential, 439 const Scalar& tau, 440 Scalar* workspace); 441 442 ///////// Jacobi module ///////// 443 444 template<typename OtherScalar> 445 void applyOnTheLeft(Index p, Index q, const JacobiRotation<OtherScalar>& j); 446 template<typename OtherScalar> 447 void applyOnTheRight(Index p, Index q, const JacobiRotation<OtherScalar>& j); 448 449 ///////// MatrixFunctions module ///////// 450 451 typedef typename internal::stem_function<Scalar>::type StemFunction; 452 const MatrixExponentialReturnValue<Derived> exp() const; 453 const MatrixFunctionReturnValue<Derived> matrixFunction(StemFunction f) const; 454 const MatrixFunctionReturnValue<Derived> cosh() const; 455 const MatrixFunctionReturnValue<Derived> sinh() const; 456 const MatrixFunctionReturnValue<Derived> cos() const; 457 const MatrixFunctionReturnValue<Derived> sin() const; 458 const MatrixSquareRootReturnValue<Derived> sqrt() const; 459 const MatrixLogarithmReturnValue<Derived> log() const; 460 const MatrixPowerReturnValue<Derived> pow(const RealScalar& p) const; 461 462 #ifdef EIGEN2_SUPPORT 463 template<typename ProductDerived, typename Lhs, typename Rhs> 464 Derived& operator+=(const Flagged<ProductBase<ProductDerived, Lhs,Rhs>, 0, 465 EvalBeforeAssigningBit>& other); 466 467 template<typename ProductDerived, typename Lhs, typename Rhs> 468 Derived& operator-=(const Flagged<ProductBase<ProductDerived, Lhs,Rhs>, 0, 469 EvalBeforeAssigningBit>& other); 470 471 /** \deprecated because .lazy() is deprecated 472 * Overloaded for cache friendly product evaluation */ 473 template<typename OtherDerived> 474 Derived& lazyAssign(const Flagged<OtherDerived, 0, EvalBeforeAssigningBit>& other) 475 { return lazyAssign(other._expression()); } 476 477 template<unsigned int Added> 478 const Flagged<Derived, Added, 0> marked() const; 479 const Flagged<Derived, 0, EvalBeforeAssigningBit> lazy() const; 480 481 inline const Cwise<Derived> cwise() const; 482 inline Cwise<Derived> cwise(); 483 484 VectorBlock<Derived> start(Index size); 485 const VectorBlock<const Derived> start(Index size) const; 486 VectorBlock<Derived> end(Index size); 487 const VectorBlock<const Derived> end(Index size) const; 488 template<int Size> VectorBlock<Derived,Size> start(); 489 template<int Size> const VectorBlock<const Derived,Size> start() const; 490 template<int Size> VectorBlock<Derived,Size> end(); 491 template<int Size> const VectorBlock<const Derived,Size> end() const; 492 493 Minor<Derived> minor(Index row, Index col); 494 const Minor<Derived> minor(Index row, Index col) const; 495 #endif 496 497 protected: 498 MatrixBase() : Base() {} 499 500 private: 501 explicit MatrixBase(int); 502 MatrixBase(int,int); 503 template<typename OtherDerived> explicit MatrixBase(const MatrixBase<OtherDerived>&); 504 protected: 505 // mixing arrays and matrices is not legal 506 template<typename OtherDerived> Derived& operator+=(const ArrayBase<OtherDerived>& ) 507 {EIGEN_STATIC_ASSERT(std::ptrdiff_t(sizeof(typename OtherDerived::Scalar))==-1,YOU_CANNOT_MIX_ARRAYS_AND_MATRICES); return *this;} 508 // mixing arrays and matrices is not legal 509 template<typename OtherDerived> Derived& operator-=(const ArrayBase<OtherDerived>& ) 510 {EIGEN_STATIC_ASSERT(std::ptrdiff_t(sizeof(typename OtherDerived::Scalar))==-1,YOU_CANNOT_MIX_ARRAYS_AND_MATRICES); return *this;} 511 }; 512 513 514 /*************************************************************************** 515 * Implementation of matrix base methods 516 ***************************************************************************/ 517 518 /** replaces \c *this by \c *this * \a other. 519 * 520 * \returns a reference to \c *this 521 * 522 * Example: \include MatrixBase_applyOnTheRight.cpp 523 * Output: \verbinclude MatrixBase_applyOnTheRight.out 524 */ 525 template<typename Derived> 526 template<typename OtherDerived> 527 inline Derived& 528 MatrixBase<Derived>::operator*=(const EigenBase<OtherDerived> &other) 529 { 530 other.derived().applyThisOnTheRight(derived()); 531 return derived(); 532 } 533 534 /** replaces \c *this by \c *this * \a other. It is equivalent to MatrixBase::operator*=(). 535 * 536 * Example: \include MatrixBase_applyOnTheRight.cpp 537 * Output: \verbinclude MatrixBase_applyOnTheRight.out 538 */ 539 template<typename Derived> 540 template<typename OtherDerived> 541 inline void MatrixBase<Derived>::applyOnTheRight(const EigenBase<OtherDerived> &other) 542 { 543 other.derived().applyThisOnTheRight(derived()); 544 } 545 546 /** replaces \c *this by \a other * \c *this. 547 * 548 * Example: \include MatrixBase_applyOnTheLeft.cpp 549 * Output: \verbinclude MatrixBase_applyOnTheLeft.out 550 */ 551 template<typename Derived> 552 template<typename OtherDerived> 553 inline void MatrixBase<Derived>::applyOnTheLeft(const EigenBase<OtherDerived> &other) 554 { 555 other.derived().applyThisOnTheLeft(derived()); 556 } 557 558 } // end namespace Eigen 559 560 #endif // EIGEN_MATRIXBASE_H 561