Public Types | Public Member Functions | Friends
SelfAdjointView< MatrixType, UpLo > Class Template Reference

Expression of a selfadjoint matrix from a triangular part of a dense matrix. More...

Inherits TriangularBase< Derived >.

List of all members.

Public Types

typedef Matrix< RealScalar,
internal::traits< MatrixType >
::ColsAtCompileTime, 1 > 
EigenvaluesReturnType
typedef NumTraits< Scalar >::Real RealScalar
typedef internal::traits
< SelfAdjointView >::Scalar 
Scalar
 The type of coefficients in this matrix.

Public Member Functions

Scalar coeff (Index row, Index col) const
ScalarcoeffRef (Index row, Index col)
Index cols () const
EigenvaluesReturnType eigenvalues () const
 Computes the eigenvalues of a matrix.
const LDLT< PlainObject, UpLo > ldlt () const
const LLT< PlainObject, UpLo > llt () const
template<typename OtherDerived >
SelfadjointProductMatrix
< MatrixType, Mode, false,
OtherDerived,
0, OtherDerived::IsVectorAtCompileTime > 
operator* (const MatrixBase< OtherDerived > &rhs) const
RealScalar operatorNorm () const
 Computes the L2 operator norm.
template<typename DerivedU , typename DerivedV >
SelfAdjointViewrankUpdate (const MatrixBase< DerivedU > &u, const MatrixBase< DerivedV > &v, Scalar alpha=Scalar(1))
template<typename DerivedU >
SelfAdjointViewrankUpdate (const MatrixBase< DerivedU > &u, Scalar alpha=Scalar(1))
Index rows () const

Friends

template<typename OtherDerived >
SelfadjointProductMatrix
< OtherDerived,
0, OtherDerived::IsVectorAtCompileTime,
MatrixType, Mode, false > 
operator* (const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs)

Detailed Description

template<typename MatrixType, unsigned int UpLo>
class Eigen::SelfAdjointView< MatrixType, UpLo >

Expression of a selfadjoint matrix from a triangular part of a dense matrix.

Parameters:
MatrixTypethe type of the dense matrix storing the coefficients
TriangularPartcan be either Lower or Upper

This class is an expression of a sefladjoint matrix from a triangular part of a matrix with given dense storage of the coefficients. It is the return type of MatrixBase::selfadjointView() and most of the time this is the only way that it is used.

See also:
class TriangularBase, MatrixBase::selfadjointView()

Member Typedef Documentation

typedef Matrix<RealScalar, internal::traits<MatrixType>::ColsAtCompileTime, 1> EigenvaluesReturnType

Return type of eigenvalues()

typedef NumTraits<Scalar>::Real RealScalar

Real part of Scalar


Member Function Documentation

Scalar coeff ( Index  row,
Index  col 
) const [inline]
See also:
MatrixBase::coeff()
Warning:
the coordinates must fit into the referenced triangular part
Scalar& coeffRef ( Index  row,
Index  col 
) [inline]
See also:
MatrixBase::coeffRef()
Warning:
the coordinates must fit into the referenced triangular part
Index cols ( void  ) const [inline]
Returns:
the number of columns.
See also:
rows(), ColsAtCompileTime

Reimplemented from EigenBase< Derived >.

SelfAdjointView< MatrixType, UpLo >::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 SelfAdjointEigenSolver class. The eigenvalues are repeated according to their algebraic multiplicity, so there are as many eigenvalues as rows in the matrix.

Example:

MatrixXd ones = MatrixXd::Ones(3,3);
VectorXd eivals = ones.selfadjointView<Lower>().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.39e-16
8.66e-17
3
See also:
SelfAdjointEigenSolver::eigenvalues(), MatrixBase::eigenvalues()
const LDLT< typename SelfAdjointView< MatrixType, UpLo >::PlainObject, UpLo > 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 SelfAdjointView< MatrixType, UpLo >::PlainObject, UpLo > llt ( ) const [inline]

This is defined in the Cholesky module.

 #include <Eigen/Cholesky> 
Returns:
the LLT decomposition of *this
SelfadjointProductMatrix<MatrixType,Mode,false,OtherDerived,0,OtherDerived::IsVectorAtCompileTime> operator* ( const MatrixBase< OtherDerived > &  rhs) const [inline]

Efficient self-adjoint matrix times vector/matrix product

SelfAdjointView< MatrixType, UpLo >::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 self-adjoint matrix. For a self-adjoint matrix, the operator norm is the largest eigenvalue.

The current implementation uses the eigenvalues of the matrix, as computed by eigenvalues(), to compute the operator norm of the matrix.

Example:

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

Output:

The operator norm of the 3x3 matrix of ones is 3
See also:
eigenvalues(), MatrixBase::operatorNorm()
SelfAdjointView< MatrixType, UpLo > & rankUpdate ( const MatrixBase< DerivedU > &  u,
const MatrixBase< DerivedV > &  v,
Scalar  alpha = Scalar(1) 
)

Perform a symmetric rank 2 update of the selfadjoint matrix *this: $ this = this + \alpha u v^* + conj(\alpha) v u^* $

Returns:
a reference to *this

The vectors u and v must be column vectors, however they can be a adjoint expression without any overhead. Only the meaningful triangular part of the matrix is updated, the rest is left unchanged.

See also:
rankUpdate(const MatrixBase<DerivedU>&, Scalar)
SelfAdjointView< MatrixType, UpLo > & rankUpdate ( const MatrixBase< DerivedU > &  u,
Scalar  alpha = Scalar(1) 
)

Perform a symmetric rank K update of the selfadjoint matrix *this: $ this = this + \alpha ( u u^* ) $ where u is a vector or matrix.

Returns:
a reference to *this

Note that to perform $ this = this + \alpha ( u^* u ) $ you can simply call this function with u.adjoint().

See also:
rankUpdate(const MatrixBase<DerivedU>&, const MatrixBase<DerivedV>&, Scalar)
Index rows ( void  ) const [inline]
Returns:
the number of rows.
See also:
cols(), RowsAtCompileTime

Reimplemented from EigenBase< Derived >.


Friends And Related Function Documentation

SelfadjointProductMatrix<OtherDerived,0,OtherDerived::IsVectorAtCompileTime,MatrixType,Mode,false> operator* ( const MatrixBase< OtherDerived > &  lhs,
const SelfAdjointView< MatrixType, UpLo > &  rhs 
) [friend]

Efficient vector/matrix times self-adjoint matrix product


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