Public Types | Public Member Functions | Static Public Member Functions
MatrixBase< Derived > Class Template Reference

Base class for all dense matrices, vectors, and expressions. More...

+ Inheritance diagram for MatrixBase< Derived >:

List of all members.

Public Types

typedef Matrix< typename
internal::traits< Derived >
::Scalar, internal::traits
< Derived >::RowsAtCompileTime,
internal::traits< Derived >
::ColsAtCompileTime, AutoAlign|(internal::traits
< Derived >::Flags
&RowMajorBit?RowMajor:ColMajor),
internal::traits< Derived >
::MaxRowsAtCompileTime,
internal::traits< Derived >
::MaxColsAtCompileTime
PlainObject
 The plain matrix type corresponding to this expression.
- Public Types inherited from DenseBase< Derived >
enum  {
  RowsAtCompileTime,
  ColsAtCompileTime,
  SizeAtCompileTime,
  MaxRowsAtCompileTime,
  MaxColsAtCompileTime,
  MaxSizeAtCompileTime,
  IsVectorAtCompileTime,
  Flags,
  IsRowMajor ,
  CoeffReadCost
}
typedef internal::traits
< Derived >::Index 
Index
 The type of indices.

Public Member Functions

const AdjointReturnType adjoint () const
void adjointInPlace ()
template<typename OtherDerived >
void applyOnTheLeft (const EigenBase< OtherDerived > &other)
template<typename OtherScalar >
void applyOnTheLeft (Index p, Index q, const JacobiRotation< OtherScalar > &j)
template<typename OtherDerived >
void applyOnTheRight (const EigenBase< OtherDerived > &other)
template<typename OtherScalar >
void applyOnTheRight (Index p, Index q, const JacobiRotation< OtherScalar > &j)
ArrayWrapper< Derived > array ()
const DiagonalWrapper< const
Derived > 
asDiagonal () const
template<typename CustomBinaryOp , typename OtherDerived >
const CwiseBinaryOp
< CustomBinaryOp, const
Derived, const OtherDerived > 
binaryExpr (const Eigen::MatrixBase< OtherDerived > &other, const CustomBinaryOp &func=CustomBinaryOp()) const
RealScalar blueNorm () const
template<typename NewType >
internal::cast_return_type
< Derived, const CwiseUnaryOp
< internal::scalar_cast_op
< typename internal::traits
< Derived >::Scalar, NewType >
, const Derived > >::type 
cast () const
const ColPivHouseholderQR
< PlainObject
colPivHouseholderQr () const
template<typename ResultType >
void computeInverseAndDetWithCheck (ResultType &inverse, typename ResultType::Scalar &determinant, bool &invertible, const RealScalar &absDeterminantThreshold=NumTraits< Scalar >::dummy_precision()) const
template<typename ResultType >
void computeInverseWithCheck (ResultType &inverse, bool &invertible, const RealScalar &absDeterminantThreshold=NumTraits< Scalar >::dummy_precision()) const
ConjugateReturnType conjugate () const
template<typename OtherDerived >
cross_product_return_type
< OtherDerived >::type 
cross (const MatrixBase< OtherDerived > &other) const
template<typename OtherDerived >
PlainObject cross3 (const MatrixBase< OtherDerived > &other) const
const CwiseUnaryOp
< internal::scalar_abs_op
< Scalar >, const Derived > 
cwiseAbs () const
const CwiseUnaryOp
< internal::scalar_abs2_op
< Scalar >, const Derived > 
cwiseAbs2 () const
template<typename OtherDerived >
const CwiseBinaryOp (operator-)(const Eigen
template<typename OtherDerived >
const CwiseBinaryOp (operator+)(const Eigen
template<typename OtherDerived >
const CwiseBinaryOp
< std::equal_to< Scalar >
, const Derived, const
OtherDerived > 
cwiseEqual (const Eigen::MatrixBase< OtherDerived > &other) const
const CwiseUnaryOp
< std::binder1st
< std::equal_to< Scalar >
>, const Derived > 
cwiseEqual (const Scalar &s) const
const CwiseUnaryOp
< internal::scalar_inverse_op
< Scalar >, const Derived > 
cwiseInverse () const
template<typename OtherDerived >
const CwiseBinaryOp
< internal::scalar_max_op
< Scalar >, const Derived,
const OtherDerived > 
cwiseMax (const Eigen::MatrixBase< OtherDerived > &other) const
template<typename OtherDerived >
const CwiseBinaryOp
< internal::scalar_min_op
< Scalar >, const Derived,
const OtherDerived > 
cwiseMin (const Eigen::MatrixBase< OtherDerived > &other) const
template<typename OtherDerived >
const CwiseBinaryOp
< std::not_equal_to< Scalar >
, const Derived, const
OtherDerived > 
cwiseNotEqual (const Eigen::MatrixBase< OtherDerived > &other) const
template<typename OtherDerived >
const CwiseBinaryOp
< internal::scalar_product_op
< typename internal::traits
< Derived >::Scalar, typename
internal::traits< OtherDerived >
::Scalar >, const Derived,
const OtherDerived > 
cwiseProduct (const Eigen::MatrixBase< OtherDerived > &other) const
template<typename OtherDerived >
const CwiseBinaryOp
< internal::scalar_quotient_op
< Scalar >, const Derived,
const OtherDerived > 
cwiseQuotient (const Eigen::MatrixBase< OtherDerived > &other) const
const CwiseUnaryOp
< internal::scalar_sqrt_op
< Scalar >, const Derived > 
cwiseSqrt () const
Scalar determinant () const
DiagonalReturnType diagonal ()
const ConstDiagonalReturnType diagonal () const
DiagonalIndexReturnType
< Dynamic >::Type 
diagonal (Index index)
ConstDiagonalIndexReturnType
< Dynamic >::Type 
diagonal (Index index) const
Index diagonalSize () const
template<typename OtherDerived >
internal::scalar_product_traits
< typename internal::traits
< Derived >::Scalar, typename
internal::traits< OtherDerived >
::Scalar >::ReturnType 
dot (const MatrixBase< OtherDerived > &other) const
EigenvaluesReturnType eigenvalues () const
 Computes the eigenvalues of a matrix.
Matrix< Scalar, 3, 1 > eulerAngles (Index a0, Index a1, Index a2) const
const ForceAlignedAccess< Derived > forceAlignedAccess () const
ForceAlignedAccess< Derived > forceAlignedAccess ()
template<bool Enable>
internal::add_const_on_value_type
< typename
internal::conditional< Enable,
ForceAlignedAccess< Derived >
, Derived & >::type >::type 
forceAlignedAccessIf () const
template<bool Enable>
internal::conditional< Enable,
ForceAlignedAccess< Derived >
, Derived & >::type 
forceAlignedAccessIf ()
const FullPivHouseholderQR
< PlainObject
fullPivHouseholderQr () const
const FullPivLU< PlainObjectfullPivLu () const
const HNormalizedReturnType hnormalized () const
HomogeneousReturnType homogeneous () const
const HouseholderQR< PlainObjecthouseholderQr () const
RealScalar hypotNorm () const
const ImagReturnType imag () const
NonConstImagReturnType imag ()
const internal::inverse_impl
< Derived > 
inverse () const
bool isDiagonal (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isIdentity (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isLowerTriangular (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
template<typename OtherDerived >
bool isOrthogonal (const MatrixBase< OtherDerived > &other, RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isUnitary (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isUpperTriangular (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
JacobiSVD< PlainObjectjacobiSvd (unsigned int computationOptions=0) const
template<typename OtherDerived >
const LazyProductReturnType
< Derived, OtherDerived >
::Type 
lazyProduct (const MatrixBase< OtherDerived > &other) const
const LDLT< PlainObjectldlt () const
const LLT< PlainObjectllt () const
template<int p>
RealScalar lpNorm () const
const PartialPivLU< PlainObjectlu () const
template<typename EssentialPart >
void makeHouseholder (EssentialPart &essential, Scalar &tau, RealScalar &beta) const
NoAlias< Derived,
Eigen::MatrixBase
noalias ()
RealScalar norm () const
void normalize ()
const PlainObject normalized () const
template<typename OtherDerived >
bool operator!= (const MatrixBase< OtherDerived > &other) const
const ScalarMultipleReturnType operator* (const Scalar &scalar) const
const CwiseUnaryOp
< internal::scalar_multiple2_op
< Scalar, std::complex< Scalar >
>, const Derived > 
operator* (const std::complex< Scalar > &scalar) const
template<typename OtherDerived >
const ProductReturnType
< Derived, OtherDerived >
::Type 
operator* (const MatrixBase< OtherDerived > &other) const
template<typename DiagonalDerived >
const DiagonalProduct< Derived,
DiagonalDerived, OnTheRight
operator* (const DiagonalBase< DiagonalDerived > &diagonal) const
ScalarMultipleReturnType operator* (const UniformScaling< Scalar > &s) const
template<typename OtherDerived >
Derived & operator*= (const EigenBase< OtherDerived > &other)
template<typename OtherDerived >
Derived & operator+= (const MatrixBase< OtherDerived > &other)
const CwiseUnaryOp
< internal::scalar_opposite_op
< typename internal::traits
< Derived >::Scalar >, const
Derived > 
operator- () const
template<typename OtherDerived >
Derived & operator-= (const MatrixBase< OtherDerived > &other)
const CwiseUnaryOp
< internal::scalar_quotient1_op
< typename internal::traits
< Derived >::Scalar >, const
Derived > 
operator/ (const Scalar &scalar) const
Derived & operator= (const MatrixBase &other)
template<typename OtherDerived >
Derived & operator= (const DenseBase< OtherDerived > &other)
template<typename OtherDerived >
Derived & operator= (const EigenBase< OtherDerived > &other)
 Copies the generic expression other into *this.
template<typename OtherDerived >
bool operator== (const MatrixBase< OtherDerived > &other) const
RealScalar operatorNorm () const
 Computes the L2 operator norm.
const PartialPivLU< PlainObjectpartialPivLu () const
RealReturnType real () const
NonConstRealReturnType real ()
Derived & setIdentity ()
Derived & setIdentity (Index rows, Index cols)
 Resizes to the given size, and writes the identity expression (not necessarily square) into *this.
RealScalar squaredNorm () const
RealScalar stableNorm () const
Scalar trace () const
template<unsigned int Mode>
TriangularViewReturnType< Mode >
::Type 
triangularView ()
template<unsigned int Mode>
ConstTriangularViewReturnType
< Mode >::Type 
triangularView () const
template<typename CustomUnaryOp >
const CwiseUnaryOp
< CustomUnaryOp, const Derived > 
unaryExpr (const CustomUnaryOp &func=CustomUnaryOp()) const
 Apply a unary operator coefficient-wise.
template<typename CustomViewOp >
const CwiseUnaryView
< CustomViewOp, const Derived > 
unaryViewExpr (const CustomViewOp &func=CustomViewOp()) const
PlainObject unitOrthogonal (void) const
- Public Member Functions inherited from DenseBase< Derived >
bool all (void) const
bool any (void) const
Block< Derived > block (Index startRow, Index startCol, Index blockRows, Index blockCols)
const Block< const Derived > block (Index startRow, Index startCol, Index blockRows, Index blockCols) const
template<int BlockRows, int BlockCols>
Block< Derived, BlockRows,
BlockCols > 
block (Index startRow, Index startCol)
template<int BlockRows, int BlockCols>
const Block< const Derived,
BlockRows, BlockCols > 
block (Index startRow, Index startCol) const
Block< Derived > bottomLeftCorner (Index cRows, Index cCols)
const Block< const Derived > bottomLeftCorner (Index cRows, Index cCols) const
template<int CRows, int CCols>
Block< Derived, CRows, CCols > bottomLeftCorner ()
template<int CRows, int CCols>
const Block< const Derived,
CRows, CCols > 
bottomLeftCorner () const
Block< Derived > bottomRightCorner (Index cRows, Index cCols)
const Block< const Derived > bottomRightCorner (Index cRows, Index cCols) const
template<int CRows, int CCols>
Block< Derived, CRows, CCols > bottomRightCorner ()
template<int CRows, int CCols>
const Block< const Derived,
CRows, CCols > 
bottomRightCorner () const
RowsBlockXpr bottomRows (Index n)
ConstRowsBlockXpr bottomRows (Index n) const
template<int N>
NRowsBlockXpr< N >::Type bottomRows ()
template<int N>
ConstNRowsBlockXpr< N >::Type bottomRows () const
ColXpr col (Index i)
ConstColXpr col (Index i) const
ConstColwiseReturnType colwise () const
ColwiseReturnType colwise ()
Index count () const
const internal::eval< Derived >
::type 
eval () const
void fill (const Scalar &value)
template<unsigned int Added, unsigned int Removed>
const Flagged< Derived, Added,
Removed > 
flagged () const
const WithFormat< Derived > format (const IOFormat &fmt) const
SegmentReturnType head (Index size)
DenseBase::ConstSegmentReturnType head (Index size) const
template<int Size>
FixedSegmentReturnType< Size >
::Type 
head ()
template<int Size>
ConstFixedSegmentReturnType
< Size >::Type 
head () const
Index innerSize () const
template<typename OtherDerived >
bool isApprox (const DenseBase< OtherDerived > &other, RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isApproxToConstant (const Scalar &value, RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isConstant (const Scalar &value, RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
template<typename Derived >
bool isMuchSmallerThan (const typename NumTraits< Scalar >::Real &other, RealScalar prec) const
template<typename OtherDerived >
bool isMuchSmallerThan (const DenseBase< OtherDerived > &other, RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isOnes (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
bool isZero (RealScalar prec=NumTraits< Scalar >::dummy_precision()) const
ColsBlockXpr leftCols (Index n)
ConstColsBlockXpr leftCols (Index n) const
template<int N>
NColsBlockXpr< N >::Type leftCols ()
template<int N>
ConstNColsBlockXpr< N >::Type leftCols () const
internal::traits< Derived >::Scalar maxCoeff () const
template<typename IndexType >
internal::traits< Derived >::Scalar maxCoeff (IndexType *row, IndexType *col) const
template<typename IndexType >
internal::traits< Derived >::Scalar maxCoeff (IndexType *index) const
Scalar mean () const
ColsBlockXpr middleCols (Index startCol, Index numCols)
ConstColsBlockXpr middleCols (Index startCol, Index numCols) const
template<int N>
NColsBlockXpr< N >::Type middleCols (Index startCol)
template<int N>
ConstNColsBlockXpr< N >::Type middleCols (Index startCol) const
RowsBlockXpr middleRows (Index startRow, Index numRows)
ConstRowsBlockXpr middleRows (Index startRow, Index numRows) const
template<int N>
NRowsBlockXpr< N >::Type middleRows (Index startRow)
template<int N>
ConstNRowsBlockXpr< N >::Type middleRows (Index startRow) const
internal::traits< Derived >::Scalar minCoeff () const
template<typename IndexType >
internal::traits< Derived >::Scalar minCoeff (IndexType *row, IndexType *col) const
template<typename IndexType >
internal::traits< Derived >::Scalar minCoeff (IndexType *index) const
const NestByValue< Derived > nestByValue () const
Index nonZeros () const
CommaInitializer< Derived > operator<< (const Scalar &s)
template<typename OtherDerived >
CommaInitializer< Derived > operator<< (const DenseBase< OtherDerived > &other)
Derived & operator= (const DenseBase &other)
Index outerSize () const
Scalar prod () const
template<typename Func >
internal::result_of< Func(typename
internal::traits< Derived >
::Scalar)>::type 
redux (const Func &func) const
template<int RowFactor, int ColFactor>
const Replicate< Derived,
RowFactor, ColFactor > 
replicate () const
const Replicate< Derived,
Dynamic, Dynamic
replicate (Index rowFacor, Index colFactor) const
void resize (Index size)
void resize (Index rows, Index cols)
ReverseReturnType reverse ()
ConstReverseReturnType reverse () const
void reverseInPlace ()
ColsBlockXpr rightCols (Index n)
ConstColsBlockXpr rightCols (Index n) const
template<int N>
NColsBlockXpr< N >::Type rightCols ()
template<int N>
ConstNColsBlockXpr< N >::Type rightCols () const
RowXpr row (Index i)
ConstRowXpr row (Index i) const
ConstRowwiseReturnType rowwise () const
RowwiseReturnType rowwise ()
SegmentReturnType segment (Index start, Index size)
DenseBase::ConstSegmentReturnType segment (Index start, Index size) const
template<int Size>
FixedSegmentReturnType< Size >
::Type 
segment (Index start)
template<int Size>
ConstFixedSegmentReturnType
< Size >::Type 
segment (Index start) const
template<typename ThenDerived , typename ElseDerived >
const Select< Derived,
ThenDerived, ElseDerived > 
select (const DenseBase< ThenDerived > &thenMatrix, const DenseBase< ElseDerived > &elseMatrix) const
template<typename ThenDerived >
const Select< Derived,
ThenDerived, typename
ThenDerived::ConstantReturnType > 
select (const DenseBase< ThenDerived > &thenMatrix, typename ThenDerived::Scalar elseScalar) const
template<typename ElseDerived >
const Select< Derived,
typename
ElseDerived::ConstantReturnType,
ElseDerived > 
select (typename ElseDerived::Scalar thenScalar, const DenseBase< ElseDerived > &elseMatrix) const
Derived & setConstant (const Scalar &value)
Derived & setLinSpaced (Index size, const Scalar &low, const Scalar &high)
 Sets a linearly space vector.
Derived & setOnes ()
Derived & setRandom ()
Derived & setZero ()
Scalar sum () const
template<typename OtherDerived >
void swap (const DenseBase< OtherDerived > &other, int=OtherDerived::ThisConstantIsPrivateInPlainObjectBase)
template<typename OtherDerived >
void swap (PlainObjectBase< OtherDerived > &other)
SegmentReturnType tail (Index size)
DenseBase::ConstSegmentReturnType tail (Index size) const
template<int Size>
FixedSegmentReturnType< Size >
::Type 
tail ()
template<int Size>
ConstFixedSegmentReturnType
< Size >::Type 
tail () const
Block< Derived > topLeftCorner (Index cRows, Index cCols)
const Block< const Derived > topLeftCorner (Index cRows, Index cCols) const
template<int CRows, int CCols>
Block< Derived, CRows, CCols > topLeftCorner ()
template<int CRows, int CCols>
const Block< const Derived,
CRows, CCols > 
topLeftCorner () const
Block< Derived > topRightCorner (Index cRows, Index cCols)
const Block< const Derived > topRightCorner (Index cRows, Index cCols) const
template<int CRows, int CCols>
Block< Derived, CRows, CCols > topRightCorner ()
template<int CRows, int CCols>
const Block< const Derived,
CRows, CCols > 
topRightCorner () const
RowsBlockXpr topRows (Index n)
ConstRowsBlockXpr topRows (Index n) const
template<int N>
NRowsBlockXpr< N >::Type topRows ()
template<int N>
ConstNRowsBlockXpr< N >::Type topRows () const
Eigen::Transpose< Derived > transpose ()
ConstTransposeReturnType transpose () const
void transposeInPlace ()
CoeffReturnType value () const
template<typename Visitor >
void visit (Visitor &func) const

Static Public Member Functions

static const IdentityReturnType Identity ()
static const IdentityReturnType Identity (Index rows, Index cols)
static const BasisReturnType Unit (Index size, Index i)
static const BasisReturnType Unit (Index i)
static const BasisReturnType UnitW ()
static const BasisReturnType UnitX ()
static const BasisReturnType UnitY ()
static const BasisReturnType UnitZ ()
- Static Public Member Functions inherited from DenseBase< Derived >
static const ConstantReturnType Constant (Index rows, Index cols, const Scalar &value)
static const ConstantReturnType Constant (Index size, const Scalar &value)
static const ConstantReturnType Constant (const Scalar &value)
static const
SequentialLinSpacedReturnType 
LinSpaced (Sequential_t, Index size, const Scalar &low, const Scalar &high)
 Sets a linearly space vector.
static const
RandomAccessLinSpacedReturnType 
LinSpaced (Index size, const Scalar &low, const Scalar &high)
 Sets a linearly space vector.
static const
SequentialLinSpacedReturnType 
LinSpaced (Sequential_t, const Scalar &low, const Scalar &high)
 Sets a linearly space vector.
static const
RandomAccessLinSpacedReturnType 
LinSpaced (const Scalar &low, const Scalar &high)
 Sets a linearly space vector.
template<typename CustomNullaryOp >
static const CwiseNullaryOp
< CustomNullaryOp, Derived > 
NullaryExpr (Index rows, Index cols, const CustomNullaryOp &func)
template<typename CustomNullaryOp >
static const CwiseNullaryOp
< CustomNullaryOp, Derived > 
NullaryExpr (Index size, const CustomNullaryOp &func)
template<typename CustomNullaryOp >
static const CwiseNullaryOp
< CustomNullaryOp, Derived > 
NullaryExpr (const CustomNullaryOp &func)
static const ConstantReturnType Ones (Index rows, Index cols)
static const ConstantReturnType Ones (Index size)
static const ConstantReturnType Ones ()
static const CwiseNullaryOp
< internal::scalar_random_op
< Scalar >, Derived > 
Random (Index rows, Index cols)
static const CwiseNullaryOp
< internal::scalar_random_op
< Scalar >, Derived > 
Random (Index size)
static const CwiseNullaryOp
< internal::scalar_random_op
< Scalar >, Derived > 
Random ()
static const ConstantReturnType Zero (Index rows, Index cols)
static const ConstantReturnType Zero (Index size)
static const ConstantReturnType Zero ()

Additional Inherited Members

- Protected Member Functions inherited from DenseBase< Derived >
 DenseBase ()
std::ostream & operator<< (std::ostream &s, const DenseBase< Derived > &m)

Detailed Description

template<typename Derived>
class Eigen::MatrixBase< Derived >

Base class for all dense matrices, vectors, and expressions.

This class is the base that is inherited by all matrix, vector, and related expression types. Most of the Eigen API is contained in this class, and its base classes. Other important classes for the Eigen API are Matrix, and VectorwiseOp.

Note that some methods are defined in other modules such as the LU module LU module for all functions related to matrix inversions.

Template Parameters:
Derivedis the derived type, e.g. a matrix type, or an expression, etc.

When writing a function taking Eigen objects as argument, if you want your function to take as argument any matrix, vector, or expression, just let it take a MatrixBase argument. As an example, here is a function printFirstRow which, given a matrix, vector, or expression x, prints the first row of x.

template<typename Derived>
void printFirstRow(const Eigen::MatrixBase<Derived>& x)
{
cout << x.row(0) << endl;
}

This class can be extended with the help of the plugin mechanism described on the page Customizing/Extending Eigen by defining the preprocessor symbol EIGEN_MATRIXBASE_PLUGIN.

See also:
The class hierarchy

Member Typedef Documentation

typedef Matrix<typename internal::traits<Derived>::Scalar, internal::traits<Derived>::RowsAtCompileTime, internal::traits<Derived>::ColsAtCompileTime, AutoAlign | (internal::traits<Derived>::Flags&RowMajorBit ? RowMajor : ColMajor), internal::traits<Derived>::MaxRowsAtCompileTime, internal::traits<Derived>::MaxColsAtCompileTime > PlainObject

The plain matrix type corresponding to this expression.

This is not necessarily exactly the return type of eval(). In the case of plain matrices, the return type of eval() is a const reference to a matrix, not a matrix! It is however guaranteed that the return type of eval() is either PlainObject or const PlainObject&.


Member Function Documentation

const MatrixBase< Derived >::AdjointReturnType adjoint ( ) const
inline
Returns:
an expression of the adjoint (i.e. conjugate transpose) of *this.

Example:

cout << "Here is the 2x2 complex matrix m:" << endl << m << endl;
cout << "Here is the adjoint of m:" << endl << m.adjoint() << endl;

Output:

Here is the 2x2 complex matrix m:
 (-0.211,0.68) (-0.605,0.823)
 (0.597,0.566)  (0.536,-0.33)
Here is the adjoint of m:
(-0.211,-0.68) (0.597,-0.566)
(-0.605,-0.823)   (0.536,0.33)
Warning:
If you want to replace a matrix by its own adjoint, do NOT do this:
m = m.adjoint(); // bug!!! caused by aliasing effect
Instead, use the adjointInPlace() method:
m.adjointInPlace();
which gives Eigen good opportunities for optimization, or alternatively you can also do:
m = m.adjoint().eval();
See also:
adjointInPlace(), transpose(), conjugate(), class Transpose, class internal::scalar_conjugate_op
void adjointInPlace ( )
inline

This is the "in place" version of adjoint(): it replaces *this by its own transpose. Thus, doing

m.adjointInPlace();

has the same effect on m as doing

m = m.adjoint().eval();

and is faster and also safer because in the latter line of code, forgetting the eval() results in a bug caused by aliasing.

Notice however that this method is only useful if you want to replace a matrix by its own adjoint. If you just need the adjoint of a matrix, use adjoint().

Note:
if the matrix is not square, then *this must be a resizable matrix.
See also:
transpose(), adjoint(), transposeInPlace()
void applyOnTheLeft ( const EigenBase< OtherDerived > &  other)
inline

replaces *this by *this * other.

void applyOnTheLeft ( Index  p,
Index  q,
const JacobiRotation< OtherScalar > &  j 
)
inline

This is defined in the Jacobi module.

#include <Eigen/Jacobi>

Applies the rotation in the plane j to the rows p and q of *this, i.e., it computes B = J * B, with $ B = \left ( \begin{array}{cc} \text{*this.row}(p) \\ \text{*this.row}(q) \end{array} \right ) $.

See also:
class JacobiRotation, MatrixBase::applyOnTheRight(), internal::apply_rotation_in_the_plane()
void applyOnTheRight ( const EigenBase< OtherDerived > &  other)
inline

replaces this by *this * other. It is equivalent to MatrixBase::operator</em>=()

ArrayWrapper<Derived> array ( )
inline
Returns:
an Array expression of this matrix
See also:
ArrayBase::matrix()
const DiagonalWrapper< const Derived > asDiagonal ( ) const
inline
Returns:
a pseudo-expression of a diagonal matrix with *this as vector of diagonal coefficients

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

Example:

cout << Matrix3i(Vector3i(2,5,6).asDiagonal()) << endl;

Output:

2 0 0
0 5 0
0 0 6
See also:
class DiagonalWrapper, class DiagonalMatrix, diagonal(), isDiagonal()
const CwiseBinaryOp<CustomBinaryOp, const Derived, const OtherDerived> binaryExpr ( const Eigen::MatrixBase< OtherDerived > &  other,
const CustomBinaryOp &  func = CustomBinaryOp() 
) const
inline
Returns:
an expression of a custom coefficient-wise operator func of *this and other

The template parameter CustomBinaryOp is the type of the functor of the custom operator (see class CwiseBinaryOp for an example)

Here is an example illustrating the use of custom functors:

#include <Eigen/Core>
#include <iostream>
using namespace Eigen;
using namespace std;
// define a custom template binary functor
template<typename Scalar> struct MakeComplexOp {
EIGEN_EMPTY_STRUCT_CTOR(MakeComplexOp)
typedef complex<Scalar> result_type;
complex<Scalar> operator()(const Scalar& a, const Scalar& b) const { return complex<Scalar>(a,b); }
};
int main(int, char**)
{
cout << m1.binaryExpr(m2, MakeComplexOp<double>()) << endl;
return 0;
}

Output:

   (0.68,0.271)  (0.823,-0.967) (-0.444,-0.687)   (-0.27,0.998)
 (-0.211,0.435) (-0.605,-0.514)  (0.108,-0.198) (0.0268,-0.563)
 (0.566,-0.717)  (-0.33,-0.726) (-0.0452,-0.74)  (0.904,0.0259)
  (0.597,0.214)   (0.536,0.608)  (0.258,-0.782)   (0.832,0.678)
See also:
class CwiseBinaryOp, operator+(), operator-(), cwiseProduct()
NumTraits< typename internal::traits< Derived >::Scalar >::Real blueNorm ( ) const
inline
Returns:
the l2 norm of *this using the Blue's algorithm. A Portable Fortran Program to Find the Euclidean Norm of a Vector, ACM TOMS, Vol 4, Issue 1, 1978.

For architecture/scalar types without vectorization, this version is much faster than stableNorm(). Otherwise the stableNorm() is faster.

See also:
norm(), stableNorm(), hypotNorm()
internal::cast_return_type<Derived,const CwiseUnaryOp<internal::scalar_cast_op<typename internal::traits<Derived>::Scalar, NewType>, const Derived> >::type cast ( ) const
inline
Returns:
an expression of *this with the Scalar type casted to NewScalar.

The template parameter NewScalar is the type we are casting the scalars to.

See also:
class CwiseUnaryOp
const ColPivHouseholderQR< typename MatrixBase< Derived >::PlainObject > colPivHouseholderQr ( ) const
Returns:
the column-pivoting Householder QR decomposition of *this.
See also:
class ColPivHouseholderQR
void computeInverseAndDetWithCheck ( ResultType &  inverse,
typename ResultType::Scalar &  determinant,
bool &  invertible,
const RealScalar &  absDeterminantThreshold = NumTraits<Scalar>::dummy_precision() 
) const
inline

This is defined in the LU module.

#include <Eigen/LU>

Computation of matrix inverse and determinant, with invertibility check.

This is only for fixed-size square matrices of size up to 4x4.

Parameters:
inverseReference to the matrix in which to store the inverse.
determinantReference to the variable in which to store the inverse.
invertibleReference to the bool variable in which to store whether the matrix is invertible.
absDeterminantThresholdOptional parameter controlling the invertibility check. The matrix will be declared invertible if the absolute value of its determinant is greater than this threshold.

Example:

cout << "Here is the matrix m:" << endl << m << endl;
bool invertible;
double determinant;
m.computeInverseAndDetWithCheck(inverse,determinant,invertible);
cout << "Its determinant is " << determinant << endl;
if(invertible) {
cout << "It is invertible, and its inverse is:" << endl << inverse << endl;
}
else {
cout << "It is not invertible." << endl;
}

Output:

Here is the matrix m:
  0.68  0.597  -0.33
-0.211  0.823  0.536
 0.566 -0.605 -0.444
Its determinant is 0.209
It is invertible, and its inverse is:
-0.199   2.23   2.84
  1.01 -0.555  -1.42
 -1.62   3.59   3.29
See also:
inverse(), computeInverseWithCheck()
void computeInverseWithCheck ( ResultType &  inverse,
bool &  invertible,
const RealScalar &  absDeterminantThreshold = NumTraits<Scalar>::dummy_precision() 
) const
inline

This is defined in the LU module.

#include <Eigen/LU>

Computation of matrix inverse, with invertibility check.

This is only for fixed-size square matrices of size up to 4x4.

Parameters:
inverseReference to the matrix in which to store the inverse.
invertibleReference to the bool variable in which to store whether the matrix is invertible.
absDeterminantThresholdOptional parameter controlling the invertibility check. The matrix will be declared invertible if the absolute value of its determinant is greater than this threshold.

Example:

cout << "Here is the matrix m:" << endl << m << endl;
bool invertible;
m.computeInverseWithCheck(inverse,invertible);
if(invertible) {
cout << "It is invertible, and its inverse is:" << endl << inverse << endl;
}
else {
cout << "It is not invertible." << endl;
}

Output:

Here is the matrix m:
  0.68  0.597  -0.33
-0.211  0.823  0.536
 0.566 -0.605 -0.444
It is invertible, and its inverse is:
-0.199   2.23   2.84
  1.01 -0.555  -1.42
 -1.62   3.59   3.29
See also:
inverse(), computeInverseAndDetWithCheck()
ConjugateReturnType conjugate ( ) const
inline
Returns:
an expression of the complex conjugate of *this.
See also:
adjoint()
MatrixBase< Derived >::template cross_product_return_type< OtherDerived >::type cross ( const MatrixBase< OtherDerived > &  other) const
inline

This is defined in the Geometry module.

#include <Eigen/Geometry>
Returns:
the cross product of *this and other

Here is a very good explanation of cross-product: http://xkcd.com/199/

See also:
MatrixBase::cross3()
MatrixBase< Derived >::PlainObject cross3 ( const MatrixBase< OtherDerived > &  other) const
inline

This is defined in the Geometry module.

#include <Eigen/Geometry>
Returns:
the cross product of *this and other using only the x, y, and z coefficients

The size of *this and other must be four. This function is especially useful when using 4D vectors instead of 3D ones to get advantage of SSE/AltiVec vectorization.

See also:
MatrixBase::cross()
const CwiseUnaryOp<internal::scalar_abs_op<Scalar>, const Derived> cwiseAbs ( ) const
inline
Returns:
an expression of the coefficient-wise absolute value of *this

Example:

MatrixXd m(2,3);
m << 2, -4, 6,
-5, 1, 0;
cout << m.cwiseAbs() << endl;

Output:

2 4 6
5 1 0
See also:
cwiseAbs2()
const CwiseUnaryOp<internal::scalar_abs2_op<Scalar>, const Derived> cwiseAbs2 ( ) const
inline
Returns:
an expression of the coefficient-wise squared absolute value of *this

Example:

MatrixXd m(2,3);
m << 2, -4, 6,
-5, 1, 0;
cout << m.cwiseAbs2() << endl;

Output:

 4 16 36
25  1  0
See also:
cwiseAbs()
const CwiseBinaryOp ( operator-  ) const
inline
Returns:
an expression of the difference of *this and other
Note:
If you want to substract a given scalar from all coefficients, see Cwise::operator-().
See also:
class CwiseBinaryOp, operator-=()
const CwiseBinaryOp ( operator+  ) const
inline
Returns:
an expression of the sum of *this and other
Note:
If you want to add a given scalar to all coefficients, see Cwise::operator+().
See also:
class CwiseBinaryOp, operator+=()
const CwiseBinaryOp<std::equal_to<Scalar>, const Derived, const OtherDerived> cwiseEqual ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the coefficient-wise == operator of *this and other
Warning:
this performs an exact comparison, which is generally a bad idea with floating-point types. In order to check for equality between two vectors or matrices with floating-point coefficients, it is generally a far better idea to use a fuzzy comparison as provided by isApprox() and isMuchSmallerThan().

Example:

MatrixXi m(2,2);
m << 1, 0,
1, 1;
cout << "Comparing m with identity matrix:" << endl;
cout << m.cwiseEqual(MatrixXi::Identity(2,2)) << endl;
int count = m.cwiseEqual(MatrixXi::Identity(2,2)).count();
cout << "Number of coefficients that are equal: " << count << endl;

Output:

Comparing m with identity matrix:
1 1
0 1
Number of coefficients that are equal: 3
See also:
cwiseNotEqual(), isApprox(), isMuchSmallerThan()
const CwiseUnaryOp<std::binder1st<std::equal_to<Scalar> >, const Derived> cwiseEqual ( const Scalar &  s) const
inline
Returns:
an expression of the coefficient-wise == operator of *this and a scalar s
Warning:
this performs an exact comparison, which is generally a bad idea with floating-point types. In order to check for equality between two vectors or matrices with floating-point coefficients, it is generally a far better idea to use a fuzzy comparison as provided by isApprox() and isMuchSmallerThan().
See also:
cwiseEqual(const MatrixBase<OtherDerived> &) const
const CwiseUnaryOp<internal::scalar_inverse_op<Scalar>, const Derived> cwiseInverse ( ) const
inline
Returns:
an expression of the coefficient-wise inverse of *this.

Example:

MatrixXd m(2,3);
m << 2, 0.5, 1,
3, 0.25, 1;
cout << m.cwiseInverse() << endl;

Output:

0.5 2 1
0.333 4 1
See also:
cwiseProduct()
const CwiseBinaryOp<internal::scalar_max_op<Scalar>, const Derived, const OtherDerived> cwiseMax ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the coefficient-wise max of *this and other

Example:

Vector3d v(2,3,4), w(4,2,3);
cout << v.cwiseMax(w) << endl;

Output:

4
3
4
See also:
class CwiseBinaryOp, min()
const CwiseBinaryOp<internal::scalar_min_op<Scalar>, const Derived, const OtherDerived> cwiseMin ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the coefficient-wise min of *this and other

Example:

Vector3d v(2,3,4), w(4,2,3);
cout << v.cwiseMin(w) << endl;

Output:

2
2
3
See also:
class CwiseBinaryOp, max()
const CwiseBinaryOp<std::not_equal_to<Scalar>, const Derived, const OtherDerived> cwiseNotEqual ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the coefficient-wise != operator of *this and other
Warning:
this performs an exact comparison, which is generally a bad idea with floating-point types. In order to check for equality between two vectors or matrices with floating-point coefficients, it is generally a far better idea to use a fuzzy comparison as provided by isApprox() and isMuchSmallerThan().

Example:

MatrixXi m(2,2);
m << 1, 0,
1, 1;
cout << "Comparing m with identity matrix:" << endl;
cout << m.cwiseNotEqual(MatrixXi::Identity(2,2)) << endl;
int count = m.cwiseNotEqual(MatrixXi::Identity(2,2)).count();
cout << "Number of coefficients that are not equal: " << count << endl;

Output:

Comparing m with identity matrix:
0 0
1 0
Number of coefficients that are not equal: 1
See also:
cwiseEqual(), isApprox(), isMuchSmallerThan()
const CwiseBinaryOp< internal::scalar_product_op< typename internal::traits< Derived >::Scalar, typename internal::traits< OtherDerived >::Scalar >, const Derived , const OtherDerived > cwiseProduct ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the Schur product (coefficient wise product) of *this and other

Example:

Matrix3i c = a.cwiseProduct(b);
cout << "a:\n" << a << "\nb:\n" << b << "\nc:\n" << c << endl;

Output:

a:
 7  6 -3
-2  9  6
 6 -6 -5
b:
 1 -3  9
 0  0  3
 3  9  5
c:
  7 -18 -27
  0   0  18
 18 -54 -25
See also:
class CwiseBinaryOp, cwiseAbs2
const CwiseBinaryOp<internal::scalar_quotient_op<Scalar>, const Derived, const OtherDerived> cwiseQuotient ( const Eigen::MatrixBase< OtherDerived > &  other) const
inline
Returns:
an expression of the coefficient-wise quotient of *this and other

Example:

Vector3d v(2,3,4), w(4,2,3);
cout << v.cwiseQuotient(w) << endl;

Output:

0.5
1.5
1.33
See also:
class CwiseBinaryOp, cwiseProduct(), cwiseInverse()
const CwiseUnaryOp<internal::scalar_sqrt_op<Scalar>, const Derived> cwiseSqrt ( ) const
inline
Returns:
an expression of the coefficient-wise square root of *this.

Example:

Vector3d v(1,2,4);
cout << v.cwiseSqrt() << endl;

Output:

1
1.41
2
See also:
cwisePow(), cwiseSquare()
internal::traits< Derived >::Scalar determinant ( ) const
inline

This is defined in the LU module.

#include <Eigen/LU>
Returns:
the determinant of this matrix
MatrixBase< Derived >::template DiagonalIndexReturnType< Index >::Type diagonal ( )
inline
Returns:
an expression of the main diagonal of the matrix *this

*this is not required to be square.

Example:

cout << "Here is the matrix m:" << endl << m << endl;
cout << "Here are the coefficients on the main diagonal of m:" << endl
<< m.diagonal() << endl;

Output:

Here is the matrix m:
 7  6 -3
-2  9  6
 6 -6 -5
Here are the coefficients on the main diagonal of m:
7
9
-5
See also:
class Diagonal
Returns:
an expression of the DiagIndex-th sub or super diagonal of the matrix *this

*this is not required to be square.

The template parameter DiagIndex represent a super diagonal if DiagIndex > 0 and a sub diagonal otherwise. DiagIndex == 0 is equivalent to the main diagonal.

Example:

cout << "Here is the matrix m:" << endl << m << endl;
cout << "Here are the coefficients on the 1st super-diagonal and 2nd sub-diagonal of m:" << endl
<< m.diagonal<1>().transpose() << endl
<< m.diagonal<-2>().transpose() << endl;

Output:

Here is the matrix m:
 7  9 -5 -3
-2 -6  1  0
 6 -3  0  9
 6  6  3  9
Here are the coefficients on the 1st super-diagonal and 2nd sub-diagonal of m:
9 1 9
6 6
See also:
MatrixBase::diagonal(), class Diagonal
MatrixBase< Derived >::template ConstDiagonalIndexReturnType< Index >::Type diagonal ( ) const
inline

This is the const version of diagonal().

This is the const version of diagonal<int>().

MatrixBase< Derived >::template DiagonalIndexReturnType< Dynamic >::Type diagonal ( Index  index)
inline
Returns:
an expression of the DiagIndex-th sub or super diagonal of the matrix *this

*this is not required to be square.

The template parameter DiagIndex represent a super diagonal if DiagIndex > 0 and a sub diagonal otherwise. DiagIndex == 0 is equivalent to the main diagonal.

Example:

cout << "Here is the matrix m:" << endl << m << endl;
cout << "Here are the coefficients on the 1st super-diagonal and 2nd sub-diagonal of m:" << endl
<< m.diagonal(1).transpose() << endl
<< m.diagonal(-2).transpose() << endl;

Output:

Here is the matrix m:
 7  9 -5 -3
-2 -6  1  0
 6 -3  0  9
 6  6  3  9
Here are the coefficients on the 1st super-diagonal and 2nd sub-diagonal of m:
9 1 9
6 6
See also:
MatrixBase::diagonal(), class Diagonal
MatrixBase< Derived >::template ConstDiagonalIndexReturnType< Dynamic >::Type diagonal ( Index  index) const
inline

This is the const version of diagonal(Index).

Index diagonalSize ( ) const
inline
Returns:
the size of the main diagonal, which is min(rows(),cols()).
See also:
rows(), cols(), SizeAtCompileTime.
internal::scalar_product_traits< typename internal::traits< Derived >::Scalar, typename internal::traits< OtherDerived >::Scalar >::ReturnType dot ( const MatrixBase< OtherDerived > &  other) const
Returns:
the dot product of *this with other.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

Note:
If the scalar type is complex numbers, then this function returns the hermitian (sesquilinear) dot product, conjugate-linear in the first variable and linear in the second variable.
See also:
squaredNorm(), norm()
MatrixBase< Derived >::EigenvaluesReturnType eigenvalues ( ) const
inline

Computes the eigenvalues of a matrix.

Returns:
Column vector containing the eigenvalues.

This is defined in the Eigenvalues module.

#include <Eigen/Eigenvalues>

This function computes the eigenvalues with the help of the EigenSolver class (for real matrices) or the ComplexEigenSolver class (for complex matrices).

The eigenvalues are repeated according to their algebraic multiplicity, so there are as many eigenvalues as rows in the matrix.

The SelfAdjointView class provides a better algorithm for selfadjoint matrices.

Example:

MatrixXd ones = MatrixXd::Ones(3,3);
VectorXcd eivals = ones.eigenvalues();
cout << "The eigenvalues of the 3x3 matrix of ones are:" << endl << eivals << endl;

Output:

The eigenvalues of the 3x3 matrix of ones are:
(-2.98e-17,0)
(3,0)
(1.81e-32,0)
See also:
EigenSolver::eigenvalues(), ComplexEigenSolver::eigenvalues(), SelfAdjointView::eigenvalues()
const ForceAlignedAccess< Derived > forceAlignedAccess ( ) const
inline
Returns:
an expression of *this with forced aligned access
See also:
forceAlignedAccessIf(),class ForceAlignedAccess

Reimplemented from DenseBase< Derived >.

ForceAlignedAccess< Derived > forceAlignedAccess ( )
inline
Returns:
an expression of *this with forced aligned access
See also:
forceAlignedAccessIf(), class ForceAlignedAccess

Reimplemented from DenseBase< Derived >.

internal::add_const_on_value_type< typename internal::conditional< Enable, ForceAlignedAccess< Derived >, Derived & >::type >::type forceAlignedAccessIf ( ) const
inline
Returns:
an expression of *this with forced aligned access if Enable is true.
See also:
forceAlignedAccess(), class ForceAlignedAccess

Reimplemented from DenseBase< Derived >.

internal::conditional< Enable, ForceAlignedAccess< Derived >, Derived & >::type forceAlignedAccessIf ( )
inline
Returns:
an expression of *this with forced aligned access if Enable is true.
See also:
forceAlignedAccess(), class ForceAlignedAccess

Reimplemented from DenseBase< Derived >.

const FullPivHouseholderQR< typename MatrixBase< Derived >::PlainObject > fullPivHouseholderQr ( ) const
Returns:
the full-pivoting Householder QR decomposition of *this.
See also:
class FullPivHouseholderQR
const FullPivLU< typename MatrixBase< Derived >::PlainObject > fullPivLu ( ) const
inline

This is defined in the LU module.

#include <Eigen/LU>
Returns:
the full-pivoting LU decomposition of *this.
See also:
class FullPivLU
const MatrixBase< Derived >::HNormalizedReturnType hnormalized ( ) const
inline

This is defined in the Geometry module.

#include <Eigen/Geometry>
Returns:
an expression of the homogeneous normalized vector of *this

Example:

Output:

See also:
VectorwiseOp::hnormalized()
MatrixBase< Derived >::HomogeneousReturnType homogeneous ( ) const
inline

This is defined in the Geometry module.

#include <Eigen/Geometry>
Returns:
an expression of the equivalent homogeneous vector

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

Example:

Output:

See also:
class Homogeneous
const HouseholderQR< typename MatrixBase< Derived >::PlainObject > householderQr ( ) const
Returns:
the Householder QR decomposition of *this.
See also:
class HouseholderQR
NumTraits< typename internal::traits< Derived >::Scalar >::Real hypotNorm ( ) const
inline
Returns:
the l2 norm of *this avoiding undeflow and overflow. This version use a concatenation of hypot() calls, and it is very slow.
See also:
norm(), stableNorm()
const MatrixBase< Derived >::IdentityReturnType Identity ( )
inlinestatic
Returns:
an expression of the identity matrix (not necessarily square).

This variant is only for fixed-size MatrixBase types. For dynamic-size types, you need to use the variant taking size arguments.

Example:

cout << Matrix<double, 3, 4>::Identity() << endl;

Output:

1 0 0 0
0 1 0 0
0 0 1 0
See also:
Identity(Index,Index), setIdentity(), isIdentity()
const MatrixBase< Derived >::IdentityReturnType Identity ( Index  rows,
Index  cols 
)
inlinestatic
Returns:
an expression of the identity matrix (not necessarily square).

The parameters rows and cols are the number of rows and of columns of the returned matrix. Must be compatible with this MatrixBase type.

This variant is meant to be used for dynamic-size matrix types. For fixed-size types, it is redundant to pass rows and cols as arguments, so Identity() should be used instead.

Example:

cout << MatrixXd::Identity(4, 3) << endl;

Output:

1 0 0
0 1 0
0 0 1
0 0 0
See also:
Identity(), setIdentity(), isIdentity()
const ImagReturnType imag ( ) const
inline
Returns:
an read-only expression of the imaginary part of *this.
See also:
real()
NonConstImagReturnType imag ( )
inline
Returns:
a non const expression of the imaginary part of *this.
See also:
real()
const internal::inverse_impl< Derived > inverse ( ) const
inline

This is defined in the LU module.

#include <Eigen/LU>
Returns:
the matrix inverse of this matrix.

For small fixed sizes up to 4x4, this method uses cofactors. In the general case, this method uses class PartialPivLU.

Note:
This matrix must be invertible, otherwise the result is undefined. If you need an invertibility check, do the following: Example:
cout << "Here is the matrix m:" << endl << m << endl;
cout << "Its inverse is:" << endl << m.inverse() << endl;
Output:
Here is the matrix m:
  0.68  0.597  -0.33
-0.211  0.823  0.536
 0.566 -0.605 -0.444
Its inverse is:
-0.199   2.23   2.84
  1.01 -0.555  -1.42
 -1.62   3.59   3.29
See also:
computeInverseAndDetWithCheck()
bool isDiagonal ( RealScalar  prec = NumTraits<Scalar>::dummy_precision()) const
Returns:
true if *this is approximately equal to a diagonal matrix, within the precision given by prec.

Example:

m(0,2) = 1;
cout << "Here's the matrix m:" << endl << m << endl;
cout << "m.isDiagonal() returns: " << m.isDiagonal() << endl;
cout << "m.isDiagonal(1e-3) returns: " << m.isDiagonal(1e-3) << endl;

Output:

Here's the matrix m:
1e+04     0     1
    0 1e+04     0
    0     0 1e+04
m.isDiagonal() returns: 0
m.isDiagonal(1e-3) returns: 1
See also:
asDiagonal()
bool isIdentity ( RealScalar  prec = NumTraits<Scalar>::dummy_precision()) const
Returns:
true if *this is approximately equal to the identity matrix (not necessarily square), within the precision given by prec.

Example:

m(0,2) = 1e-4;
cout << "Here's the matrix m:" << endl << m << endl;
cout << "m.isIdentity() returns: " << m.isIdentity() << endl;
cout << "m.isIdentity(1e-3) returns: " << m.isIdentity(1e-3) << endl;

Output:

Here's the matrix m:
     1      0 0.0001
     0      1      0
     0      0      1
m.isIdentity() returns: 0
m.isIdentity(1e-3) returns: 1
See also:
class CwiseNullaryOp, Identity(), Identity(Index,Index), setIdentity()
bool isLowerTriangular ( RealScalar  prec = NumTraits<Scalar>::dummy_precision()) const
Returns:
true if *this is approximately equal to a lower triangular matrix, within the precision given by prec.
See also:
isUpperTriangular()
bool isOrthogonal ( const MatrixBase< OtherDerived > &  other,
RealScalar  prec = NumTraits<Scalar>::dummy_precision() 
) const
Returns:
true if *this is approximately orthogonal to other, within the precision given by prec.

Example:

Vector3d v(1,0,0);
Vector3d w(1e-4,0,1);
cout << "Here's the vector v:" << endl << v << endl;
cout << "Here's the vector w:" << endl << w << endl;
cout << "v.isOrthogonal(w) returns: " << v.isOrthogonal(w) << endl;
cout << "v.isOrthogonal(w,1e-3) returns: " << v.isOrthogonal(w,1e-3) << endl;

Output:

Here's the vector v:
1
0
0
Here's the vector w:
0.0001
0
1
v.isOrthogonal(w) returns: 0
v.isOrthogonal(w,1e-3) returns: 1
bool isUnitary ( RealScalar  prec = NumTraits<Scalar>::dummy_precision()) const
Returns:
true if *this is approximately an unitary matrix, within the precision given by prec. In the case where the Scalar type is real numbers, a unitary matrix is an orthogonal matrix, whence the name.
Note:
This can be used to check whether a family of vectors forms an orthonormal basis. Indeed, m.isUnitary() returns true if and only if the columns (equivalently, the rows) of m form an orthonormal basis.

Example:

m(0,2) = 1e-4;
cout << "Here's the matrix m:" << endl << m << endl;
cout << "m.isUnitary() returns: " << m.isUnitary() << endl;
cout << "m.isUnitary(1e-3) returns: " << m.isUnitary(1e-3) << endl;

Output:

Here's the matrix m:
     1      0 0.0001
     0      1      0
     0      0      1
m.isUnitary() returns: 0
m.isUnitary(1e-3) returns: 1
bool isUpperTriangular ( RealScalar  prec = NumTraits<Scalar>::dummy_precision()) const
Returns:
true if *this is approximately equal to an upper triangular matrix, within the precision given by prec.
See also:
isLowerTriangular()
JacobiSVD< typename MatrixBase< Derived >::PlainObject > jacobiSvd ( unsigned int  computationOptions = 0) const

This is defined in the SVD module.

#include <Eigen/SVD>
Returns:
the singular value decomposition of *this computed by two-sided Jacobi transformations.
See also:
class JacobiSVD
const LazyProductReturnType< Derived, OtherDerived >::Type lazyProduct ( const MatrixBase< OtherDerived > &  other) const
Returns:
an expression of the matrix product of *this and other without implicit evaluation.

The returned product will behave like any other expressions: the coefficients of the product will be computed once at a time as requested. This might be useful in some extremely rare cases when only a small and no coherent fraction of the result's coefficients have to be computed.

Warning:
This version of the matrix product can be much much slower. So use it only if you know what you are doing and that you measured a true speed improvement.
See also:
operator*(const MatrixBase&)
const LDLT< typename MatrixBase< Derived >::PlainObject > ldlt ( ) const
inline

This is defined in the Cholesky module.

#include <Eigen/Cholesky>
Returns:
the Cholesky decomposition with full pivoting without square root of *this
const LLT< typename MatrixBase< Derived >::PlainObject > llt ( ) const
inline

This is defined in the Cholesky module.

#include <Eigen/Cholesky>
Returns:
the LLT decomposition of *this
NumTraits< typename internal::traits< Derived >::Scalar >::Real lpNorm ( ) const
inline
Returns:
the $ \ell^p $ norm of *this, that is, returns the p-th root of the sum of the p-th powers of the absolute values of the coefficients of *this. If p is the special value Eigen::Infinity, this function returns the $ \ell^\infty $ norm, that is the maximum of the absolute values of the coefficients of *this.
See also:
norm()

Reimplemented from DenseBase< Derived >.

const PartialPivLU< typename MatrixBase< Derived >::PlainObject > lu ( ) const
inline

This is defined in the LU module.

#include <Eigen/LU>

Synonym of partialPivLu().

Returns:
the partial-pivoting LU decomposition of *this.
See also:
class PartialPivLU
void makeHouseholder ( EssentialPart &  essential,
Scalar &  tau,
RealScalar &  beta 
) const

Computes the elementary reflector H such that: $ H *this = [ beta 0 ... 0]^T $ where the transformation H is: $ H = I - tau v v^*$ and the vector v is: $ v^T = [1 essential^T] $

On output:

Parameters:
essentialthe essential part of the vector v
tauthe scaling factor of the householder transformation
betathe result of H * *this
See also:
MatrixBase::makeHouseholderInPlace(), MatrixBase::applyHouseholderOnTheLeft(), MatrixBase::applyHouseholderOnTheRight()
NoAlias< Derived, MatrixBase > noalias ( )
Returns:
a pseudo expression of *this with an operator= assuming no aliasing between *this and the source expression.

More precisely, noalias() allows to bypass the EvalBeforeAssignBit flag. Currently, even though several expressions may alias, only product expressions have this flag. Therefore, noalias() is only usefull when the source expression contains a matrix product.

Here are some examples where noalias is usefull:

D.noalias() = A * B;
D.noalias() += A.transpose() * B;
D.noalias() -= 2 * A * B.adjoint();

On the other hand the following example will lead to a wrong result:

A.noalias() = A * B;

because the result matrix A is also an operand of the matrix product. Therefore, there is no alternative than evaluating A * B in a temporary, that is the default behavior when you write:

A = A * B;
See also:
class NoAlias
NumTraits< typename internal::traits< Derived >::Scalar >::Real norm ( ) const
inline
Returns:
, for vectors, the l2 norm of *this, and for matrices the Frobenius norm. In both cases, it consists in the square root of the sum of the square of all the matrix entries. For vectors, this is also equals to the square root of the dot product of *this with itself.
See also:
dot(), squaredNorm()
void normalize ( )
inline

Normalizes the vector, i.e. divides it by its own norm.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
norm(), normalized()
const MatrixBase< Derived >::PlainObject normalized ( ) const
inline
Returns:
an expression of the quotient of *this by its own norm.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
norm(), normalize()
bool operator!= ( const MatrixBase< OtherDerived > &  other) const
inline
Returns:
true if at least one pair of coefficients of *this and other are not exactly equal to each other.
Warning:
When using floating point scalar values you probably should rather use a fuzzy comparison such as isApprox()
See also:
isApprox(), operator==
const ScalarMultipleReturnType operator* ( const Scalar &  scalar) const
inline
Returns:
an expression of *this scaled by the scalar factor scalar
const CwiseUnaryOp<internal::scalar_multiple2_op<Scalar,std::complex<Scalar> >, const Derived> operator* ( const std::complex< Scalar > &  scalar) const
inline

Overloaded for efficient real matrix times complex scalar value

const ProductReturnType< Derived, OtherDerived >::Type operator* ( const MatrixBase< OtherDerived > &  other) const
inline
Returns:
the matrix product of *this and other.
Note:
If instead of the matrix product you want the coefficient-wise product, see Cwise::operator*().
See also:
lazyProduct(), operator*=(const MatrixBase&), Cwise::operator*()
const DiagonalProduct< Derived, DiagonalDerived, OnTheRight > operator* ( const DiagonalBase< DiagonalDerived > &  diagonal) const
inline
Returns:
the diagonal matrix product of *this by the diagonal matrix diagonal.
MatrixBase< Derived >::ScalarMultipleReturnType operator* ( const UniformScaling< Scalar > &  s) const

Concatenates a linear transformation matrix and a uniform scaling

Derived & operator*= ( const EigenBase< OtherDerived > &  other)
inline

replaces *this by *this * other.

Returns:
a reference to *this
Derived & operator+= ( const MatrixBase< OtherDerived > &  other)
inline

replaces *this by *this + other.

Returns:
a reference to *this
const CwiseUnaryOp<internal::scalar_opposite_op<typename internal::traits<Derived>::Scalar>, const Derived> operator- ( ) const
inline
Returns:
an expression of the opposite of *this
Derived & operator-= ( const MatrixBase< OtherDerived > &  other)
inline

replaces *this by *this - other.

Returns:
a reference to *this
const CwiseUnaryOp<internal::scalar_quotient1_op<typename internal::traits<Derived>::Scalar>, const Derived> operator/ ( const Scalar &  scalar) const
inline
Returns:
an expression of *this divided by the scalar value scalar
Derived& operator= ( const MatrixBase< Derived > &  other)

Special case of the template operator=, in order to prevent the compiler from generating a default operator= (issue hit with g++ 4.1)

Derived& operator= ( const DenseBase< OtherDerived > &  other)

Copies other into *this.

Returns:
a reference to *this.

Reimplemented from DenseBase< Derived >.

Derived& operator= ( const EigenBase< OtherDerived > &  other)

Copies the generic expression other into *this.

The expression must provide a (templated) evalTo(Derived& dst) const function which does the actual job. In practice, this allows any user to write its own special matrix without having to modify MatrixBase

Returns:
a reference to *this.

Reimplemented from DenseBase< Derived >.

bool operator== ( const MatrixBase< OtherDerived > &  other) const
inline
Returns:
true if each coefficients of *this and other are all exactly equal.
Warning:
When using floating point scalar values you probably should rather use a fuzzy comparison such as isApprox()
See also:
isApprox(), operator!=
MatrixBase< Derived >::RealScalar operatorNorm ( ) const
inline

Computes the L2 operator norm.

Returns:
Operator norm of the matrix.

This is defined in the Eigenvalues module.

#include <Eigen/Eigenvalues>

This function computes the L2 operator norm of a matrix, which is also known as the spectral norm. The norm of a matrix $ A $ is defined to be

\[ \|A\|_2 = \max_x \frac{\|Ax\|_2}{\|x\|_2} \]

where the maximum is over all vectors and the norm on the right is the Euclidean vector norm. The norm equals the largest singular value, which is the square root of the largest eigenvalue of the positive semi-definite matrix $ A^*A $.

The current implementation uses the eigenvalues of $ A^*A $, as computed by SelfAdjointView::eigenvalues(), to compute the operator norm of a matrix. The SelfAdjointView class provides a better algorithm for selfadjoint matrices.

Example:

MatrixXd ones = MatrixXd::Ones(3,3);
cout << "The operator norm of the 3x3 matrix of ones is "
<< ones.operatorNorm() << endl;

Output:

The operator norm of the 3x3 matrix of ones is 3
See also:
SelfAdjointView::eigenvalues(), SelfAdjointView::operatorNorm()
const PartialPivLU< typename MatrixBase< Derived >::PlainObject > partialPivLu ( ) const
inline

This is defined in the LU module.

#include <Eigen/LU>
Returns:
the partial-pivoting LU decomposition of *this.
See also:
class PartialPivLU
RealReturnType real ( ) const
inline
Returns:
a read-only expression of the real part of *this.
See also:
imag()
NonConstRealReturnType real ( )
inline
Returns:
a non const expression of the real part of *this.
See also:
imag()
Derived & setIdentity ( )
inline

Writes the identity expression (not necessarily square) into *this.

Example:

m.block<3,3>(1,0).setIdentity();
cout << m << endl;

Output:

0 0 0 0
1 0 0 0
0 1 0 0
0 0 1 0
See also:
class CwiseNullaryOp, Identity(), Identity(Index,Index), isIdentity()
Derived & setIdentity ( Index  rows,
Index  cols 
)
inline

Resizes to the given size, and writes the identity expression (not necessarily square) into *this.

Parameters:
rowsthe new number of rows
colsthe new number of columns

Example:

m.setIdentity(3, 3);
cout << m << endl;

Output:

1 0 0
0 1 0
0 0 1
See also:
MatrixBase::setIdentity(), class CwiseNullaryOp, MatrixBase::Identity()
NumTraits< typename internal::traits< Derived >::Scalar >::Real squaredNorm ( ) const
inline
Returns:
, for vectors, the squared l2 norm of *this, and for matrices the Frobenius norm. In both cases, it consists in the sum of the square of all the matrix entries. For vectors, this is also equals to the dot product of *this with itself.
See also:
dot(), norm()
NumTraits< typename internal::traits< Derived >::Scalar >::Real stableNorm ( ) const
inline
Returns:
the l2 norm of *this avoiding underflow and overflow. This version use a blockwise two passes algorithm: 1 - find the absolute largest coefficient s 2 - compute $ s \Vert \frac{*this}{s} \Vert $ in a standard way

For architecture/scalar types supporting vectorization, this version is faster than blueNorm(). Otherwise the blueNorm() is much faster.

See also:
norm(), blueNorm(), hypotNorm()
internal::traits< Derived >::Scalar trace ( ) const
inline
Returns:
the trace of *this, i.e. the sum of the coefficients on the main diagonal.

*this can be any matrix, not necessarily square.

See also:
diagonal(), sum()

Reimplemented from DenseBase< Derived >.

MatrixBase< Derived >::template TriangularViewReturnType< Mode >::Type triangularView ( )
Returns:
an expression of a triangular view extracted from the current matrix

The parameter Mode can have the following values: Upper, StrictlyUpper, UnitUpper, Lower, StrictlyLower, UnitLower.

Example:

#ifndef _MSC_VER
#warning deprecated
#endif
/* deprecated
Matrix3i m = Matrix3i::Random();
cout << "Here is the matrix m:" << endl << m << endl;
cout << "Here is the upper-triangular matrix extracted from m:" << endl
<< m.part<Eigen::UpperTriangular>() << endl;
cout << "Here is the strictly-upper-triangular matrix extracted from m:" << endl
<< m.part<Eigen::StrictlyUpperTriangular>() << endl;
cout << "Here is the unit-lower-triangular matrix extracted from m:" << endl
<< m.part<Eigen::UnitLowerTriangular>() << endl;
*/

Output:

See also:
class TriangularView
MatrixBase< Derived >::template ConstTriangularViewReturnType< Mode >::Type triangularView ( ) const

This is the const version of MatrixBase::triangularView()

const CwiseUnaryOp<CustomUnaryOp, const Derived> unaryExpr ( const CustomUnaryOp &  func = CustomUnaryOp()) const
inline

Apply a unary operator coefficient-wise.

Parameters:
[in]funcFunctor implementing the unary operator
Template Parameters:
CustomUnaryOpType of func
Returns:
An expression of a custom coefficient-wise unary operator func of *this

The function ptr_fun() from the C++ standard library can be used to make functors out of normal functions.

Example:

#include <Eigen/Core>
#include <iostream>
using namespace Eigen;
using namespace std;
// define function to be applied coefficient-wise
double ramp(double x)
{
if (x > 0)
return x;
else
return 0;
}
int main(int, char**)
{
cout << m1 << endl << "becomes: " << endl << m1.unaryExpr(ptr_fun(ramp)) << endl;
return 0;
}

Output:

   0.68   0.823  -0.444   -0.27
 -0.211  -0.605   0.108  0.0268
  0.566   -0.33 -0.0452   0.904
  0.597   0.536   0.258   0.832
becomes: 
  0.68  0.823      0      0
     0      0  0.108 0.0268
 0.566      0      0  0.904
 0.597  0.536  0.258  0.832

Genuine functors allow for more possibilities, for instance it may contain a state.

Example:

#include <Eigen/Core>
#include <iostream>
using namespace Eigen;
using namespace std;
// define a custom template unary functor
template<typename Scalar>
struct CwiseClampOp {
CwiseClampOp(const Scalar& inf, const Scalar& sup) : m_inf(inf), m_sup(sup) {}
const Scalar operator()(const Scalar& x) const { return x<m_inf ? m_inf : (x>m_sup ? m_sup : x); }
Scalar m_inf, m_sup;
};
int main(int, char**)
{
cout << m1 << endl << "becomes: " << endl << m1.unaryExpr(CwiseClampOp<double>(-0.5,0.5)) << endl;
return 0;
}

Output:

   0.68   0.823  -0.444   -0.27
 -0.211  -0.605   0.108  0.0268
  0.566   -0.33 -0.0452   0.904
  0.597   0.536   0.258   0.832
becomes: 
    0.5     0.5  -0.444   -0.27
 -0.211    -0.5   0.108  0.0268
    0.5   -0.33 -0.0452     0.5
    0.5     0.5   0.258     0.5
See also:
class CwiseUnaryOp, class CwiseBinaryOp
const CwiseUnaryView<CustomViewOp, const Derived> unaryViewExpr ( const CustomViewOp &  func = CustomViewOp()) const
inline
Returns:
an expression of a custom coefficient-wise unary operator func of *this

The template parameter CustomUnaryOp is the type of the functor of the custom unary operator.

Example:

#include <Eigen/Core>
#include <iostream>
using namespace Eigen;
using namespace std;
// define a custom template unary functor
template<typename Scalar>
struct CwiseClampOp {
CwiseClampOp(const Scalar& inf, const Scalar& sup) : m_inf(inf), m_sup(sup) {}
const Scalar operator()(const Scalar& x) const { return x<m_inf ? m_inf : (x>m_sup ? m_sup : x); }
Scalar m_inf, m_sup;
};
int main(int, char**)
{
cout << m1 << endl << "becomes: " << endl << m1.unaryExpr(CwiseClampOp<double>(-0.5,0.5)) << endl;
return 0;
}

Output:

   0.68   0.823  -0.444   -0.27
 -0.211  -0.605   0.108  0.0268
  0.566   -0.33 -0.0452   0.904
  0.597   0.536   0.258   0.832
becomes: 
    0.5     0.5  -0.444   -0.27
 -0.211    -0.5   0.108  0.0268
    0.5   -0.33 -0.0452     0.5
    0.5     0.5   0.258     0.5
See also:
class CwiseUnaryOp, class CwiseBinaryOp
const MatrixBase< Derived >::BasisReturnType Unit ( Index  size,
Index  i 
)
inlinestatic
Returns:
an expression of the i-th unit (basis) vector.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
MatrixBase::Unit(Index), MatrixBase::UnitX(), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()
const MatrixBase< Derived >::BasisReturnType Unit ( Index  i)
inlinestatic
Returns:
an expression of the i-th unit (basis) vector.

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

This variant is for fixed-size vector only.

See also:
MatrixBase::Unit(Index,Index), MatrixBase::UnitX(), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()
MatrixBase< Derived >::PlainObject unitOrthogonal ( void  ) const
Returns:
a unit vector which is orthogonal to *this

The size of *this must be at least 2. If the size is exactly 2, then the returned vector is a counter clock wise rotation of *this, i.e., (-y,x).normalized().

See also:
cross()
const MatrixBase< Derived >::BasisReturnType UnitW ( )
inlinestatic
Returns:
an expression of the W axis unit vector (0,0,0,1)

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
MatrixBase::Unit(Index,Index), MatrixBase::Unit(Index), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()
const MatrixBase< Derived >::BasisReturnType UnitX ( )
inlinestatic
Returns:
an expression of the X axis unit vector (1{,0}^*)

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
MatrixBase::Unit(Index,Index), MatrixBase::Unit(Index), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()
const MatrixBase< Derived >::BasisReturnType UnitY ( )
inlinestatic
Returns:
an expression of the Y axis unit vector (0,1{,0}^*)

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
MatrixBase::Unit(Index,Index), MatrixBase::Unit(Index), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()
const MatrixBase< Derived >::BasisReturnType UnitZ ( )
inlinestatic
Returns:
an expression of the Z axis unit vector (0,0,1{,0}^*)

This is only for vectors (either row-vectors or column-vectors), i.e. matrices which are known at compile-time to have either one row or one column.

See also:
MatrixBase::Unit(Index,Index), MatrixBase::Unit(Index), MatrixBase::UnitY(), MatrixBase::UnitZ(), MatrixBase::UnitW()

The documentation for this class was generated from the following files: