From d90a0534bed2503112286425b2538b85fe1252ec Mon Sep 17 00:00:00 2001 From: Charles Schlosser Date: Mon, 5 Jan 2026 05:19:49 +0000 Subject: [PATCH] fix polynomialsolver test failures libeigen/eigen!2104 --- .../Eigen/src/Polynomials/PolynomialSolver.h | 22 +++++++++---------- unsupported/test/polynomialsolver.cpp | 11 ++++++---- 2 files changed, 18 insertions(+), 15 deletions(-) diff --git a/unsupported/Eigen/src/Polynomials/PolynomialSolver.h b/unsupported/Eigen/src/Polynomials/PolynomialSolver.h index aa357a41a..a79a49ec7 100644 --- a/unsupported/Eigen/src/Polynomials/PolynomialSolver.h +++ b/unsupported/Eigen/src/Polynomials/PolynomialSolver.h @@ -82,8 +82,8 @@ class PolynomialSolverBase { } protected: - template - inline const RootType& selectComplexRoot_withRespectToNorm(squaredNormBinaryPredicate& pred) const { + template + inline const RootType& selectComplexRoot_withRespectToNorm(Predicate& pred) const { Index res = 0; RealScalar norm2 = numext::abs2(m_roots[0]); for (Index i = 1; i < m_roots.size(); ++i) { @@ -114,25 +114,25 @@ class PolynomialSolverBase { } protected: - template + template inline const RealScalar& selectRealRoot_withRespectToAbsRealPart( - squaredRealPartBinaryPredicate& pred, bool& hasArealRoot, + Predicate& pred, bool& hasArealRoot, const RealScalar& absImaginaryThreshold = NumTraits::dummy_precision()) const { using std::abs; hasArealRoot = false; Index res = 0; - RealScalar abs2(0); + RealScalar val(0); for (Index i = 0; i < m_roots.size(); ++i) { if (abs(m_roots[i].imag()) <= absImaginaryThreshold) { if (!hasArealRoot) { hasArealRoot = true; res = i; - abs2 = m_roots[i].real() * m_roots[i].real(); + val = abs(m_roots[i].real()); } else { - const RealScalar currAbs2 = m_roots[i].real() * m_roots[i].real(); - if (pred(currAbs2, abs2)) { - abs2 = currAbs2; + const RealScalar curr = abs(m_roots[i].real()); + if (pred(curr, val)) { + val = curr; res = i; } } @@ -145,9 +145,9 @@ class PolynomialSolverBase { return numext::real_ref(m_roots[res]); } - template + template inline const RealScalar& selectRealRoot_withRespectToRealPart( - RealPartBinaryPredicate& pred, bool& hasArealRoot, + Predicate& pred, bool& hasArealRoot, const RealScalar& absImaginaryThreshold = NumTraits::dummy_precision()) const { using std::abs; hasArealRoot = false; diff --git a/unsupported/test/polynomialsolver.cpp b/unsupported/test/polynomialsolver.cpp index 12d137d24..fc66c4240 100644 --- a/unsupported/test/polynomialsolver.cpp +++ b/unsupported/test/polynomialsolver.cpp @@ -131,28 +131,28 @@ void evalSolverSugarFunction(const POLYNOMIAL& pols, const ROOTS& roots, const R bool hasRealRoot; // Test absGreatestRealRoot - RealScalar r = psolve.absGreatestRealRoot(hasRealRoot); + RealScalar r = psolve.absGreatestRealRoot(hasRealRoot, test_precision()); VERIFY(hasRealRoot == (real_roots.size() > 0)); if (hasRealRoot) { VERIFY(internal::isApprox(real_roots.array().abs().maxCoeff(), abs(r), psPrec)); } // Test absSmallestRealRoot - r = psolve.absSmallestRealRoot(hasRealRoot); + r = psolve.absSmallestRealRoot(hasRealRoot, test_precision()); VERIFY(hasRealRoot == (real_roots.size() > 0)); if (hasRealRoot) { VERIFY(internal::isApprox(real_roots.array().abs().minCoeff(), abs(r), psPrec)); } // Test greatestRealRoot - r = psolve.greatestRealRoot(hasRealRoot); + r = psolve.greatestRealRoot(hasRealRoot, test_precision()); VERIFY(hasRealRoot == (real_roots.size() > 0)); if (hasRealRoot) { VERIFY(internal::isApprox(real_roots.array().maxCoeff(), r, psPrec)); } // Test smallestRealRoot - r = psolve.smallestRealRoot(hasRealRoot); + r = psolve.smallestRealRoot(hasRealRoot, test_precision()); VERIFY(hasRealRoot == (real_roots.size() > 0)); if (hasRealRoot) { VERIFY(internal::isApprox(real_roots.array().minCoeff(), r, psPrec)); @@ -180,6 +180,9 @@ void polynomialsolver(int deg) { cout << "Test sugar" << endl; RealRootsType realRoots = RealRootsType::Random(deg); + // sort by ascending absolute value to mitigate precision lost during polynomial expansion + std::sort(realRoots.begin(), realRoots.end(), + [](RealScalar a, RealScalar b) { return numext::abs(a) < numext::abs(b); }); roots_to_monicPolynomial(realRoots, pols); evalSolverSugarFunction(pols, realRoots.template cast >().eval(), realRoots); }