aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src/Eigenvalues
diff options
context:
space:
mode:
Diffstat (limited to 'Eigen/src/Eigenvalues')
-rw-r--r--Eigen/src/Eigenvalues/ComplexEigenSolver.h10
-rw-r--r--Eigen/src/Eigenvalues/ComplexSchur.h61
-rw-r--r--Eigen/src/Eigenvalues/EigenSolver.h58
-rw-r--r--Eigen/src/Eigenvalues/GeneralizedSelfAdjointEigenSolver.h4
-rw-r--r--Eigen/src/Eigenvalues/HessenbergDecomposition.h14
-rw-r--r--Eigen/src/Eigenvalues/MatrixBaseEigenvalues.h14
-rw-r--r--Eigen/src/Eigenvalues/RealSchur.h24
-rw-r--r--Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h42
-rw-r--r--Eigen/src/Eigenvalues/Tridiagonalization.h64
9 files changed, 159 insertions, 132 deletions
diff --git a/Eigen/src/Eigenvalues/ComplexEigenSolver.h b/Eigen/src/Eigenvalues/ComplexEigenSolver.h
index 7bf1d140e..d737d4cf0 100644
--- a/Eigen/src/Eigenvalues/ComplexEigenSolver.h
+++ b/Eigen/src/Eigenvalues/ComplexEigenSolver.h
@@ -169,8 +169,8 @@ template<typename _MatrixType> class ComplexEigenSolver
*/
const EigenvectorType& eigenvectors() const
{
- ei_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
return m_eivec;
}
@@ -193,7 +193,7 @@ template<typename _MatrixType> class ComplexEigenSolver
*/
const EigenvalueType& eigenvalues() const
{
- ei_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
return m_eivalues;
}
@@ -229,7 +229,7 @@ template<typename _MatrixType> class ComplexEigenSolver
*/
ComputationInfo info() const
{
- ei_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
return m_schur.info();
}
@@ -293,7 +293,7 @@ void ComplexEigenSolver<MatrixType>::doComputeEigenvectors(RealScalar matrixnorm
{
// If the i-th and k-th eigenvalue are equal, then z equals 0.
// Use a small value instead, to prevent division by zero.
- ei_real_ref(z) = NumTraits<RealScalar>::epsilon() * matrixnorm;
+ internal::real_ref(z) = NumTraits<RealScalar>::epsilon() * matrixnorm;
}
m_matX.coeffRef(i,k) = m_matX.coeff(i,k) / z;
}
diff --git a/Eigen/src/Eigenvalues/ComplexSchur.h b/Eigen/src/Eigenvalues/ComplexSchur.h
index 612573895..b1830f642 100644
--- a/Eigen/src/Eigenvalues/ComplexSchur.h
+++ b/Eigen/src/Eigenvalues/ComplexSchur.h
@@ -30,7 +30,9 @@
#include "./EigenvaluesCommon.h"
#include "./HessenbergDecomposition.h"
-template<typename MatrixType, bool IsComplex> struct ei_complex_schur_reduce_to_hessenberg;
+namespace internal {
+template<typename MatrixType, bool IsComplex> struct complex_schur_reduce_to_hessenberg;
+}
/** \eigenvalues_module \ingroup Eigenvalues_Module
*
@@ -146,8 +148,8 @@ template<typename _MatrixType> class ComplexSchur
*/
const ComplexMatrixType& matrixU() const
{
- ei_assert(m_isInitialized && "ComplexSchur is not initialized.");
- ei_assert(m_matUisUptodate && "The matrix U has not been computed during the ComplexSchur decomposition.");
+ eigen_assert(m_isInitialized && "ComplexSchur is not initialized.");
+ eigen_assert(m_matUisUptodate && "The matrix U has not been computed during the ComplexSchur decomposition.");
return m_matU;
}
@@ -170,7 +172,7 @@ template<typename _MatrixType> class ComplexSchur
*/
const ComplexMatrixType& matrixT() const
{
- ei_assert(m_isInitialized && "ComplexSchur is not initialized.");
+ eigen_assert(m_isInitialized && "ComplexSchur is not initialized.");
return m_matT;
}
@@ -201,7 +203,7 @@ template<typename _MatrixType> class ComplexSchur
*/
ComputationInfo info() const
{
- ei_assert(m_isInitialized && "RealSchur is not initialized.");
+ eigen_assert(m_isInitialized && "RealSchur is not initialized.");
return m_info;
}
@@ -222,22 +224,24 @@ template<typename _MatrixType> class ComplexSchur
bool subdiagonalEntryIsNeglegible(Index i);
ComplexScalar computeShift(Index iu, Index iter);
void reduceToTriangularForm(bool computeU);
- friend struct ei_complex_schur_reduce_to_hessenberg<MatrixType, NumTraits<Scalar>::IsComplex>;
+ friend struct internal::complex_schur_reduce_to_hessenberg<MatrixType, NumTraits<Scalar>::IsComplex>;
};
+namespace internal {
+
/** Computes the principal value of the square root of the complex \a z. */
template<typename RealScalar>
-std::complex<RealScalar> ei_sqrt(const std::complex<RealScalar> &z)
+std::complex<RealScalar> sqrt(const std::complex<RealScalar> &z)
{
RealScalar t, tre, tim;
- t = ei_abs(z);
+ t = abs(z);
- if (ei_abs(ei_real(z)) <= ei_abs(ei_imag(z)))
+ if (abs(real(z)) <= abs(imag(z)))
{
// No cancellation in these formulas
- tre = ei_sqrt(RealScalar(0.5)*(t + ei_real(z)));
- tim = ei_sqrt(RealScalar(0.5)*(t - ei_real(z)));
+ tre = sqrt(RealScalar(0.5)*(t + real(z)));
+ tim = sqrt(RealScalar(0.5)*(t - real(z)));
}
else
{
@@ -245,14 +249,14 @@ std::complex<RealScalar> ei_sqrt(const std::complex<RealScalar> &z)
if (z.real() > RealScalar(0))
{
tre = t + z.real();
- tim = ei_abs(ei_imag(z))*ei_sqrt(RealScalar(0.5)/tre);
- tre = ei_sqrt(RealScalar(0.5)*tre);
+ tim = abs(imag(z))*sqrt(RealScalar(0.5)/tre);
+ tre = sqrt(RealScalar(0.5)*tre);
}
else
{
tim = t - z.real();
- tre = ei_abs(ei_imag(z))*ei_sqrt(RealScalar(0.5)/tim);
- tim = ei_sqrt(RealScalar(0.5)*tim);
+ tre = abs(imag(z))*sqrt(RealScalar(0.5)/tim);
+ tim = sqrt(RealScalar(0.5)*tim);
}
}
if(z.imag() < RealScalar(0))
@@ -260,6 +264,7 @@ std::complex<RealScalar> ei_sqrt(const std::complex<RealScalar> &z)
return (std::complex<RealScalar>(tre,tim));
}
+} // end namespace internal
/** If m_matT(i+1,i) is neglegible in floating point arithmetic
@@ -268,9 +273,9 @@ std::complex<RealScalar> ei_sqrt(const std::complex<RealScalar> &z)
template<typename MatrixType>
inline bool ComplexSchur<MatrixType>::subdiagonalEntryIsNeglegible(Index i)
{
- RealScalar d = ei_norm1(m_matT.coeff(i,i)) + ei_norm1(m_matT.coeff(i+1,i+1));
- RealScalar sd = ei_norm1(m_matT.coeff(i+1,i));
- if (ei_isMuchSmallerThan(sd, d, NumTraits<RealScalar>::epsilon()))
+ RealScalar d = internal::norm1(m_matT.coeff(i,i)) + internal::norm1(m_matT.coeff(i+1,i+1));
+ RealScalar sd = internal::norm1(m_matT.coeff(i+1,i));
+ if (internal::isMuchSmallerThan(sd, d, NumTraits<RealScalar>::epsilon()))
{
m_matT.coeffRef(i+1,i) = ComplexScalar(0);
return true;
@@ -286,7 +291,7 @@ typename ComplexSchur<MatrixType>::ComplexScalar ComplexSchur<MatrixType>::compu
if (iter == 10 || iter == 20)
{
// exceptional shift, taken from http://www.netlib.org/eispack/comqr.f
- return ei_abs(ei_real(m_matT.coeff(iu,iu-1))) + ei_abs(ei_real(m_matT.coeff(iu-1,iu-2)));
+ return internal::abs(internal::real(m_matT.coeff(iu,iu-1))) + internal::abs(internal::real(m_matT.coeff(iu-1,iu-2)));
}
// compute the shift as one of the eigenvalues of t, the 2x2
@@ -297,19 +302,19 @@ typename ComplexSchur<MatrixType>::ComplexScalar ComplexSchur<MatrixType>::compu
ComplexScalar b = t.coeff(0,1) * t.coeff(1,0);
ComplexScalar c = t.coeff(0,0) - t.coeff(1,1);
- ComplexScalar disc = ei_sqrt(c*c + RealScalar(4)*b);
+ ComplexScalar disc = internal::sqrt(c*c + RealScalar(4)*b);
ComplexScalar det = t.coeff(0,0) * t.coeff(1,1) - b;
ComplexScalar trace = t.coeff(0,0) + t.coeff(1,1);
ComplexScalar eival1 = (trace + disc) / RealScalar(2);
ComplexScalar eival2 = (trace - disc) / RealScalar(2);
- if(ei_norm1(eival1) > ei_norm1(eival2))
+ if(internal::norm1(eival1) > internal::norm1(eival2))
eival2 = det / eival1;
else
eival1 = det / eival2;
// choose the eigenvalue closest to the bottom entry of the diagonal
- if(ei_norm1(eival1-t.coeff(1,1)) < ei_norm1(eival2-t.coeff(1,1)))
+ if(internal::norm1(eival1-t.coeff(1,1)) < internal::norm1(eival2-t.coeff(1,1)))
return normt * eival1;
else
return normt * eival2;
@@ -320,7 +325,7 @@ template<typename MatrixType>
ComplexSchur<MatrixType>& ComplexSchur<MatrixType>::compute(const MatrixType& matrix, bool computeU)
{
m_matUisUptodate = false;
- ei_assert(matrix.cols() == matrix.rows());
+ eigen_assert(matrix.cols() == matrix.rows());
if(matrix.cols() == 1)
{
@@ -332,14 +337,16 @@ ComplexSchur<MatrixType>& ComplexSchur<MatrixType>::compute(const MatrixType& ma
return *this;
}
- ei_complex_schur_reduce_to_hessenberg<MatrixType, NumTraits<Scalar>::IsComplex>::run(*this, matrix, computeU);
+ internal::complex_schur_reduce_to_hessenberg<MatrixType, NumTraits<Scalar>::IsComplex>::run(*this, matrix, computeU);
reduceToTriangularForm(computeU);
return *this;
}
+namespace internal {
+
/* Reduce given matrix to Hessenberg form */
template<typename MatrixType, bool IsComplex>
-struct ei_complex_schur_reduce_to_hessenberg
+struct complex_schur_reduce_to_hessenberg
{
// this is the implementation for the case IsComplex = true
static void run(ComplexSchur<MatrixType>& _this, const MatrixType& matrix, bool computeU)
@@ -351,7 +358,7 @@ struct ei_complex_schur_reduce_to_hessenberg
};
template<typename MatrixType>
-struct ei_complex_schur_reduce_to_hessenberg<MatrixType, false>
+struct complex_schur_reduce_to_hessenberg<MatrixType, false>
{
static void run(ComplexSchur<MatrixType>& _this, const MatrixType& matrix, bool computeU)
{
@@ -370,6 +377,8 @@ struct ei_complex_schur_reduce_to_hessenberg<MatrixType, false>
}
};
+} // end namespace internal
+
// Reduce the Hessenberg matrix m_matT to triangular form by QR iteration.
template<typename MatrixType>
void ComplexSchur<MatrixType>::reduceToTriangularForm(bool computeU)
diff --git a/Eigen/src/Eigenvalues/EigenSolver.h b/Eigen/src/Eigenvalues/EigenSolver.h
index 0a5faec52..5a59ccbd4 100644
--- a/Eigen/src/Eigenvalues/EigenSolver.h
+++ b/Eigen/src/Eigenvalues/EigenSolver.h
@@ -211,8 +211,8 @@ template<typename _MatrixType> class EigenSolver
*/
const MatrixType& pseudoEigenvectors() const
{
- ei_assert(m_isInitialized && "EigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
return m_eivec;
}
@@ -254,7 +254,7 @@ template<typename _MatrixType> class EigenSolver
*/
const EigenvalueType& eigenvalues() const
{
- ei_assert(m_isInitialized && "EigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
return m_eivalues;
}
@@ -289,7 +289,7 @@ template<typename _MatrixType> class EigenSolver
ComputationInfo info() const
{
- ei_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "ComplexEigenSolver is not initialized.");
return m_realSchur.info();
}
@@ -311,17 +311,17 @@ template<typename _MatrixType> class EigenSolver
template<typename MatrixType>
MatrixType EigenSolver<MatrixType>::pseudoEigenvalueMatrix() const
{
- ei_assert(m_isInitialized && "EigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
Index n = m_eivalues.rows();
MatrixType matD = MatrixType::Zero(n,n);
for (Index i=0; i<n; ++i)
{
- if (ei_isMuchSmallerThan(ei_imag(m_eivalues.coeff(i)), ei_real(m_eivalues.coeff(i))))
- matD.coeffRef(i,i) = ei_real(m_eivalues.coeff(i));
+ if (internal::isMuchSmallerThan(internal::imag(m_eivalues.coeff(i)), internal::real(m_eivalues.coeff(i))))
+ matD.coeffRef(i,i) = internal::real(m_eivalues.coeff(i));
else
{
- matD.template block<2,2>(i,i) << ei_real(m_eivalues.coeff(i)), ei_imag(m_eivalues.coeff(i)),
- -ei_imag(m_eivalues.coeff(i)), ei_real(m_eivalues.coeff(i));
+ matD.template block<2,2>(i,i) << internal::real(m_eivalues.coeff(i)), internal::imag(m_eivalues.coeff(i)),
+ -internal::imag(m_eivalues.coeff(i)), internal::real(m_eivalues.coeff(i));
++i;
}
}
@@ -331,13 +331,13 @@ MatrixType EigenSolver<MatrixType>::pseudoEigenvalueMatrix() const
template<typename MatrixType>
typename EigenSolver<MatrixType>::EigenvectorsType EigenSolver<MatrixType>::eigenvectors() const
{
- ei_assert(m_isInitialized && "EigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "EigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
Index n = m_eivec.cols();
EigenvectorsType matV(n,n);
for (Index j=0; j<n; ++j)
{
- if (ei_isMuchSmallerThan(ei_imag(m_eivalues.coeff(j)), ei_real(m_eivalues.coeff(j))))
+ if (internal::isMuchSmallerThan(internal::imag(m_eivalues.coeff(j)), internal::real(m_eivalues.coeff(j))))
{
// we have a real eigen value
matV.col(j) = m_eivec.col(j).template cast<ComplexScalar>();
@@ -384,7 +384,7 @@ EigenSolver<MatrixType>& EigenSolver<MatrixType>::compute(const MatrixType& matr
else
{
Scalar p = Scalar(0.5) * (m_matT.coeff(i, i) - m_matT.coeff(i+1, i+1));
- Scalar z = ei_sqrt(ei_abs(p * p + m_matT.coeff(i+1, i) * m_matT.coeff(i, i+1)));
+ Scalar z = internal::sqrt(internal::abs(p * p + m_matT.coeff(i+1, i) * m_matT.coeff(i, i+1)));
m_eivalues.coeffRef(i) = ComplexScalar(m_matT.coeff(i+1, i+1) + p, z);
m_eivalues.coeffRef(i+1) = ComplexScalar(m_matT.coeff(i+1, i+1) + p, -z);
i += 2;
@@ -407,7 +407,7 @@ template<typename Scalar>
std::complex<Scalar> cdiv(Scalar xr, Scalar xi, Scalar yr, Scalar yi)
{
Scalar r,d;
- if (ei_abs(yr) > ei_abs(yi))
+ if (internal::abs(yr) > internal::abs(yi))
{
r = yi/yr;
d = yr + r*yi;
@@ -480,14 +480,14 @@ void EigenSolver<MatrixType>::doComputeEigenvectors()
Scalar denom = (m_eivalues.coeff(i).real() - p) * (m_eivalues.coeff(i).real() - p) + m_eivalues.coeff(i).imag() * m_eivalues.coeff(i).imag();
Scalar t = (x * lastr - lastw * r) / denom;
m_matT.coeffRef(i,n) = t;
- if (ei_abs(x) > ei_abs(lastw))
+ if (internal::abs(x) > internal::abs(lastw))
m_matT.coeffRef(i+1,n) = (-r - w * t) / x;
else
m_matT.coeffRef(i+1,n) = (-lastr - y * t) / lastw;
}
// Overflow control
- Scalar t = ei_abs(m_matT.coeff(i,n));
+ Scalar t = internal::abs(m_matT.coeff(i,n));
if ((eps * t) * t > 1)
m_matT.col(n).tail(size-i) /= t;
}
@@ -499,16 +499,16 @@ void EigenSolver<MatrixType>::doComputeEigenvectors()
Index l = n-1;
// Last vector component imaginary so matrix is triangular
- if (ei_abs(m_matT.coeff(n,n-1)) > ei_abs(m_matT.coeff(n-1,n)))
+ if (internal::abs(m_matT.coeff(n,n-1)) > internal::abs(m_matT.coeff(n-1,n)))
{
m_matT.coeffRef(n-1,n-1) = q / m_matT.coeff(n,n-1);
m_matT.coeffRef(n-1,n) = -(m_matT.coeff(n,n) - p) / m_matT.coeff(n,n-1);
}
else
{
- std::complex<Scalar> cc = cdiv<Scalar>(0.0,-m_matT.coeff(n-1,n),m_matT.coeff(n-1,n-1)-p,q);
- m_matT.coeffRef(n-1,n-1) = ei_real(cc);
- m_matT.coeffRef(n-1,n) = ei_imag(cc);
+ std::complex<Scalar> cc = cdiv<Scalar>(0.0,-m_matT.coeff(n-1,n),m_matT.coeff(n-1,n-1)-p,q);
+ m_matT.coeffRef(n-1,n-1) = internal::real(cc);
+ m_matT.coeffRef(n-1,n) = internal::imag(cc);
}
m_matT.coeffRef(n,n-1) = 0.0;
m_matT.coeffRef(n,n) = 1.0;
@@ -530,8 +530,8 @@ void EigenSolver<MatrixType>::doComputeEigenvectors()
if (m_eivalues.coeff(i).imag() == 0)
{
std::complex<Scalar> cc = cdiv(-ra,-sa,w,q);
- m_matT.coeffRef(i,n-1) = ei_real(cc);
- m_matT.coeffRef(i,n) = ei_imag(cc);
+ m_matT.coeffRef(i,n-1) = internal::real(cc);
+ m_matT.coeffRef(i,n) = internal::imag(cc);
}
else
{
@@ -541,12 +541,12 @@ void EigenSolver<MatrixType>::doComputeEigenvectors()
Scalar vr = (m_eivalues.coeff(i).real() - p) * (m_eivalues.coeff(i).real() - p) + m_eivalues.coeff(i).imag() * m_eivalues.coeff(i).imag() - q * q;
Scalar vi = (m_eivalues.coeff(i).real() - p) * Scalar(2) * q;
if ((vr == 0.0) && (vi == 0.0))
- vr = eps * norm * (ei_abs(w) + ei_abs(q) + ei_abs(x) + ei_abs(y) + ei_abs(lastw));
+ vr = eps * norm * (internal::abs(w) + internal::abs(q) + internal::abs(x) + internal::abs(y) + internal::abs(lastw));
std::complex<Scalar> cc = cdiv(x*lastra-lastw*ra+q*sa,x*lastsa-lastw*sa-q*ra,vr,vi);
- m_matT.coeffRef(i,n-1) = ei_real(cc);
- m_matT.coeffRef(i,n) = ei_imag(cc);
- if (ei_abs(x) > (ei_abs(lastw) + ei_abs(q)))
+ m_matT.coeffRef(i,n-1) = internal::real(cc);
+ m_matT.coeffRef(i,n) = internal::imag(cc);
+ if (internal::abs(x) > (internal::abs(lastw) + internal::abs(q)))
{
m_matT.coeffRef(i+1,n-1) = (-ra - w * m_matT.coeff(i,n-1) + q * m_matT.coeff(i,n)) / x;
m_matT.coeffRef(i+1,n) = (-sa - w * m_matT.coeff(i,n) - q * m_matT.coeff(i,n-1)) / x;
@@ -554,13 +554,13 @@ void EigenSolver<MatrixType>::doComputeEigenvectors()
else
{
cc = cdiv(-lastra-y*m_matT.coeff(i,n-1),-lastsa-y*m_matT.coeff(i,n),lastw,q);
- m_matT.coeffRef(i+1,n-1) = ei_real(cc);
- m_matT.coeffRef(i+1,n) = ei_imag(cc);
+ m_matT.coeffRef(i+1,n-1) = internal::real(cc);
+ m_matT.coeffRef(i+1,n) = internal::imag(cc);
}
}
// Overflow control
- Scalar t = std::max(ei_abs(m_matT.coeff(i,n-1)),ei_abs(m_matT.coeff(i,n)));
+ Scalar t = std::max(internal::abs(m_matT.coeff(i,n-1)),internal::abs(m_matT.coeff(i,n)));
if ((eps * t) * t > 1)
m_matT.block(i, n-1, size-i, 2) /= t;
diff --git a/Eigen/src/Eigenvalues/GeneralizedSelfAdjointEigenSolver.h b/Eigen/src/Eigenvalues/GeneralizedSelfAdjointEigenSolver.h
index a47a3dc2e..a19e8cf24 100644
--- a/Eigen/src/Eigenvalues/GeneralizedSelfAdjointEigenSolver.h
+++ b/Eigen/src/Eigenvalues/GeneralizedSelfAdjointEigenSolver.h
@@ -182,8 +182,8 @@ template<typename MatrixType>
GeneralizedSelfAdjointEigenSolver<MatrixType>& GeneralizedSelfAdjointEigenSolver<MatrixType>::
compute(const MatrixType& matA, const MatrixType& matB, int options)
{
- ei_assert(matA.cols()==matA.rows() && matB.rows()==matA.rows() && matB.cols()==matB.rows());
- ei_assert((options&~(EigVecMask|GenEigMask))==0
+ eigen_assert(matA.cols()==matA.rows() && matB.rows()==matA.rows() && matB.cols()==matB.rows());
+ eigen_assert((options&~(EigVecMask|GenEigMask))==0
&& (options&EigVecMask)!=EigVecMask
&& ((options&GenEigMask)==0 || (options&GenEigMask)==Ax_lBx
|| (options&GenEigMask)==ABx_lx || (options&GenEigMask)==BAx_lx)
diff --git a/Eigen/src/Eigenvalues/HessenbergDecomposition.h b/Eigen/src/Eigenvalues/HessenbergDecomposition.h
index 79554187a..b7da136a8 100644
--- a/Eigen/src/Eigenvalues/HessenbergDecomposition.h
+++ b/Eigen/src/Eigenvalues/HessenbergDecomposition.h
@@ -28,11 +28,13 @@
template<typename MatrixType> struct HessenbergDecompositionMatrixHReturnType;
+namespace internal {
template<typename MatrixType>
-struct ei_traits<HessenbergDecompositionMatrixHReturnType<MatrixType> >
+struct traits<HessenbergDecompositionMatrixHReturnType<MatrixType> >
{
typedef MatrixType ReturnType;
};
+}
/** \eigenvalues_module \ingroup Eigenvalues_Module
*
@@ -184,7 +186,7 @@ template<typename _MatrixType> class HessenbergDecomposition
*/
const CoeffVectorType& householderCoefficients() const
{
- ei_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
+ eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
return m_hCoeffs;
}
@@ -219,7 +221,7 @@ template<typename _MatrixType> class HessenbergDecomposition
*/
const MatrixType& packedMatrix() const
{
- ei_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
+ eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
return m_matrix;
}
@@ -239,7 +241,7 @@ template<typename _MatrixType> class HessenbergDecomposition
*/
HouseholderSequenceType matrixQ() const
{
- ei_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
+ eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
return HouseholderSequenceType(m_matrix, m_hCoeffs.conjugate(), false, m_matrix.rows() - 1, 1);
}
@@ -265,7 +267,7 @@ template<typename _MatrixType> class HessenbergDecomposition
*/
HessenbergDecompositionMatrixHReturnType<MatrixType> matrixH() const
{
- ei_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
+ eigen_assert(m_isInitialized && "HessenbergDecomposition is not initialized.");
return HessenbergDecompositionMatrixHReturnType<MatrixType>(*this);
}
@@ -319,7 +321,7 @@ void HessenbergDecomposition<MatrixType>::_compute(MatrixType& matA, CoeffVector
// A = A H'
matA.rightCols(remainingSize)
- .applyHouseholderOnTheRight(matA.col(i).tail(remainingSize-1).conjugate(), ei_conj(h), &temp.coeffRef(0));
+ .applyHouseholderOnTheRight(matA.col(i).tail(remainingSize-1).conjugate(), internal::conj(h), &temp.coeffRef(0));
}
}
diff --git a/Eigen/src/Eigenvalues/MatrixBaseEigenvalues.h b/Eigen/src/Eigenvalues/MatrixBaseEigenvalues.h
index e517b6e5a..5591519fb 100644
--- a/Eigen/src/Eigenvalues/MatrixBaseEigenvalues.h
+++ b/Eigen/src/Eigenvalues/MatrixBaseEigenvalues.h
@@ -26,10 +26,10 @@
#ifndef EIGEN_MATRIXBASEEIGENVALUES_H
#define EIGEN_MATRIXBASEEIGENVALUES_H
-
+namespace internal {
template<typename Derived, bool IsComplex>
-struct ei_eigenvalues_selector
+struct eigenvalues_selector
{
// this is the implementation for the case IsComplex = true
static inline typename MatrixBase<Derived>::EigenvaluesReturnType const
@@ -42,7 +42,7 @@ struct ei_eigenvalues_selector
};
template<typename Derived>
-struct ei_eigenvalues_selector<Derived, false>
+struct eigenvalues_selector<Derived, false>
{
static inline typename MatrixBase<Derived>::EigenvaluesReturnType const
run(const MatrixBase<Derived>& m)
@@ -53,6 +53,8 @@ struct ei_eigenvalues_selector<Derived, false>
}
};
+} // end namespace internal
+
/** \brief Computes the eigenvalues of a matrix
* \returns Column vector containing the eigenvalues.
*
@@ -77,8 +79,8 @@ template<typename Derived>
inline typename MatrixBase<Derived>::EigenvaluesReturnType
MatrixBase<Derived>::eigenvalues() const
{
- typedef typename ei_traits<Derived>::Scalar Scalar;
- return ei_eigenvalues_selector<Derived, NumTraits<Scalar>::IsComplex>::run(derived());
+ typedef typename internal::traits<Derived>::Scalar Scalar;
+ return internal::eigenvalues_selector<Derived, NumTraits<Scalar>::IsComplex>::run(derived());
}
/** \brief Computes the eigenvalues of a matrix
@@ -135,7 +137,7 @@ MatrixBase<Derived>::operatorNorm() const
typename Derived::PlainObject m_eval(derived());
// FIXME if it is really guaranteed that the eigenvalues are already sorted,
// then we don't need to compute a maxCoeff() here, comparing the 1st and last ones is enough.
- return ei_sqrt((m_eval*m_eval.adjoint())
+ return internal::sqrt((m_eval*m_eval.adjoint())
.eval()
.template selfadjointView<Lower>()
.eigenvalues()
diff --git a/Eigen/src/Eigenvalues/RealSchur.h b/Eigen/src/Eigenvalues/RealSchur.h
index 1d299d5f6..e250d9a86 100644
--- a/Eigen/src/Eigenvalues/RealSchur.h
+++ b/Eigen/src/Eigenvalues/RealSchur.h
@@ -137,8 +137,8 @@ template<typename _MatrixType> class RealSchur
*/
const MatrixType& matrixU() const
{
- ei_assert(m_isInitialized && "RealSchur is not initialized.");
- ei_assert(m_matUisUptodate && "The matrix U has not been computed during the RealSchur decomposition.");
+ eigen_assert(m_isInitialized && "RealSchur is not initialized.");
+ eigen_assert(m_matUisUptodate && "The matrix U has not been computed during the RealSchur decomposition.");
return m_matU;
}
@@ -154,7 +154,7 @@ template<typename _MatrixType> class RealSchur
*/
const MatrixType& matrixT() const
{
- ei_assert(m_isInitialized && "RealSchur is not initialized.");
+ eigen_assert(m_isInitialized && "RealSchur is not initialized.");
return m_matT;
}
@@ -183,7 +183,7 @@ template<typename _MatrixType> class RealSchur
*/
ComputationInfo info() const
{
- ei_assert(m_isInitialized && "RealSchur is not initialized.");
+ eigen_assert(m_isInitialized && "RealSchur is not initialized.");
return m_info;
}
@@ -300,10 +300,10 @@ inline typename MatrixType::Index RealSchur<MatrixType>::findSmallSubdiagEntry(I
Index res = iu;
while (res > 0)
{
- Scalar s = ei_abs(m_matT.coeff(res-1,res-1)) + ei_abs(m_matT.coeff(res,res));
+ Scalar s = internal::abs(m_matT.coeff(res-1,res-1)) + internal::abs(m_matT.coeff(res,res));
if (s == 0.0)
s = norm;
- if (ei_abs(m_matT.coeff(res,res-1)) < NumTraits<Scalar>::epsilon() * s)
+ if (internal::abs(m_matT.coeff(res,res-1)) < NumTraits<Scalar>::epsilon() * s)
break;
res--;
}
@@ -325,7 +325,7 @@ inline void RealSchur<MatrixType>::splitOffTwoRows(Index iu, bool computeU, Scal
if (q >= 0) // Two real eigenvalues
{
- Scalar z = ei_sqrt(ei_abs(q));
+ Scalar z = internal::sqrt(internal::abs(q));
JacobiRotation<Scalar> rot;
if (p >= 0)
rot.makeGivens(p + z, m_matT.coeff(iu, iu-1));
@@ -357,7 +357,7 @@ inline void RealSchur<MatrixType>::computeShift(Index iu, Index iter, Scalar& ex
exshift += shiftInfo.coeff(0);
for (Index i = 0; i <= iu; ++i)
m_matT.coeffRef(i,i) -= shiftInfo.coeff(0);
- Scalar s = ei_abs(m_matT.coeff(iu,iu-1)) + ei_abs(m_matT.coeff(iu-1,iu-2));
+ Scalar s = internal::abs(m_matT.coeff(iu,iu-1)) + internal::abs(m_matT.coeff(iu-1,iu-2));
shiftInfo.coeffRef(0) = Scalar(0.75) * s;
shiftInfo.coeffRef(1) = Scalar(0.75) * s;
shiftInfo.coeffRef(2) = Scalar(-0.4375) * s * s;
@@ -370,7 +370,7 @@ inline void RealSchur<MatrixType>::computeShift(Index iu, Index iter, Scalar& ex
s = s * s + shiftInfo.coeff(2);
if (s > 0)
{
- s = ei_sqrt(s);
+ s = internal::sqrt(s);
if (shiftInfo.coeff(1) < shiftInfo.coeff(0))
s = -s;
s = s + (shiftInfo.coeff(1) - shiftInfo.coeff(0)) / Scalar(2.0);
@@ -400,9 +400,9 @@ inline void RealSchur<MatrixType>::initFrancisQRStep(Index il, Index iu, const V
if (im == il) {
break;
}
- const Scalar lhs = m_matT.coeff(im,im-1) * (ei_abs(v.coeff(1)) + ei_abs(v.coeff(2)));
- const Scalar rhs = v.coeff(0) * (ei_abs(m_matT.coeff(im-1,im-1)) + ei_abs(Tmm) + ei_abs(m_matT.coeff(im+1,im+1)));
- if (ei_abs(lhs) < NumTraits<Scalar>::epsilon() * rhs)
+ const Scalar lhs = m_matT.coeff(im,im-1) * (internal::abs(v.coeff(1)) + internal::abs(v.coeff(2)));
+ const Scalar rhs = v.coeff(0) * (internal::abs(m_matT.coeff(im-1,im-1)) + internal::abs(Tmm) + internal::abs(m_matT.coeff(im+1,im+1)));
+ if (internal::abs(lhs) < NumTraits<Scalar>::epsilon() * rhs)
{
break;
}
diff --git a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
index d32e223fb..7cef9e529 100644
--- a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
+++ b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
@@ -101,7 +101,7 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
* This is a column vector with entries of type #RealScalar.
* The length of the vector is the size of \p _MatrixType.
*/
- typedef typename ei_plain_col_type<MatrixType, RealScalar>::type RealVectorType;
+ typedef typename internal::plain_col_type<MatrixType, RealScalar>::type RealVectorType;
typedef Tridiagonalization<MatrixType> TridiagonalizationType;
/** \brief Default constructor for fixed-size matrices.
@@ -219,8 +219,8 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
*/
const MatrixType& eigenvectors() const
{
- ei_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
return m_eivec;
}
@@ -240,7 +240,7 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
*/
const RealVectorType& eigenvalues() const
{
- ei_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
return m_eivalues;
}
@@ -264,8 +264,8 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
*/
MatrixType operatorSqrt() const
{
- ei_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
return m_eivec * m_eivalues.cwiseSqrt().asDiagonal() * m_eivec.adjoint();
}
@@ -289,8 +289,8 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
*/
MatrixType operatorInverseSqrt() const
{
- ei_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
- ei_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
+ eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+ eigen_assert(m_eigenvectorsOk && "The eigenvectors have not been computed together with the eigenvalues.");
return m_eivec * m_eivalues.cwiseInverse().cwiseSqrt().asDiagonal() * m_eivec.adjoint();
}
@@ -300,7 +300,7 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
*/
ComputationInfo info() const
{
- ei_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
+ eigen_assert(m_isInitialized && "SelfAdjointEigenSolver is not initialized.");
return m_info;
}
@@ -335,15 +335,17 @@ template<typename _MatrixType> class SelfAdjointEigenSolver
* Implemented from Golub's "Matrix Computations", algorithm 8.3.2:
* "implicit symmetric QR step with Wilkinson shift"
*/
+namespace internal {
template<typename RealScalar, typename Scalar, typename Index>
-static void ei_tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n);
+static void tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n);
+}
template<typename MatrixType>
SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
::compute(const MatrixType& matrix, int options)
{
- ei_assert(matrix.cols() == matrix.rows());
- ei_assert((options&~(EigVecMask|GenEigMask))==0
+ eigen_assert(matrix.cols() == matrix.rows());
+ eigen_assert((options&~(EigVecMask|GenEigMask))==0
&& (options&EigVecMask)!=EigVecMask
&& "invalid option parameter");
bool computeEigenvectors = (options&ComputeEigenvectors)==ComputeEigenvectors;
@@ -352,7 +354,7 @@ SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
if(n==1)
{
- m_eivalues.coeffRef(0,0) = ei_real(matrix.coeff(0,0));
+ m_eivalues.coeffRef(0,0) = internal::real(matrix.coeff(0,0));
if(computeEigenvectors)
m_eivec.setOnes();
m_info = Success;
@@ -367,7 +369,7 @@ SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
mat = matrix;
m_subdiag.resize(n-1);
- ei_tridiagonalization_inplace(mat, diag, m_subdiag, computeEigenvectors);
+ internal::tridiagonalization_inplace(mat, diag, m_subdiag, computeEigenvectors);
Index end = n-1;
Index start = 0;
@@ -376,7 +378,7 @@ SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
while (end>0)
{
for (Index i = start; i<end; ++i)
- if (ei_isMuchSmallerThan(ei_abs(m_subdiag[i]),(ei_abs(diag[i])+ei_abs(diag[i+1]))))
+ if (internal::isMuchSmallerThan(internal::abs(m_subdiag[i]),(internal::abs(diag[i])+internal::abs(diag[i+1]))))
m_subdiag[i] = 0;
// find the largest unreduced block
@@ -396,7 +398,7 @@ SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
while (start>0 && m_subdiag[start-1]!=0)
start--;
- ei_tridiagonal_qr_step(diag.data(), m_subdiag.data(), start, end, computeEigenvectors ? m_eivec.data() : (Scalar*)0, n);
+ internal::tridiagonal_qr_step(diag.data(), m_subdiag.data(), start, end, computeEigenvectors ? m_eivec.data() : (Scalar*)0, n);
}
if (iter <= m_maxIterations)
@@ -427,12 +429,13 @@ SelfAdjointEigenSolver<MatrixType>& SelfAdjointEigenSolver<MatrixType>
return *this;
}
+namespace internal {
template<typename RealScalar, typename Scalar, typename Index>
-static void ei_tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n)
+static void tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index start, Index end, Scalar* matrixQ, Index n)
{
RealScalar td = (diag[end-1] - diag[end])*RealScalar(0.5);
- RealScalar e2 = ei_abs2(subdiag[end-1]);
- RealScalar mu = diag[end] - e2 / (td + (td>0 ? 1 : -1) * ei_sqrt(td*td + e2));
+ RealScalar e2 = abs2(subdiag[end-1]);
+ RealScalar mu = diag[end] - e2 / (td + (td>0 ? 1 : -1) * sqrt(td*td + e2));
RealScalar x = diag[start] - mu;
RealScalar z = subdiag[start];
@@ -468,5 +471,6 @@ static void ei_tridiagonal_qr_step(RealScalar* diag, RealScalar* subdiag, Index
}
}
}
+} // end namespace internal
#endif // EIGEN_SELFADJOINTEIGENSOLVER_H
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
generated by cgit v1.2.3 (git 2.39.1) at 2024-06-20 19:06:10 +0000