![]() |
Eigen
3.3.3
|
Base class for a triangular part in a dense matrix.
This class is an abstract base class of class TriangularView, and objects of type TriangularViewImpl cannot be instantiated. It extends class TriangularView with additional methods which available for dense expressions only.
Public Member Functions | |
Scalar | coeff (Index row, Index col) const |
Scalar & | coeffRef (Index row, Index col) |
void | fill (const Scalar &value) |
Index | innerStride () const |
template<typename OtherDerived > | |
const Product < TriangularViewType, OtherDerived > | operator* (const MatrixBase< OtherDerived > &rhs) const |
TriangularViewType & | operator*= (const typename internal::traits< MatrixType >::Scalar &other) |
template<typename Other > | |
TriangularViewType & | operator+= (const DenseBase< Other > &other) |
template<typename Other > | |
TriangularViewType & | operator-= (const DenseBase< Other > &other) |
TriangularViewType & | operator/= (const typename internal::traits< MatrixType >::Scalar &other) |
template<typename OtherDerived > | |
TriangularViewType & | operator= (const TriangularBase< OtherDerived > &other) |
template<typename OtherDerived > | |
TriangularViewType & | operator= (const MatrixBase< OtherDerived > &other) |
Index | outerStride () const |
TriangularViewType & | setConstant (const Scalar &value) |
TriangularViewType & | setOnes () |
TriangularViewType & | setZero () |
template<int Side, typename Other > | |
const internal::triangular_solve_retval < Side, TriangularViewType, Other > | solve (const MatrixBase< Other > &other) const |
template<int Side, typename OtherDerived > | |
void | solveInPlace (const MatrixBase< OtherDerived > &other) const |
template<typename OtherDerived > | |
void | swap (TriangularBase< OtherDerived > &other) |
template<typename OtherDerived > | |
void | swap (MatrixBase< OtherDerived > const &other) |
Friends | |
template<typename OtherDerived > | |
const Product< OtherDerived, TriangularViewType > | operator* (const MatrixBase< OtherDerived > &lhs, const TriangularViewImpl &rhs) |
Scalar Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::coeff | ( | Index | row, |
Index | col | ||
) | const [inline] |
Reimplemented from Eigen::TriangularBase< TriangularView< _MatrixType, _Mode > >.
Scalar& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::coeffRef | ( | Index | row, |
Index | col | ||
) | [inline] |
Reimplemented from Eigen::TriangularBase< TriangularView< _MatrixType, _Mode > >.
void Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::fill | ( | const Scalar & | value | ) | [inline] |
Index Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::innerStride | ( | ) | const [inline] |
Reimplemented from Eigen::TriangularBase< TriangularView< _MatrixType, _Mode > >.
const Product<TriangularViewType,OtherDerived> Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator* | ( | const MatrixBase< OtherDerived > & | rhs | ) | const [inline] |
Efficient triangular matrix times vector/matrix product
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator*= | ( | const typename internal::traits< MatrixType >::Scalar & | other | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator+= | ( | const DenseBase< Other > & | other | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator-= | ( | const DenseBase< Other > & | other | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator/= | ( | const typename internal::traits< MatrixType >::Scalar & | other | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator= | ( | const TriangularBase< OtherDerived > & | other | ) |
Assigns a triangular matrix to a triangular part of a dense matrix
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::operator= | ( | const MatrixBase< OtherDerived > & | other | ) |
Shortcut for
*this = other.other.triangularView<(*this)::Mode>()
Index Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::outerStride | ( | ) | const [inline] |
Reimplemented from Eigen::TriangularBase< TriangularView< _MatrixType, _Mode > >.
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::setConstant | ( | const Scalar & | value | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::setOnes | ( | ) | [inline] |
TriangularViewType& Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::setZero | ( | ) | [inline] |
const internal::triangular_solve_retval<Side,TriangularViewType, Other> Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::solve | ( | const MatrixBase< Other > & | other | ) | const [inline] |
*this
with other, *this being triangular.This function computes the inverse-matrix matrix product inverse(*this
) * other if Side==OnTheLeft (the default), or the right-inverse-multiply other * inverse(*this
) if Side==OnTheRight.
Note that the template parameter Side
can be ommitted, in which case Side==OnTheLeft
The matrix *this
must be triangular and invertible (i.e., all the coefficients of the diagonal must be non zero). It works as a forward (resp. backward) substitution if *this
is an upper (resp. lower) triangular matrix.
Example:
Matrix3d m = Matrix3d::Zero(); m.triangularView<Eigen::Upper>().setOnes(); cout << "Here is the matrix m:\n" << m << endl; Matrix3d n = Matrix3d::Ones(); n.triangularView<Eigen::Lower>() *= 2; cout << "Here is the matrix n:\n" << n << endl; cout << "And now here is m.inverse()*n, taking advantage of the fact that" " m is upper-triangular:\n" << m.triangularView<Eigen::Upper>().solve(n) << endl; cout << "And this is n*m.inverse():\n" << m.triangularView<Eigen::Upper>().solve<Eigen::OnTheRight>(n);
Output:
Here is the matrix m: 1 1 1 0 1 1 0 0 1 Here is the matrix n: 2 1 1 2 2 1 2 2 2 And now here is m.inverse()*n, taking advantage of the fact that m is upper-triangular: 0 -1 0 0 0 -1 2 2 2 And this is n*m.inverse(): 2 -1 0 2 0 -1 2 0 0
This function returns an expression of the inverse-multiply and can works in-place if it is assigned to the same matrix or vector other.
For users coming from BLAS, this function (and more specifically solveInPlace()) offer all the operations supported by the *TRSV
and *TRSM
BLAS routines.
void Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::solveInPlace | ( | const MatrixBase< OtherDerived > & | other | ) | const |
"in-place" version of TriangularView::solve() where the result is written in other
Note that the template parameter Side
can be ommitted, in which case Side==OnTheLeft
See TriangularView:solve() for the details.
void Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::swap | ( | TriangularBase< OtherDerived > & | other | ) | [inline] |
Swaps the coefficients of the common triangular parts of two matrices
void Eigen::TriangularViewImpl< _MatrixType, _Mode, Dense >::swap | ( | MatrixBase< OtherDerived > const & | other | ) | [inline] |
(*this).swap(other.triangularView<(*this)::Mode>())
const Product<OtherDerived,TriangularViewType> operator* | ( | const MatrixBase< OtherDerived > & | lhs, |
const TriangularViewImpl< _MatrixType, _Mode, Dense > & | rhs | ||
) | [friend] |
Efficient vector/matrix times triangular matrix product