More numerical tweaks.

2 files changed, 12 insertions, 1 deletions

diff --git a/libdimension/dimension/geometry.h b/libdimension/dimension/geometry.h
index 4c4f8b8..6ff9941 100644
--- a/libdimension/dimension/geometry.h
+++ b/libdimension/dimension/geometry.h
@@ -63,7 +63,7 @@ typedef struct dmnsn_bounding_box { dmnsn_vector min, max; } dmnsn_bounding_box;
 
 /* Constants */
 
-#define dmnsn_epsilon 1.0e-9
+#define dmnsn_epsilon 1.0e-10
 
 static const dmnsn_vector dmnsn_zero = { 0.0, 0.0, 0.0 };
 static const dmnsn_vector dmnsn_x    = { 1.0, 0.0, 0.0 };
diff --git a/libdimension/polynomial.c b/libdimension/polynomial.c
index 6e2aed2..ab80f69 100644
--- a/libdimension/polynomial.c
+++ b/libdimension/polynomial.c
@@ -57,6 +57,15 @@ dmnsn_real_degree(const double poly[], size_t degree)
   return 0;
 }
 
+/* Divide each coefficient by the leading coefficient */
+static inline void
+dmnsn_normalize_polynomial(double poly[], size_t degree)
+{
+  for (size_t i = 0; i <= degree; ++i) {
+    poly[i] /= poly[degree];
+  }
+}
+
 /* Eliminate trivial zero roots from poly[] */
 static inline void
 dmnsn_eliminate_zero_roots(double poly[], size_t *degree)
@@ -285,6 +294,8 @@ dmnsn_solve_polynomial(const double poly[], size_t degree, double x[])
 
   /* Account for leading zero coefficients */
   degree = dmnsn_real_degree(p, degree);
+  /* Normalize the leading coefficient to 1.0 */
+  dmnsn_normalize_polynomial(p, degree);
   /* Eliminate simple zero roots */
   dmnsn_eliminate_zero_roots(p, &degree);