Fix the polynomial solver when the bounds are exact.

1 files changed, 24 insertions, 0 deletions

diff --git a/libdimension/tests/polynomial.c b/libdimension/tests/polynomial.c
index 208edf9..335fd27 100644
--- a/libdimension/tests/polynomial.c
+++ b/libdimension/tests/polynomial.c
@@ -23,6 +23,7 @@
  */
 
 #include "tests.h"
+#include "../polynomial.c"
 
 /* poly[] = 2*(x + 1)*(x - 1.2345)*(x - 2.3456)*(x - 5)*(x - 100) */
 static const double poly[6] = {
@@ -66,3 +67,26 @@ DMNSN_TEST(polynomial, accurate_roots)
     ck_assert(fabs(ev) < DMNSN_CLOSE_ENOUGH);
   }
 }
+
+/* repeated_root[] = (x - 1)^6 */
+static const double repeated_root[7] = {
+  [6] =  1.0,
+  [5] = -6.0,
+  [4] =  15.0,
+  [3] = -20.0,
+  [2] =  15.0,
+  [1] = -6.0,
+  [0] =  1.0,
+};
+
+DMNSN_TEST(stability, equal_bounds)
+{
+  double root = dmnsn_bisect_root(repeated_root, 6, 1.0, 1.0);
+  ck_assert_msg(root == 1.0, "root == %.17g", root);
+}
+
+DMNSN_TEST(stability, equal_values_at_bounds)
+{
+  double root = dmnsn_bisect_root(repeated_root, 6, 0.9, 1.1);
+  ck_assert_msg(fabs(root - 1.0) < dmnsn_epsilon, "root == %.17g", root);
+}