summaryrefslogtreecommitdiffstats
path: root/src/continued-fraction.c
diff options
context:
space:
mode:
Diffstat (limited to 'src/continued-fraction.c')
-rw-r--r--src/continued-fraction.c69
1 files changed, 69 insertions, 0 deletions
diff --git a/src/continued-fraction.c b/src/continued-fraction.c
new file mode 100644
index 0000000..6cbf8a5
--- /dev/null
+++ b/src/continued-fraction.c
@@ -0,0 +1,69 @@
+/*
+ * Stolen from gnumeric
+ *
+ * Original code from Morten Welinder <terra@gnome.org>
+ */
+
+#include <stdio.h>
+#include <limits.h>
+
+#include "continued-fraction.h"
+
+void sa_continued_fraction(double val, int max_denom, int *res_num, int *res_denom) {
+ int n1, n2, d1, d2;
+ double x, y;
+
+ if (val < 0) {
+ sa_continued_fraction (-val, max_denom, res_num, res_denom);
+ *res_num = -*res_num;
+ return;
+ }
+
+ n1 = 0; d1 = 1;
+ n2 = 1; d2 = 0;
+
+ x = val;
+ y = 1;
+
+ do {
+ int a = (int) (x / y);
+ double newy = x - a * y;
+ int n3, d3;
+
+ if ((n2 && a > (INT_MAX - n1) / n2) ||
+ (d2 && a > (INT_MAX - d1) / d2) ||
+ a * d2 + d1 > max_denom) {
+ *res_num = n2;
+ *res_denom = d2;
+ return;
+ }
+
+ n3 = a * n2 + n1;
+ d3 = a * d2 + d1;
+
+ x = y;
+ y = newy;
+
+ n1 = n2; n2 = n3;
+ d1 = d2; d2 = d3;
+ } while (y > 1e-10);
+
+ *res_num = n2;
+ *res_denom = d2;
+}
+
+#if 0
+
+int main(int argc, char*argv[]) {
+
+ float value = 3.1417235;
+ int res_num, res_denom;
+
+ sa_continued_fraction(value, 0xFFFF, &res_num, &res_denom);
+
+ printf("%g ~= %i/%i = %g\n", value, res_num, res_denom, (float) res_num / res_denom);
+
+ return 0;
+}
+
+#endif