aboutsummaryrefslogtreecommitdiffhomepage
path: root/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
diff options
context:
space:
mode:
authorGravatar Gael Guennebaud <g.gael@free.fr>2015-05-07 15:55:12 +0200
committerGravatar Gael Guennebaud <g.gael@free.fr>2015-05-07 15:55:12 +0200
commit4a936974a5b895f61d0f80a61b27ea642037cb11 (patch)
tree311f4b521eecb64b95d54494b64d07ae2561d943 /Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
parentc2107d30ce9b9f30ff1f2d436667f3d09a4d9bd5 (diff)
bug #1013: fix 2x2 direct eigensolver for identical eiegnvalues
Diffstat (limited to 'Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h')
-rw-r--r--Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h30
1 files changed, 19 insertions, 11 deletions
diff --git a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
index 1dcfacf0b..8240b9465 100644
--- a/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
+++ b/Eigen/src/Eigenvalues/SelfAdjointEigenSolver.h
@@ -739,6 +739,7 @@ struct direct_selfadjoint_eigenvalues<SolverType,2,false>
static inline void run(SolverType& solver, const MatrixType& mat, int options)
{
EIGEN_USING_STD_MATH(sqrt);
+ EIGEN_USING_STD_MATH(abs);
eigen_assert(mat.cols() == 2 && mat.cols() == mat.rows());
eigen_assert((options&~(EigVecMask|GenEigMask))==0
@@ -760,22 +761,29 @@ struct direct_selfadjoint_eigenvalues<SolverType,2,false>
// compute the eigen vectors
if(computeEigenvectors)
{
- scaledMat.diagonal().array () -= eivals(1);
- Scalar a2 = numext::abs2(scaledMat(0,0));
- Scalar c2 = numext::abs2(scaledMat(1,1));
- Scalar b2 = numext::abs2(scaledMat(1,0));
- if(a2>c2)
+ if((eivals(1)-eivals(0))<=abs(eivals(1))*Eigen::NumTraits<Scalar>::epsilon())
{
- eivecs.col(1) << -scaledMat(1,0), scaledMat(0,0);
- eivecs.col(1) /= sqrt(a2+b2);
+ eivecs.setIdentity();
}
else
{
- eivecs.col(1) << -scaledMat(1,1), scaledMat(1,0);
- eivecs.col(1) /= sqrt(c2+b2);
- }
+ scaledMat.diagonal().array () -= eivals(1);
+ Scalar a2 = numext::abs2(scaledMat(0,0));
+ Scalar c2 = numext::abs2(scaledMat(1,1));
+ Scalar b2 = numext::abs2(scaledMat(1,0));
+ if(a2>c2)
+ {
+ eivecs.col(1) << -scaledMat(1,0), scaledMat(0,0);
+ eivecs.col(1) /= sqrt(a2+b2);
+ }
+ else
+ {
+ eivecs.col(1) << -scaledMat(1,1), scaledMat(1,0);
+ eivecs.col(1) /= sqrt(c2+b2);
+ }
- eivecs.col(0) << eivecs.col(1).unitOrthogonal();
+ eivecs.col(0) << eivecs.col(1).unitOrthogonal();
+ }
}
// Rescale back to the original size.