![]() |
Eigen
3.3.3
|
Expression of a selfadjoint matrix from a triangular part of a dense matrix.
MatrixType | the type of the dense matrix storing the coefficients |
TriangularPart | can 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.
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 | |
const AdjointReturnType | adjoint () const |
Scalar | coeff (Index row, Index col) const |
Scalar & | coeffRef (Index row, Index col) |
Index | cols () const |
const ConjugateReturnType | conjugate () const |
MatrixType::ConstDiagonalReturnType | diagonal () 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 > | |
const Product< SelfAdjointView, OtherDerived > | operator* (const MatrixBase< OtherDerived > &rhs) const |
RealScalar | operatorNorm () const |
Computes the L2 operator norm. | |
template<typename DerivedU , typename DerivedV > | |
SelfAdjointView & | rankUpdate (const MatrixBase< DerivedU > &u, const MatrixBase< DerivedV > &v, const Scalar &alpha=Scalar(1)) |
template<typename DerivedU > | |
SelfAdjointView & | rankUpdate (const MatrixBase< DerivedU > &u, const Scalar &alpha=Scalar(1)) |
Index | rows () const |
TransposeReturnType | transpose () |
const ConstTransposeReturnType | transpose () const |
Public Attributes | |
internal::conditional<(TriMode &(Upper|Lower))==(UpLo &(Upper|Lower)), TriangularView< MatrixType, TriMode >, TriangularView < typename MatrixType::AdjointReturnType, TriMode > >::type | triangularView () const |
Friends | |
template<typename OtherDerived > | |
const Product< OtherDerived, SelfAdjointView > | operator* (const MatrixBase< OtherDerived > &lhs, const SelfAdjointView &rhs) |
typedef Matrix<RealScalar, internal::traits<MatrixType>::ColsAtCompileTime, 1> Eigen::SelfAdjointView< _MatrixType, UpLo >::EigenvaluesReturnType |
Return type of eigenvalues()
typedef NumTraits<Scalar>::Real Eigen::SelfAdjointView< _MatrixType, UpLo >::RealScalar |
Real part of Scalar
const AdjointReturnType Eigen::SelfAdjointView< _MatrixType, UpLo >::adjoint | ( | ) | const [inline] |
Scalar Eigen::SelfAdjointView< _MatrixType, UpLo >::coeff | ( | Index | row, |
Index | col | ||
) | const [inline] |
Reimplemented from Eigen::TriangularBase< SelfAdjointView< _MatrixType, UpLo > >.
Scalar& Eigen::SelfAdjointView< _MatrixType, UpLo >::coeffRef | ( | Index | row, |
Index | col | ||
) | [inline] |
Reimplemented from Eigen::TriangularBase< SelfAdjointView< _MatrixType, UpLo > >.
Index Eigen::SelfAdjointView< _MatrixType, UpLo >::cols | ( | void | ) | const [inline] |
Reimplemented from Eigen::TriangularBase< SelfAdjointView< _MatrixType, UpLo > >.
const ConjugateReturnType Eigen::SelfAdjointView< _MatrixType, UpLo >::conjugate | ( | ) | const [inline] |
MatrixType::ConstDiagonalReturnType Eigen::SelfAdjointView< _MatrixType, UpLo >::diagonal | ( | ) | const [inline] |
*this
This method simply returns the diagonal of the nested expression, thus by-passing the SelfAdjointView decorator.
SelfAdjointView< MatrixType, UpLo >::EigenvaluesReturnType Eigen::SelfAdjointView< MatrixType, UpLo >::eigenvalues | ( | ) | const [inline] |
Computes the eigenvalues of a matrix.
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: -3.09e-16 0 3
const LDLT< typename SelfAdjointView< MatrixType, UpLo >::PlainObject, UpLo > Eigen::SelfAdjointView< MatrixType, UpLo >::ldlt | ( | ) | const [inline] |
This is defined in the Cholesky module.
#include <Eigen/Cholesky>
*this
const LLT< typename SelfAdjointView< MatrixType, UpLo >::PlainObject, UpLo > Eigen::SelfAdjointView< MatrixType, UpLo >::llt | ( | ) | const [inline] |
This is defined in the Cholesky module.
#include <Eigen/Cholesky>
*this
const Product<SelfAdjointView,OtherDerived> Eigen::SelfAdjointView< _MatrixType, UpLo >::operator* | ( | const MatrixBase< OtherDerived > & | rhs | ) | const [inline] |
Efficient triangular matrix times vector/matrix product
SelfAdjointView< MatrixType, UpLo >::RealScalar Eigen::SelfAdjointView< MatrixType, UpLo >::operatorNorm | ( | ) | const [inline] |
Computes the L2 operator norm.
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
SelfAdjointView& Eigen::SelfAdjointView< _MatrixType, UpLo >::rankUpdate | ( | const MatrixBase< DerivedU > & | u, |
const MatrixBase< DerivedV > & | v, | ||
const Scalar & | alpha = Scalar(1) |
||
) |
Perform a symmetric rank 2 update of the selfadjoint matrix *this
:
*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.
SelfAdjointView& Eigen::SelfAdjointView< _MatrixType, UpLo >::rankUpdate | ( | const MatrixBase< DerivedU > & | u, |
const Scalar & | alpha = Scalar(1) |
||
) |
Perform a symmetric rank K update of the selfadjoint matrix *this
: where u is a vector or matrix.
*this
Note that to perform you can simply call this function with u.adjoint().
Index Eigen::SelfAdjointView< _MatrixType, UpLo >::rows | ( | void | ) | const [inline] |
Reimplemented from Eigen::TriangularBase< SelfAdjointView< _MatrixType, UpLo > >.
TransposeReturnType Eigen::SelfAdjointView< _MatrixType, UpLo >::transpose | ( | ) | [inline] |
const ConstTransposeReturnType Eigen::SelfAdjointView< _MatrixType, UpLo >::transpose | ( | ) | const [inline] |
const Product<OtherDerived,SelfAdjointView> operator* | ( | const MatrixBase< OtherDerived > & | lhs, |
const SelfAdjointView< _MatrixType, UpLo > & | rhs | ||
) | [friend] |
Efficient vector/matrix times triangular matrix product
internal::conditional<(TriMode&(Upper|Lower))==(UpLo&(Upper|Lower)), TriangularView<MatrixType,TriMode>, TriangularView<typename MatrixType::AdjointReturnType,TriMode> >::type Eigen::SelfAdjointView< _MatrixType, UpLo >::triangularView() const [inline] |
The parameter TriMode can have the following values: #Upper
, #StrictlyUpper
, #UnitUpper
, #Lower
, #StrictlyLower
, #UnitLower
.
If TriMode
references the same triangular part than *this
, then this method simply return a TriangularView
of the nested expression, otherwise, the nested expression is first transposed, thus returning a TriangularView<Transpose<MatrixType>>
object.