11 #ifndef EIGEN_MATRIXBASE_H
12 #define EIGEN_MATRIXBASE_H
52 #ifndef EIGEN_PARSED_BY_DOXYGEN
54 typedef typename internal::traits<Derived>::StorageKind StorageKind;
55 typedef typename internal::traits<Derived>::StorageIndex
StorageIndex;
56 typedef typename internal::traits<Derived>::Scalar
Scalar;
57 typedef typename internal::packet_traits<Scalar>::type PacketScalar;
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;
71 using Base::const_cast_derived;
77 using Base::lazyAssign;
79 using Base::operator+=;
80 using Base::operator-=;
81 using Base::operator*=;
82 using Base::operator/=;
83 using Base::operator*;
84 using Base::operator/;
86 typedef typename Base::CoeffReturnType CoeffReturnType;
88 typedef typename Base::RowXpr RowXpr;
89 typedef typename Base::ColXpr ColXpr;
90 #endif // not EIGEN_PARSED_BY_DOXYGEN
94 #ifndef EIGEN_PARSED_BY_DOXYGEN
98 #endif // not EIGEN_PARSED_BY_DOXYGEN
103 inline Index
diagonalSize()
const {
return (std::min)(rows(),cols()); }
107 #ifndef EIGEN_PARSED_BY_DOXYGEN
111 typedef typename internal::conditional<NumTraits<Scalar>::IsComplex,
113 ConstTransposeReturnType
114 >::type AdjointReturnType;
121 internal::traits<Derived>::RowsAtCompileTime,
122 internal::traits<Derived>::ColsAtCompileTime> BasisReturnType;
123 #endif // not EIGEN_PARSED_BY_DOXYGEN
125 #define EIGEN_CURRENT_STORAGE_BASE_CLASS Eigen::MatrixBase
126 # include "../plugins/CommonCwiseUnaryOps.h"
127 # include "../plugins/CommonCwiseBinaryOps.h"
128 # include "../plugins/MatrixCwiseUnaryOps.h"
129 # include "../plugins/MatrixCwiseBinaryOps.h"
130 # ifdef EIGEN_MATRIXBASE_PLUGIN
131 # include EIGEN_MATRIXBASE_PLUGIN
133 #undef EIGEN_CURRENT_STORAGE_BASE_CLASS
144 template <
typename OtherDerived>
148 template <
typename OtherDerived>
152 template<
typename OtherDerived>
154 Derived&
operator=(
const ReturnByValue<OtherDerived>& other);
156 template<
typename OtherDerived>
159 template<
typename OtherDerived>
164 template<
typename OtherDerived>
171 template<
typename OtherDerived>
177 template<
typename OtherDerived>
182 template<
typename OtherDerived>
185 template<
typename OtherDerived>
188 template<
typename OtherDerived>
191 template<
typename DiagonalDerived>
196 template<
typename OtherDerived>
198 typename internal::scalar_product_traits<typename internal::traits<Derived>::Scalar,
typename internal::traits<OtherDerived>::Scalar>::ReturnType
202 EIGEN_DEVICE_FUNC RealScalar
norm()
const;
206 EIGEN_DEVICE_FUNC
const PlainObject
normalized()
const;
209 EIGEN_DEVICE_FUNC
const AdjointReturnType
adjoint()
const;
218 ConstDiagonalReturnType
diagonal()
const;
225 typename DiagonalIndexReturnType<Index>::Type
diagonal();
229 typename ConstDiagonalIndexReturnType<Index>::Type
diagonal()
const;
231 typedef Diagonal<Derived,DynamicIndex> DiagonalDynamicIndexReturnType;
232 typedef typename internal::add_const<Diagonal<const Derived,DynamicIndex> >::type ConstDiagonalDynamicIndexReturnType;
235 DiagonalDynamicIndexReturnType
diagonal(Index index);
237 ConstDiagonalDynamicIndexReturnType
diagonal(Index index)
const;
239 template<
unsigned int Mode>
struct TriangularViewReturnType {
typedef TriangularView<Derived, Mode> Type; };
240 template<
unsigned int Mode>
struct ConstTriangularViewReturnType {
typedef const TriangularView<const Derived, Mode> Type; };
242 template<
unsigned int Mode>
244 typename TriangularViewReturnType<Mode>::Type triangularView();
245 template<
unsigned int Mode>
247 typename ConstTriangularViewReturnType<Mode>::Type triangularView()
const;
249 template<
unsigned int UpLo>
struct SelfAdjointViewReturnType {
typedef SelfAdjointView<Derived, UpLo> Type; };
250 template<
unsigned int UpLo>
struct ConstSelfAdjointViewReturnType {
typedef const SelfAdjointView<const Derived, UpLo> Type; };
252 template<
unsigned int UpLo>
254 typename SelfAdjointViewReturnType<UpLo>::Type selfadjointView();
255 template<
unsigned int UpLo>
257 typename ConstSelfAdjointViewReturnType<UpLo>::Type selfadjointView()
const;
259 const SparseView<Derived> sparseView(
const Scalar& m_reference =
Scalar(0),
260 const typename NumTraits<Scalar>::Real& m_epsilon = NumTraits<Scalar>::dummy_precision())
const;
261 EIGEN_DEVICE_FUNC
static const IdentityReturnType
Identity();
262 EIGEN_DEVICE_FUNC
static const IdentityReturnType
Identity(Index rows, Index cols);
263 EIGEN_DEVICE_FUNC
static const BasisReturnType
Unit(Index size, Index i);
264 EIGEN_DEVICE_FUNC
static const BasisReturnType
Unit(Index i);
265 EIGEN_DEVICE_FUNC
static const BasisReturnType
UnitX();
266 EIGEN_DEVICE_FUNC
static const BasisReturnType
UnitY();
267 EIGEN_DEVICE_FUNC
static const BasisReturnType
UnitZ();
268 EIGEN_DEVICE_FUNC
static const BasisReturnType
UnitW();
271 const DiagonalWrapper<const Derived>
asDiagonal()
const;
272 const PermutationWrapper<const Derived> asPermutation()
const;
279 bool isIdentity(
const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
280 bool isDiagonal(
const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
282 bool isUpperTriangular(
const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
283 bool isLowerTriangular(
const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
285 template<
typename OtherDerived>
286 bool isOrthogonal(
const MatrixBase<OtherDerived>& other,
287 const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
288 bool isUnitary(
const RealScalar& prec = NumTraits<Scalar>::dummy_precision())
const;
294 template<
typename OtherDerived>
302 template<
typename OtherDerived>
312 template<
bool Enable>
inline const Derived& forceAlignedAccessIf()
const {
return derived(); }
313 template<
bool Enable>
inline Derived& forceAlignedAccessIf() {
return derived(); }
315 EIGEN_DEVICE_FUNC Scalar
trace()
const;
317 template<
int p> EIGEN_DEVICE_FUNC RealScalar lpNorm()
const;
319 EIGEN_DEVICE_FUNC MatrixBase<Derived>& matrix() {
return *
this; }
320 EIGEN_DEVICE_FUNC
const MatrixBase<Derived>& matrix()
const {
return *
this; }
342 template<
typename ResultType>
349 template<
typename ResultType>
380 #ifndef EIGEN_PARSED_BY_DOXYGEN
381 template<
typename OtherDerived>
struct cross_product_return_type {
383 typedef typename internal::scalar_product_traits<typename internal::traits<Derived>::Scalar,
typename internal::traits<OtherDerived>::Scalar>::ReturnType
Scalar;
386 #endif // EIGEN_PARSED_BY_DOXYGEN
387 template<
typename OtherDerived>
389 inline typename cross_product_return_type<OtherDerived>::type
390 cross(
const MatrixBase<OtherDerived>& other)
const;
392 template<
typename OtherDerived>
394 inline PlainObject
cross3(
const MatrixBase<OtherDerived>& other)
const;
399 inline Matrix<Scalar,3,1>
eulerAngles(Index a0, Index a1, Index a2)
const;
401 inline ScalarMultipleReturnType
operator*(
const UniformScaling<Scalar>& s)
const;
405 typedef Homogeneous<Derived, HomogeneousReturnTypeDirection> HomogeneousReturnType;
411 typedef Block<
const Derived,
412 internal::traits<Derived>::ColsAtCompileTime==1 ? SizeMinusOne : 1,
413 internal::traits<Derived>::ColsAtCompileTime==1 ? 1 : SizeMinusOne> ConstStartMinusOne;
414 typedef CwiseUnaryOp<internal::scalar_quotient1_op<typename internal::traits<Derived>::Scalar>,
415 const ConstStartMinusOne > HNormalizedReturnType;
417 inline const HNormalizedReturnType
hnormalized()
const;
422 template<
typename EssentialPart>
424 Scalar& tau, RealScalar& beta)
const;
425 template<
typename EssentialPart>
429 template<
typename EssentialPart>
436 template<
typename OtherScalar>
437 void applyOnTheLeft(Index p, Index q,
const JacobiRotation<OtherScalar>& j);
438 template<
typename OtherScalar>
439 void applyOnTheRight(Index p, Index q,
const JacobiRotation<OtherScalar>& j);
443 template<
typename OtherDerived>
444 EIGEN_STRONG_INLINE
const typename SparseMatrixBase<OtherDerived>::template CwiseProductDenseReturnType<Derived>::Type
445 cwiseProduct(
const SparseMatrixBase<OtherDerived> &other)
const
447 return other.cwiseProduct(derived());
452 typedef typename internal::stem_function<Scalar>::type StemFunction;
453 const MatrixExponentialReturnValue<Derived> exp()
const;
454 const MatrixFunctionReturnValue<Derived> matrixFunction(StemFunction f)
const;
455 const MatrixFunctionReturnValue<Derived> cosh()
const;
456 const MatrixFunctionReturnValue<Derived> sinh()
const;
457 const MatrixFunctionReturnValue<Derived> cos()
const;
458 const MatrixFunctionReturnValue<Derived> sin()
const;
459 const MatrixSquareRootReturnValue<Derived> sqrt()
const;
460 const MatrixLogarithmReturnValue<Derived> log()
const;
461 const MatrixPowerReturnValue<Derived> pow(
const RealScalar& p)
const;
462 const MatrixComplexPowerReturnValue<Derived> pow(
const std::complex<RealScalar>& p)
const;
465 EIGEN_DEVICE_FUNC MatrixBase() : Base() {}
468 EIGEN_DEVICE_FUNC
explicit MatrixBase(
int);
469 EIGEN_DEVICE_FUNC MatrixBase(
int,
int);
470 template<
typename OtherDerived> EIGEN_DEVICE_FUNC
explicit MatrixBase(
const MatrixBase<OtherDerived>&);
473 template<
typename OtherDerived> Derived&
operator+=(
const ArrayBase<OtherDerived>& )
474 {EIGEN_STATIC_ASSERT(std::ptrdiff_t(
sizeof(
typename OtherDerived::Scalar))==-1,YOU_CANNOT_MIX_ARRAYS_AND_MATRICES);
return *
this;}
476 template<
typename OtherDerived> Derived&
operator-=(
const ArrayBase<OtherDerived>& )
477 {EIGEN_STATIC_ASSERT(std::ptrdiff_t(
sizeof(
typename OtherDerived::Scalar))==-1,YOU_CANNOT_MIX_ARRAYS_AND_MATRICES);
return *
this;}
492 template<
typename Derived>
493 template<
typename OtherDerived>
497 other.
derived().applyThisOnTheRight(derived());
506 template<
typename Derived>
507 template<
typename OtherDerived>
510 other.
derived().applyThisOnTheRight(derived());
518 template<
typename Derived>
519 template<
typename OtherDerived>
522 other.
derived().applyThisOnTheLeft(derived());
527 #endif // EIGEN_MATRIXBASE_H
Generic expression of a matrix where all coefficients are defined by a functor.
Definition: CwiseNullaryOp.h:44
Robust Cholesky decomposition of a matrix with pivoting.
Definition: LDLT.h:48
Scalar determinant() const
Definition: Determinant.h:92
PlainObject unitOrthogonal(void) const
Definition: OrthoMethods.h:219
const LDLT< PlainObject > ldlt() const
Definition: LDLT.h:593
internal::traits< Derived >::Scalar Scalar
Definition: DenseBase.h:68
Expression of the product of two arbitrary matrices or vectors.
Definition: Product.h:107
bool operator!=(const MatrixBase< OtherDerived > &other) const
Definition: MatrixBase.h:303
Expression of a mathematical vector or matrix as an array object.
Definition: ArrayWrapper.h:41
Householder rank-revealing QR decomposition of a matrix with full pivoting.
Definition: ForwardDeclarations.h:254
HomogeneousReturnType homogeneous() const
Definition: Homogeneous.h:128
bool operator==(const MatrixBase< OtherDerived > &other) const
Definition: MatrixBase.h:295
void applyOnTheLeft(const EigenBase< OtherDerived > &other)
Definition: MatrixBase.h:520
Definition: Constants.h:265
Pseudo expression providing an operator = assuming no aliasing.
Definition: NoAlias.h:31
const Product< Derived, OtherDerived, LazyProduct > lazyProduct(const MatrixBase< OtherDerived > &other) const
Definition: GeneralProduct.h:434
Expression of the transpose of a matrix.
Definition: Transpose.h:53
static const BasisReturnType UnitW()
Definition: CwiseNullaryOp.h:866
RealScalar operatorNorm() const
Computes the L2 operator norm.
Definition: MatrixBaseEigenvalues.h:122
bool isLowerTriangular(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: TriangularMatrix.h:676
void applyOnTheRight(const EigenBase< OtherDerived > &other)
Definition: MatrixBase.h:508
LU decomposition of a matrix with partial pivoting, and related features.
Definition: ForwardDeclarations.h:248
Definition: Constants.h:320
const Derived & forceAlignedAccess() const
Definition: MatrixBase.h:310
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Definition: NumTraits.h:107
void makeHouseholder(EssentialPart &essential, Scalar &tau, RealScalar &beta) const
Definition: Householder.h:65
Derived & setIdentity()
Definition: CwiseNullaryOp.h:779
EigenvaluesReturnType eigenvalues() const
Computes the eigenvalues of a matrix.
Definition: MatrixBaseEigenvalues.h:67
Derived & derived()
Definition: EigenBase.h:44
const unsigned int RowMajorBit
Definition: Constants.h:61
Base class for all dense matrices, vectors, and arrays.
Definition: DenseBase.h:41
const ArrayWrapper< const Derived > array() const
Definition: MatrixBase.h:327
Matrix< Scalar, 3, 1 > eulerAngles(Index a0, Index a1, Index a2) const
Definition: EulerAngles.h:37
const PartialPivLU< PlainObject > partialPivLu() const
Definition: PartialPivLU.h:538
const HNormalizedReturnType hnormalized() const
Definition: Homogeneous.h:159
BDCSVD< PlainObject > bdcSvd(unsigned int computationOptions=0) const
Definition: BDCSVD.h:1200
internal::traits< Derived >::StorageIndex StorageIndex
The type used to store indices.
Definition: DenseBase.h:65
Definition: EigenBase.h:28
void applyHouseholderOnTheLeft(const EssentialPart &essential, const Scalar &tau, Scalar *workspace)
Definition: Householder.h:113
void normalize()
Definition: Dot.h:128
Definition: DenseBase.h:104
Expression of the inverse of another expression.
Definition: Inverse.h:43
const CwiseBinaryOp< std::equal_to< Scalar >, const Derived, const OtherDerived > cwiseEqual(const Eigen::MatrixBase< OtherDerived > &other) const
Definition: MatrixBase.h:44
RealScalar squaredNorm() const
Definition: Dot.h:88
static const BasisReturnType UnitX()
Definition: CwiseNullaryOp.h:836
Derived & operator+=(const MatrixBase< OtherDerived > &other)
Definition: CwiseBinaryOp.h:175
RealScalar blueNorm() const
Definition: StableNorm.h:200
ArrayWrapper< Derived > array()
Definition: MatrixBase.h:324
Householder rank-revealing QR decomposition of a matrix with column-pivoting.
Definition: ForwardDeclarations.h:253
Standard Cholesky decomposition (LL^T) of a matrix and associated features.
Definition: LLT.h:50
const ScalarMultipleReturnType operator*(const Scalar &scalar) const
Definition: MatrixBase.h:57
const PlainObject normalized() const
Definition: Dot.h:114
void computeInverseAndDetWithCheck(ResultType &inverse, typename ResultType::Scalar &determinant, bool &invertible, const RealScalar &absDeterminantThreshold=NumTraits< Scalar >::dummy_precision()) const
Definition: InverseImpl.h:358
PlainObject cross3(const MatrixBase< OtherDerived > &other) const
Definition: OrthoMethods.h:78
static const BasisReturnType UnitY()
Definition: CwiseNullaryOp.h:846
void adjointInPlace()
Definition: Transpose.h:309
RealScalar stableNorm() const
Definition: StableNorm.h:157
const ColPivHouseholderQR< PlainObject > colPivHouseholderQr() const
Definition: ColPivHouseholderQR.h:606
static const BasisReturnType Unit(Index size, Index i)
Definition: CwiseNullaryOp.h:808
Index diagonalSize() const
Definition: MatrixBase.h:103
bool isIdentity(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: CwiseNullaryOp.h:722
class Bidiagonal Divide and Conquer SVD
Definition: ForwardDeclarations.h:256
const CwiseBinaryOp< internal::scalar_product_op< typename Derived::Scalar, typename OtherDerived::Scalar >, const Derived, const OtherDerived > cwiseProduct(const Eigen::MatrixBase< OtherDerived > &other) const
Definition: MatrixBase.h:24
Definition: Constants.h:268
Definition: DenseBase.h:117
void computeInverseWithCheck(ResultType &inverse, bool &invertible, const RealScalar &absDeterminantThreshold=NumTraits< Scalar >::dummy_precision()) const
Definition: InverseImpl.h:397
internal::conditional< internal::is_same< typename internal::traits< Derived >::XprKind, MatrixXpr >::value, PlainMatrix, PlainArray >::type PlainObject
The plain matrix or array type corresponding to this expression.
Definition: DenseBase.h:209
bool isUnitary(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:227
static const IdentityReturnType Identity()
Definition: CwiseNullaryOp.h:705
Expression of a fixed-size or dynamic-size block.
Definition: Block.h:104
Derived & operator-=(const MatrixBase< OtherDerived > &other)
Definition: CwiseBinaryOp.h:162
RealScalar hypotNorm() const
Definition: StableNorm.h:212
LU decomposition of a matrix with complete pivoting, and related features.
Definition: ForwardDeclarations.h:247
void applyHouseholderOnTheRight(const EssentialPart &essential, const Scalar &tau, Scalar *workspace)
Definition: Householder.h:150
const FullPivLU< PlainObject > fullPivLu() const
Definition: FullPivLU.h:853
Householder QR decomposition of a matrix.
Definition: ForwardDeclarations.h:252
bool isUpperTriangular(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: TriangularMatrix.h:650
Two-sided Jacobi SVD decomposition of a rectangular matrix.
Definition: ForwardDeclarations.h:255
Scalar trace() const
Definition: Redux.h:487
Derived & operator=(const MatrixBase &other)
Definition: Assign.h:55
Derived & forceAlignedAccess()
Definition: MatrixBase.h:311
Derived & operator*=(const EigenBase< OtherDerived > &other)
Definition: MatrixBase.h:495
DiagonalReturnType diagonal()
Definition: Diagonal.h:188
NoAlias< Derived, Eigen::MatrixBase > noalias()
Definition: NoAlias.h:101
bool isOrthogonal(const MatrixBase< OtherDerived > &other, const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: Dot.h:208
Expression of a diagonal/subdiagonal/superdiagonal in a matrix.
Definition: Diagonal.h:63
Definition: DenseBase.h:110
const LLT< PlainObject > llt() const
Definition: LLT.h:473
void makeHouseholderInPlace(Scalar &tau, RealScalar &beta)
Definition: Householder.h:42
Generic expression where a coefficient-wise unary operator is applied to an expression.
Definition: CwiseUnaryOp.h:56
bool isDiagonal(const RealScalar &prec=NumTraits< Scalar >::dummy_precision()) const
Definition: DiagonalMatrix.h:292
const DiagonalWrapper< const Derived > asDiagonal() const
Definition: DiagonalMatrix.h:278
The matrix class, also used for vectors and row-vectors.
Definition: Matrix.h:178
const AdjointReturnType adjoint() const
Definition: Transpose.h:204
const Inverse< Derived > inverse() const
Definition: InverseImpl.h:331
Base class for all dense matrices, vectors, and expressions.
Definition: MatrixBase.h:48
const CwiseBinaryOp< std::not_equal_to< Scalar >, const Derived, const OtherDerived > cwiseNotEqual(const Eigen::MatrixBase< OtherDerived > &other) const
Definition: MatrixBase.h:64
RealScalar norm() const
Definition: Dot.h:100
static const BasisReturnType UnitZ()
Definition: CwiseNullaryOp.h:856
const HouseholderQR< PlainObject > householderQr() const
Definition: HouseholderQR.h:384
JacobiSVD< PlainObject > jacobiSvd(unsigned int computationOptions=0) const
Definition: JacobiSVD.h:795
const FullPivHouseholderQR< PlainObject > fullPivHouseholderQr() const
Definition: FullPivHouseholderQR.h:654
const PartialPivLU< PlainObject > lu() const
Definition: PartialPivLU.h:555