diff --git a/unsupported/test/matrix_functions.h b/unsupported/test/matrix_functions.h index 5817caef6..295da16b6 100644 --- a/unsupported/test/matrix_functions.h +++ b/unsupported/test/matrix_functions.h @@ -10,27 +10,48 @@ #include "main.h" #include +// For complex matrices, any matrix is fine. +template::Scalar>::IsComplex> +struct processTriangularMatrix +{ + static void run(MatrixType&, MatrixType&, const MatrixType&) + { } +}; + +// For real matrices, make sure none of the eigenvalues are negative. +template +struct processTriangularMatrix +{ + static void run(MatrixType& m, MatrixType& T, const MatrixType& U) + { + typedef typename MatrixType::Index Index; + const Index size = m.cols(); + + for (Index i=0; i < size; ++i) { + if (i == size - 1 || T.coeff(i+1,i) == 0) + T.coeffRef(i,i) = std::abs(T.coeff(i,i)); + else + ++i; + } + m = U * T * U.transpose(); + } +}; + template ::Scalar>::IsComplex> struct generateTestMatrix; -// for real matrices, make sure none of the eigenvalues are negative template struct generateTestMatrix { static void run(MatrixType& result, typename MatrixType::Index size) { - MatrixType mat = MatrixType::Random(size, size); - EigenSolver es(mat); - typename EigenSolver::EigenvalueType eivals = es.eigenvalues(); - for (typename MatrixType::Index i = 0; i < size; ++i) { - if (eivals(i).imag() == 0 && eivals(i).real() < 0) - eivals(i) = -eivals(i); - } - result = (es.eigenvectors() * eivals.asDiagonal() * es.eigenvectors().inverse()).real(); + result = MatrixType::Random(size, size); + RealSchur schur(result); + MatrixType T = schur.matrixT(); + processTriangularMatrix::run(result, T, schur.matrixU()); } }; -// for complex matrices, any matrix is fine template struct generateTestMatrix { diff --git a/unsupported/test/matrix_power.cpp b/unsupported/test/matrix_power.cpp index 3ee19fc56..849e4287b 100644 --- a/unsupported/test/matrix_power.cpp +++ b/unsupported/test/matrix_power.cpp @@ -96,33 +96,6 @@ void testGeneral(const MatrixType& m, double tol) } } -// For complex matrices, any matrix is fine. -template::Scalar>::IsComplex> -struct processTriangularMatrix -{ - static void run(MatrixType&, MatrixType&, const MatrixType&) - { } -}; - -// For real matrices, make sure none of the eigenvalues are negative. -template -struct processTriangularMatrix -{ - static void run(MatrixType& m, MatrixType& T, const MatrixType& U) - { - typedef typename MatrixType::Index Index; - const Index size = m.cols(); - - for (Index i=0; i < size; ++i) { - if (i == size - 1 || T.coeff(i+1,i) == 0) - T.coeffRef(i,i) = std::abs(T.coeff(i,i)); - else - ++i; - } - m = U * T * U.adjoint(); - } -}; - template void testSingular(MatrixType m, double tol) {