11 #ifndef EIGEN_TRIDIAGONALIZATION_H 12 #define EIGEN_TRIDIAGONALIZATION_H 19 template<
typename MatrixType>
21 :
public traits<typename MatrixType::PlainObject>
27 template<
typename MatrixType,
typename CoeffVectorType>
76 Size = MatrixType::RowsAtCompileTime,
78 Options = MatrixType::Options,
79 MaxSize = MatrixType::MaxRowsAtCompileTime,
80 MaxSizeMinusOne = MaxSize ==
Dynamic ?
Dynamic : (MaxSize > 1 ? MaxSize - 1 : 1)
94 typedef typename internal::conditional<NumTraits<Scalar>::IsComplex,
117 m_isInitialized(false)
130 template<
typename InputType>
132 : m_matrix(matrix.derived()),
133 m_hCoeffs(matrix.
cols() > 1 ? matrix.
cols()-1 : 1),
134 m_isInitialized(false)
137 m_isInitialized =
true;
157 template<
typename InputType>
161 m_hCoeffs.resize(matrix.
rows()-1, 1);
163 m_isInitialized =
true;
185 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
222 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
243 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
244 return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate())
245 .setLength(m_matrix.rows() - 1)
268 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
269 return MatrixTReturnType(m_matrix.real());
285 DiagonalReturnType
diagonal()
const;
297 SubDiagonalReturnType subDiagonal()
const;
306 template<
typename MatrixType>
310 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
311 return m_matrix.diagonal().real();
314 template<
typename MatrixType>
318 eigen_assert(m_isInitialized &&
"Tridiagonalization is not initialized.");
347 template<
typename MatrixType,
typename CoeffVectorType>
354 Index
n = matA.rows();
358 for (Index
i = 0;
i<n-1; ++
i)
360 Index remainingSize = n-
i-1;
363 matA.col(i).tail(remainingSize).makeHouseholderInPlace(h, beta);
367 matA.col(i).coeffRef(i+1) = 1;
369 hCoeffs.tail(n-i-1).noalias() = (matA.bottomRightCorner(remainingSize,remainingSize).template selfadjointView<Lower>()
370 * (
conj(h) * matA.col(i).tail(remainingSize)));
372 hCoeffs.tail(n-i-1) += (
conj(h)*
RealScalar(-0.5)*(hCoeffs.tail(remainingSize).dot(matA.col(i).tail(remainingSize)))) * matA.col(i).tail(n-i-1);
374 matA.bottomRightCorner(remainingSize, remainingSize).template selfadjointView<Lower>()
375 .rankUpdate(matA.col(i).tail(remainingSize), hCoeffs.tail(remainingSize),
Scalar(-1));
377 matA.col(i).coeffRef(i+1) = beta;
378 hCoeffs.coeffRef(i) =
h;
384 int Size=MatrixType::ColsAtCompileTime,
428 template<
typename MatrixType,
typename DiagonalType,
typename SubDiagonalType,
typename CoeffVectorType>
431 CoeffVectorType& hcoeffs,
bool extractQ)
433 eigen_assert(mat.cols()==mat.rows() && diag.size()==mat.rows() && subdiag.size()==mat.rows()-1);
440 template<
typename MatrixType,
int Size,
bool IsComplex>
445 template<
typename DiagonalType,
typename SubDiagonalType>
447 void run(MatrixType&
mat, DiagonalType&
diag, SubDiagonalType& subdiag, CoeffVectorType& hCoeffs,
bool extractQ)
450 diag = mat.diagonal().real();
453 mat = HouseholderSequenceType(mat, hCoeffs.conjugate())
454 .setLength(mat.rows() - 1)
463 template<
typename MatrixType>
469 template<
typename DiagonalType,
typename SubDiagonalType,
typename CoeffVectorType>
470 static void run(MatrixType&
mat, DiagonalType&
diag, SubDiagonalType& subdiag, CoeffVectorType&,
bool extractQ)
480 subdiag[0] =
mat(1,0);
481 subdiag[1] =
mat(2,1);
489 Scalar m01 =
mat(1,0) * invBeta;
490 Scalar m02 =
mat(2,0) * invBeta;
492 diag[1] =
mat(1,1) + m02*
q;
493 diag[2] =
mat(2,2) - m02*
q;
495 subdiag[1] =
mat(2,1) - m01 *
q;
509 template<
typename MatrixType,
bool IsComplex>
514 template<
typename DiagonalType,
typename SubDiagonalType,
typename CoeffVectorType>
516 void run(MatrixType&
mat, DiagonalType&
diag, SubDiagonalType&, CoeffVectorType&,
bool extractQ)
532 :
public ReturnByValue<TridiagonalizationMatrixTReturnType<MatrixType> >
541 template <
typename ResultType>
545 result.template diagonal<1>() = m_matrix.template diagonal<-1>().conjugate();
546 result.diagonal() = m_matrix.diagonal();
547 result.template
diagonal<-1>() = m_matrix.template diagonal<-1>();
561 #endif // EIGEN_TRIDIAGONALIZATION_H Tridiagonalization(Index size=Size==Dynamic ? 2 :Size)
Default constructor.
DiagonalReturnType diagonal() const
Returns the diagonal of the tridiagonal matrix T in the decomposition.
HouseholderSequence< MatrixType, typename internal::remove_all< typename CoeffVectorType::ConjugateReturnType >::type > HouseholderSequenceType
Return type of matrixQ()
Matrix< Scalar, SizeMinusOne, 1, Options &~RowMajor, MaxSizeMinusOne, 1 > CoeffVectorType
static EIGEN_DEVICE_FUNC void run(MatrixType &mat, DiagonalType &diag, SubDiagonalType &subdiag, CoeffVectorType &hCoeffs, bool extractQ)
const AutoDiffScalar< DerType > & conj(const AutoDiffScalar< DerType > &x)
Tridiagonalization< MatrixType >::HouseholderSequenceType HouseholderSequenceType
internal::remove_all< typename MatrixType::RealReturnType >::type MatrixTypeRealView
static EIGEN_DEVICE_FUNC void run(MatrixType &mat, DiagonalType &diag, SubDiagonalType &, CoeffVectorType &, bool extractQ)
Matrix diag(const std::vector< Matrix > &Hs)
internal::conditional< NumTraits< Scalar >::IsComplex, typename internal::add_const_on_value_type< typename Diagonal< const MatrixType, -1 >::RealReturnType >::type, const Diagonal< const MatrixType, -1 > >::type SubDiagonalReturnType
void diagonal(const MatrixType &m)
HouseholderSequenceType matrixQ() const
Returns the unitary matrix Q in the decomposition.
Namespace containing all symbols from the Eigen library.
Tridiagonalization & compute(const EigenBase< InputType > &matrix)
Computes tridiagonal decomposition of given matrix.
Holds information about the various numeric (i.e. scalar) types allowed by Eigen. ...
Tridiagonal decomposition of a selfadjoint matrix.
MatrixType::Nested m_matrix
AnnoyingScalar conj(const AnnoyingScalar &x)
EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
Sequence of Householder reflections acting on subspaces with decreasing size.
const MatrixType & packedMatrix() const
Returns the internal representation of the decomposition.
CoeffVectorType householderCoefficients() const
Returns the Householder coefficients.
Tridiagonalization< MatrixType >::CoeffVectorType CoeffVectorType
CoeffVectorType m_hCoeffs
EIGEN_DEFAULT_DENSE_INDEX_TYPE Index
The Index type as used for the API.
MatrixTReturnType matrixT() const
Returns an expression of the tridiagonal matrix T in the decomposition.
EIGEN_DEVICE_FUNC const Scalar & q
MatrixType::Scalar Scalar
EIGEN_DEVICE_FUNC EIGEN_CONSTEXPR Index rows() const EIGEN_NOEXCEPT
MatrixType::Scalar Scalar
Tridiagonalization(const EigenBase< InputType > &matrix)
Constructor; computes tridiagonal decomposition of given matrix.
NumTraits< Scalar >::Real RealScalar
EIGEN_CONSTEXPR Index size(const T &x)
#define EIGEN_DEVICE_FUNC
EIGEN_CONSTEXPR Index cols() const EIGEN_NOEXCEPT
internal::TridiagonalizationMatrixTReturnType< MatrixTypeRealView > MatrixTReturnType
internal::plain_col_type< MatrixType, RealScalar >::type DiagonalType
TridiagonalizationMatrixTReturnType(const MatrixType &mat)
Constructor.
MatrixType::RealScalar RealScalar
void evalTo(ResultType &result) const
internal::conditional< NumTraits< Scalar >::IsComplex, typename internal::add_const_on_value_type< typename Diagonal< const MatrixType >::RealReturnType >::type, const Diagonal< const MatrixType > >::type DiagonalReturnType
EIGEN_DEVICE_FUNC void tridiagonalization_inplace(MatrixType &matA, CoeffVectorType &hCoeffs)
Expression of a diagonal/subdiagonal/superdiagonal in a matrix.
Jet< T, N > sqrt(const Jet< T, N > &f)
MatrixType::Scalar Scalar
Matrix< RealScalar, SizeMinusOne, 1, Options &~RowMajor, MaxSizeMinusOne, 1 > SubDiagonalType
Generic expression where a coefficient-wise unary operator is applied to an expression.
static void run(MatrixType &mat, DiagonalType &diag, SubDiagonalType &subdiag, CoeffVectorType &, bool extractQ)
Map< Matrix< T, Dynamic, Dynamic, ColMajor >, 0, OuterStride<> > matrix(T *data, int rows, int cols, int stride)
_MatrixType MatrixType
Synonym for the template parameter _MatrixType.
const AutoDiffScalar< DerType > & real(const AutoDiffScalar< DerType > &x)
NumTraits< Scalar >::Real RealScalar
EIGEN_DEVICE_FUNC bool abs2(bool x)
EIGEN_DEVICE_FUNC Derived & derived()
MatrixType::PlainObject ReturnType
SubDiagonalReturnType subDiagonal() const
Returns the subdiagonal of the tridiagonal matrix T in the decomposition.