aboutsummaryrefslogtreecommitdiffstats
path: root/research/experiments
diff options
context:
space:
mode:
Diffstat (limited to 'research/experiments')
-rw-r--r--research/experiments/040-against-c/2026-10-05-newton.jsonl3
-rw-r--r--research/experiments/040-against-c/newton.c101
2 files changed, 104 insertions, 0 deletions
diff --git a/research/experiments/040-against-c/2026-10-05-newton.jsonl b/research/experiments/040-against-c/2026-10-05-newton.jsonl
new file mode 100644
index 0000000..fd56bd7
--- /dev/null
+++ b/research/experiments/040-against-c/2026-10-05-newton.jsonl
@@ -0,0 +1,3 @@
+{"moved":7.73,"kept":7.03,"checked":9.77,"int":8.69,"sub":9.58,"after":9.62,"same":1}
+{"moved":7.70,"kept":7.03,"checked":9.73,"int":8.69,"sub":9.57,"after":9.62,"same":1}
+{"moved":7.70,"kept":7.03,"checked":9.73,"int":8.69,"sub":9.59,"after":9.63,"same":1}
diff --git a/research/experiments/040-against-c/newton.c b/research/experiments/040-against-c/newton.c
new file mode 100644
index 0000000..5a01033
--- /dev/null
+++ b/research/experiments/040-against-c/newton.c
@@ -0,0 +1,101 @@
+/* roots' Newton loop (natural.h's c_roots, i to 100,000) in the forms the
+ * tier may make it, as AArch64 assembly (research 124, docs/design/023):
+ * moved next apart from r, moved to r on the way back (the tier's)
+ * kept next in r's register, the old r kept apart for the compare (gcc's)
+ * checked moved, next checked finite as numbers (fabs, fcmp with infinity)
+ * int moved, next checked finite in integers (its exponent's bits)
+ * sub moved, next - next checked 0 (not a NaN)
+ * after moved, the check after the compare
+ * Each sums r and the steps (the answers must be equal), three rounds. */
+#define _POSIX_C_SOURCE 200809L
+#include <stdint.h>
+#include <stdio.h>
+#include <time.h>
+static double now(void) { struct timespec t; clock_gettime(CLOCK_MONOTONIC, &t); return t.tv_sec * 1e3 + t.tv_nsec / 1e6; }
+/* Newton for each i: r = i; repeat next = (r + i/r) * 0.5; if next >= r break; r = next */
+static double moved(uint32_t n) { /* r and next apart: fmov on the way back */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul d2, d1, %d[half]\n fcmp d2, %d[r]\n b.ge 2f\n"
+ "add %[c], %[c], #1\n fmov %d[r], d2\n cmp %[c], #64\n b.lo 1b\n 2:\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half) : "d1", "d2", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+static double kept(uint32_t n) { /* next in r's register, the old r kept apart for the compare (as gcc) */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fmov d3, %d[r]\n fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul %d[r], d1, %d[half]\n fcmp %d[r], d3\n b.ge 2f\n"
+ "add %[c], %[c], #1\n cmp %[c], #64\n b.lo 1b\n 2: fmov %d[r], d3\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half) : "d1", "d2", "d3", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+static double checked(uint32_t n) { /* moved, plus a finite check (fabs, fcmp, b.hs) on next */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5, inf = 1.0 / 0.0;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul d2, d1, %d[half]\n fabs d4, d2\n fcmp d4, %d[inf]\n b.hs 2f\n fcmp d2, %d[r]\n b.ge 2f\n"
+ "add %[c], %[c], #1\n fmov %d[r], d2\n cmp %[c], #64\n b.lo 1b\n 2:\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half), [inf] "w"(inf) : "d1", "d2", "d4", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+
+static double checked_int(uint32_t n) { /* the finite check in integers: the exponent's bits all ones */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul d2, d1, %d[half]\n fmov x10, d2\n ubfx x10, x10, #52, #11\n cmp x10, #2047\n b.eq 2f\n fcmp d2, %d[r]\n b.ge 2f\n"
+ "add %[c], %[c], #1\n fmov %d[r], d2\n cmp %[c], #64\n b.lo 1b\n 2:\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half) : "d1", "d2", "x10", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+static double checked_sub(uint32_t n) { /* the finite check as next - next being 0 (not a NaN) */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul d2, d1, %d[half]\n fsub d4, d2, d2\n fcmp d4, #0.0\n b.ne 2f\n fcmp d2, %d[r]\n b.ge 2f\n"
+ "add %[c], %[c], #1\n fmov %d[r], d2\n cmp %[c], #64\n b.lo 1b\n 2:\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half) : "d1", "d2", "d4", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+static double checked_after(uint32_t n) { /* the finite check after the compare (off the way to the next round) */
+ double sum = 0; uint64_t steps = 0;
+ for (uint32_t i = 1; i <= n; ++i) {
+ double r = i, x = i, half = 0.5, inf = 1.0 / 0.0;
+ uint64_t c = 0;
+ __asm__ volatile(
+ "1: fdiv d1, %d[x], %d[r]\n fadd d1, %d[r], d1\n fmul d2, d1, %d[half]\n fcmp d2, %d[r]\n b.ge 3f\n"
+ "add %[c], %[c], #1\n fmov %d[r], d2\n fabs d4, d2\n fcmp d4, %d[inf]\n b.hs 2f\n cmp %[c], #64\n b.lo 1b\n b 2f\n 3: fabs d4, d2\n fcmp d4, %d[inf]\n 2:\n"
+ : [r] "+&w"(r), [c] "+r"(c) : [x] "w"(x), [half] "w"(half), [inf] "w"(inf) : "d1", "d2", "d4", "cc");
+ sum += r, steps += c;
+ }
+ return sum + steps;
+}
+int main(void) {
+ for (int k = 0; k < 3; ++k) {
+ double t0 = now(), a = moved(100000), t1 = now(), b = kept(100000), t2 = now(), c = checked(100000), t3 = now();
+ double d = checked_int(100000), t4 = now(), e = checked_sub(100000), t5 = now(), f = checked_after(100000), t6 = now();
+ printf("{\"moved\":%.2f,\"kept\":%.2f,\"checked\":%.2f,\"int\":%.2f,\"sub\":%.2f,\"after\":%.2f,\"same\":%d}\n", t1 - t0, t2 - t1, t3 - t2, t4 - t3, t5 - t4, t6 - t5, a == b && b == c && c == d && d == e && e == f); if (0) printf("moved %.2f kept %.2f checked %.2f int %.2f sub %.2f after %.2f ms (%d)\n", t1 - t0, t2 - t1, t3 - t2, t4 - t3, t5 - t4, t6 - t5, a == b && b == c && c == d && d == e && e == f);
+ }
+}