summaryrefslogtreecommitdiffstats
diff options
context:
space:
mode:
-rw-r--r--libdimension/polynomial.c37
1 files changed, 28 insertions, 9 deletions
diff --git a/libdimension/polynomial.c b/libdimension/polynomial.c
index d78361c..77e561d 100644
--- a/libdimension/polynomial.c
+++ b/libdimension/polynomial.c
@@ -148,21 +148,40 @@ dmnsn_descartes_rule(double poly[], size_t degree)
return changes;
}
+#define DMNSN_NBINOM 11
+static const double dmnsn_pascals_triangle[DMNSN_NBINOM][DMNSN_NBINOM] = {
+ { 1.0 },
+ { 1.0, 1.0 },
+ { 1.0, 2.0, 1.0 },
+ { 1.0, 3.0, 3.0, 1.0 },
+ { 1.0, 4.0, 6.0, 4.0, 1.0 },
+ { 1.0, 5.0, 10.0, 10.0, 5.0, 1.0 },
+ { 1.0, 6.0, 15.0, 20.0, 15.0, 6.0, 1.0 },
+ { 1.0, 7.0, 21.0, 35.0, 35.0, 21.0, 7.0, 1.0 },
+ { 1.0, 8.0, 28.0, 56.0, 70.0, 56.0, 28.0, 8.0, 1.0 },
+ { 1.0, 9.0, 36.0, 84.0, 126.0, 126.0, 84.0, 36.0, 9.0, 1.0 },
+ { 1.0, 10.0, 45.0, 120.0, 210.0, 252.0, 210.0, 120.0, 45.0, 10.0, 1.0 }
+};
+
/* Get the (n k) binomial coefficient */
static inline double
dmnsn_binom(size_t n, size_t k)
{
- if (k > n - k) {
- k = n - k;
- }
+ if (n < DMNSN_NBINOM && k < DMNSN_NBINOM) {
+ return dmnsn_pascals_triangle[n][k];
+ } else {
+ if (k > n - k) {
+ k = n - k;
+ }
- double ret = 1.0;
- for (size_t i = 0; i < k; ++i) {
- ret *= n - i;
- ret /= i + 1;
- }
+ double ret = 1.0;
+ for (size_t i = 0; i < k; ++i) {
+ ret *= n - i;
+ ret /= i + 1;
+ }
- return ret;
+ return ret;
+ }
}
/* Find all ranges that contain a single root, with Uspensky's algorithm */