aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src/Eigenvalues/Tridiagonalization.h
diff options
context:
space:
mode:
Diffstat (limited to 'Eigen/src/Eigenvalues/Tridiagonalization.h')
-rw-r--r--Eigen/src/Eigenvalues/Tridiagonalization.h64
1 files changed, 37 insertions, 27 deletions
diff --git a/Eigen/src/Eigenvalues/Tridiagonalization.h b/Eigen/src/Eigenvalues/Tridiagonalization.h
index 23ae748d4..fac79cf96 100644
--- a/Eigen/src/Eigenvalues/Tridiagonalization.h
+++ b/Eigen/src/Eigenvalues/Tridiagonalization.h
@@ -26,6 +26,11 @@
#ifndef EIGEN_TRIDIAGONALIZATION_H
#define EIGEN_TRIDIAGONALIZATION_H
+namespace internal {
+template<typename MatrixType, typename CoeffVectorType>
+void tridiagonalization_inplace(MatrixType& matA, CoeffVectorType& hCoeffs);
+}
+
/** \eigenvalues_module \ingroup Eigenvalues_Module
*
*
@@ -78,15 +83,15 @@ template<typename _MatrixType> class Tridiagonalization
};
typedef Matrix<Scalar, SizeMinusOne, 1, Options & ~RowMajor, MaxSizeMinusOne, 1> CoeffVectorType;
- typedef typename ei_plain_col_type<MatrixType, RealScalar>::type DiagonalType;
+ typedef typename internal::plain_col_type<MatrixType, RealScalar>::type DiagonalType;
typedef Matrix<RealScalar, SizeMinusOne, 1, Options & ~RowMajor, MaxSizeMinusOne, 1> SubDiagonalType;
- typedef typename ei_meta_if<NumTraits<Scalar>::IsComplex,
+ typedef typename internal::meta_if<NumTraits<Scalar>::IsComplex,
typename Diagonal<MatrixType,0>::RealReturnType,
Diagonal<MatrixType,0>
>::ret DiagonalReturnType;
- typedef typename ei_meta_if<NumTraits<Scalar>::IsComplex,
+ typedef typename internal::meta_if<NumTraits<Scalar>::IsComplex,
typename Diagonal<
Block<MatrixType,SizeMinusOne,SizeMinusOne>,0 >::RealReturnType,
Diagonal<
@@ -129,7 +134,7 @@ template<typename _MatrixType> class Tridiagonalization
m_hCoeffs(matrix.cols() > 1 ? matrix.cols()-1 : 1),
m_isInitialized(false)
{
- ei_tridiagonalization_inplace(m_matrix, m_hCoeffs);
+ internal::tridiagonalization_inplace(m_matrix, m_hCoeffs);
m_isInitialized = true;
}
@@ -154,7 +159,7 @@ template<typename _MatrixType> class Tridiagonalization
{
m_matrix = matrix;
m_hCoeffs.resize(matrix.rows()-1, 1);
- ei_tridiagonalization_inplace(m_matrix, m_hCoeffs);
+ internal::tridiagonalization_inplace(m_matrix, m_hCoeffs);
m_isInitialized = true;
return *this;
}
@@ -177,7 +182,7 @@ template<typename _MatrixType> class Tridiagonalization
*/
inline CoeffVectorType householderCoefficients() const
{
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
return m_hCoeffs;
}
@@ -214,7 +219,7 @@ template<typename _MatrixType> class Tridiagonalization
*/
inline const MatrixType& packedMatrix() const
{
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
return m_matrix;
}
@@ -235,7 +240,7 @@ template<typename _MatrixType> class Tridiagonalization
*/
HouseholderSequenceType matrixQ() const
{
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate(), false, m_matrix.rows() - 1, 1);
}
@@ -296,7 +301,7 @@ template<typename MatrixType>
const typename Tridiagonalization<MatrixType>::DiagonalReturnType
Tridiagonalization<MatrixType>::diagonal() const
{
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
return m_matrix.diagonal();
}
@@ -304,7 +309,7 @@ template<typename MatrixType>
const typename Tridiagonalization<MatrixType>::SubDiagonalReturnType
Tridiagonalization<MatrixType>::subDiagonal() const
{
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
Index n = m_matrix.rows();
return Block<MatrixType,SizeMinusOne,SizeMinusOne>(m_matrix, 1, 0, n-1,n-1).diagonal();
}
@@ -315,7 +320,7 @@ Tridiagonalization<MatrixType>::matrixT() const
{
// FIXME should this function (and other similar ones) rather take a matrix as argument
// and fill it ? (to avoid temporaries)
- ei_assert(m_isInitialized && "Tridiagonalization is not initialized.");
+ eigen_assert(m_isInitialized && "Tridiagonalization is not initialized.");
Index n = m_matrix.rows();
MatrixType matT = m_matrix;
matT.topRightCorner(n-1, n-1).diagonal() = subDiagonal().template cast<Scalar>().conjugate();
@@ -327,6 +332,8 @@ Tridiagonalization<MatrixType>::matrixT() const
return matT;
}
+namespace internal {
+
/** \internal
* Performs a tridiagonal decomposition of the selfadjoint matrix \a matA in-place.
*
@@ -351,10 +358,10 @@ Tridiagonalization<MatrixType>::matrixT() const
* \sa Tridiagonalization::packedMatrix()
*/
template<typename MatrixType, typename CoeffVectorType>
-void ei_tridiagonalization_inplace(MatrixType& matA, CoeffVectorType& hCoeffs)
+void tridiagonalization_inplace(MatrixType& matA, CoeffVectorType& hCoeffs)
{
- ei_assert(matA.rows()==matA.cols());
- ei_assert(matA.rows()==hCoeffs.size()+1);
+ eigen_assert(matA.rows()==matA.cols());
+ eigen_assert(matA.rows()==hCoeffs.size()+1);
typedef typename MatrixType::Index Index;
typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
@@ -371,9 +378,9 @@ void ei_tridiagonalization_inplace(MatrixType& matA, CoeffVectorType& hCoeffs)
matA.col(i).coeffRef(i+1) = 1;
hCoeffs.tail(n-i-1).noalias() = (matA.bottomRightCorner(remainingSize,remainingSize).template selfadjointView<Lower>()
- * (ei_conj(h) * matA.col(i).tail(remainingSize)));
+ * (conj(h) * matA.col(i).tail(remainingSize)));
- hCoeffs.tail(n-i-1) += (ei_conj(h)*Scalar(-0.5)*(hCoeffs.tail(remainingSize).dot(matA.col(i).tail(remainingSize)))) * matA.col(i).tail(n-i-1);
+ hCoeffs.tail(n-i-1) += (conj(h)*Scalar(-0.5)*(hCoeffs.tail(remainingSize).dot(matA.col(i).tail(remainingSize)))) * matA.col(i).tail(n-i-1);
matA.bottomRightCorner(remainingSize, remainingSize).template selfadjointView<Lower>()
.rankUpdate(matA.col(i).tail(remainingSize), hCoeffs.tail(remainingSize), -1);
@@ -387,7 +394,7 @@ void ei_tridiagonalization_inplace(MatrixType& matA, CoeffVectorType& hCoeffs)
template<typename MatrixType,
int Size=MatrixType::ColsAtCompileTime,
bool IsComplex=NumTraits<typename MatrixType::Scalar>::IsComplex>
-struct ei_tridiagonalization_inplace_selector;
+struct tridiagonalization_inplace_selector;
/** \brief Performs a full tridiagonalization in place
*
@@ -430,19 +437,19 @@ struct ei_tridiagonalization_inplace_selector;
* \sa class Tridiagonalization
*/
template<typename MatrixType, typename DiagonalType, typename SubDiagonalType>
-void ei_tridiagonalization_inplace(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
+void tridiagonalization_inplace(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
{
typedef typename MatrixType::Index Index;
//Index n = mat.rows();
- ei_assert(mat.cols()==mat.rows() && diag.size()==mat.rows() && subdiag.size()==mat.rows()-1);
- ei_tridiagonalization_inplace_selector<MatrixType>::run(mat, diag, subdiag, extractQ);
+ eigen_assert(mat.cols()==mat.rows() && diag.size()==mat.rows() && subdiag.size()==mat.rows()-1);
+ tridiagonalization_inplace_selector<MatrixType>::run(mat, diag, subdiag, extractQ);
}
/** \internal
* General full tridiagonalization
*/
template<typename MatrixType, int Size, bool IsComplex>
-struct ei_tridiagonalization_inplace_selector
+struct tridiagonalization_inplace_selector
{
typedef typename Tridiagonalization<MatrixType>::CoeffVectorType CoeffVectorType;
typedef typename Tridiagonalization<MatrixType>::HouseholderSequenceType HouseholderSequenceType;
@@ -451,7 +458,7 @@ struct ei_tridiagonalization_inplace_selector
static void run(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
{
CoeffVectorType hCoeffs(mat.cols()-1);
- ei_tridiagonalization_inplace(mat,hCoeffs);
+ tridiagonalization_inplace(mat,hCoeffs);
diag = mat.diagonal().real();
subdiag = mat.template diagonal<-1>().real();
if(extractQ)
@@ -464,7 +471,7 @@ struct ei_tridiagonalization_inplace_selector
* Especially useful for plane fitting.
*/
template<typename MatrixType>
-struct ei_tridiagonalization_inplace_selector<MatrixType,3,false>
+struct tridiagonalization_inplace_selector<MatrixType,3,false>
{
typedef typename MatrixType::Scalar Scalar;
typedef typename MatrixType::RealScalar RealScalar;
@@ -473,7 +480,7 @@ struct ei_tridiagonalization_inplace_selector<MatrixType,3,false>
static void run(MatrixType& mat, DiagonalType& diag, SubDiagonalType& subdiag, bool extractQ)
{
diag[0] = mat(0,0);
- RealScalar v1norm2 = ei_abs2(mat(2,0));
+ RealScalar v1norm2 = abs2(mat(2,0));
if(v1norm2 == RealScalar(0))
{
diag[1] = mat(1,1);
@@ -485,7 +492,7 @@ struct ei_tridiagonalization_inplace_selector<MatrixType,3,false>
}
else
{
- RealScalar beta = ei_sqrt(ei_abs2(mat(1,0)) + v1norm2);
+ RealScalar beta = sqrt(abs2(mat(1,0)) + v1norm2);
RealScalar invBeta = RealScalar(1)/beta;
Scalar m01 = mat(1,0) * invBeta;
Scalar m02 = mat(2,0) * invBeta;
@@ -508,16 +515,19 @@ struct ei_tridiagonalization_inplace_selector<MatrixType,3,false>
* Trivial specialization for 1x1 matrices
*/
template<typename MatrixType, bool IsComplex>
-struct ei_tridiagonalization_inplace_selector<MatrixType,1,IsComplex>
+struct tridiagonalization_inplace_selector<MatrixType,1,IsComplex>
{
typedef typename MatrixType::Scalar Scalar;
template<typename DiagonalType, typename SubDiagonalType>
static void run(MatrixType& mat, DiagonalType& diag, SubDiagonalType&, bool extractQ)
{
- diag(0,0) = ei_real(mat(0,0));
+ diag(0,0) = real(mat(0,0));
if(extractQ)
mat(0,0) = Scalar(1);
}
};
+
+} // end namespace internal
+
#endif // EIGEN_TRIDIAGONALIZATION_H