From c5ad724ac51240d57ebe7a8e4875374b052a3211 Mon Sep 17 00:00:00 2001 From: SatoshiTerasaki Date: Thu, 12 Dec 2024 21:12:11 +0900 Subject: [PATCH] Fix --- include/sparseir/_root.hpp | 2 +- include/sparseir/poly.hpp | 2 -- test/poly.cxx | 9 ++++++--- 3 files changed, 7 insertions(+), 6 deletions(-) diff --git a/include/sparseir/_root.hpp b/include/sparseir/_root.hpp index be666f9..8f03f04 100644 --- a/include/sparseir/_root.hpp +++ b/include/sparseir/_root.hpp @@ -235,7 +235,7 @@ std::vector discrete_extrema(F f, const std::vector &xgrid) std::vector res; for (size_t i = 0; i < a.size(); ++i) { // abs ∘ f - auto abf = [f](double x) { return std::abs(f(x)); }; + auto abf = [f](double x) { return std::fabs(f(x)); }; res.push_back( bisect_discr_extremum(abf, a[i], b[i], absf_a[i], absf_b[i])); } diff --git a/include/sparseir/poly.hpp b/include/sparseir/poly.hpp index 201abe4..0b4a316 100644 --- a/include/sparseir/poly.hpp +++ b/include/sparseir/poly.hpp @@ -250,8 +250,6 @@ class PiecewiseLegendrePoly { { Eigen::VectorXd grid = this->knots; - std::cout << "grid: " << grid.transpose() << "\n"; - Eigen::VectorXd refined_grid = refine_grid(grid, 2); auto f = [this](double x) { return this->operator()(x); }; // convert to std::vector diff --git a/test/poly.cxx b/test/poly.cxx index 5073047..b789d33 100644 --- a/test/poly.cxx +++ b/test/poly.cxx @@ -329,10 +329,13 @@ TEST_CASE("Roots") REQUIRE(roots.size() == expected_roots.size()); for(Eigen::Index i = 0; i < roots.size(); ++i) { - REQUIRE(std::abs(roots[i] - expected_roots[i]) < 1e-10); - //Verify roots are in domain + REQUIRE(std::fabs(roots[i] - expected_roots[i]) < 1e-10); + // Verifying the polynomial evaluates to zero + //at the roots with tight tolerance + REQUIRE(std::fabs(pwlp(roots[i])) < 1e-10); + // Verify roots are in domain REQUIRE(roots[i] >= knots[0]); - REQUIRE(roots[i] <= knots[knots.size()-1]); + REQUIRE(roots[i] <= knots[knots.size() - 1]); // Verify these are actually roots REQUIRE(std::abs(pwlp(roots[i])) < 1e-10); }